The Potential Impact of Climate Extremes on Cotton and Wheat Crops in Southern Punjab, Pakistan

The assessment of climate extremes’ impact on crop yield is essential to improve our understanding of agricultural resilience. In the present study, we analyzed the potential impact of climate extremes on wheat and cotton production in Southern Punjab, Pakistan using 30-year observed data from the Pakistan Meteorological Department (PMD) and the fifth-generation reanalysis data (ERA-5) from the European Center for Medium-Range Weather Forecasts (ECMWF). Cotton is a Kharif season crop that is sown in May and harvested in October, and wheat is a Rabi season crop that is planted in November and harvested in April. The agricultural data (1985–2015) that contained the crop area and crop yield were obtained from the Bureau of Statistics, Punjab for six selected districts in Southern Punjab. Three precipitation indices, namely consecutive dry days (CDD), consecutive wet days (CWD) and total precipitation of wet days (PRCPTOT), and four temperature indices, namely warm days (TX90p), warm nights (TN90p), cool days (TX10p) and cool nights (TN10p), were selected to analyze the potential impacts of climate extremes on crop production. (1) We found a potential association of TX10p, TN10p, TX90p and TN90p with crop yield in those years for which the production area remained the same. (2) In a few districts of the study area, the wheat yield losses in the Rabi season were associated with an increase in warmer days and warmer nights. (3) The grain size was suppressed due to an increase in the frequency of TX90p and TN90p, which ultimately reduced the net crop production. (4) In some districts, we found strong positive correlations between extreme temperature indices and crop yield; however, other potential factors such as the use of advanced technology, fertilizer, seeds, etc., may lead to improved net production. This study can help in adaptation planning for resilient agricultural production under the stress of climate extreme events in Southern Punjab.


Introduction
It is expected that extreme events will increase worldwide in the greenhouse gas (GHG)-induced environment. The extreme events will potentially affect agricultural productivity-for instance, in the form of heatwaves and droughts. The adverse impact of climate extremes on crop yield has implications for food security and the livelihoods of communities. The changes in the intensity and frequency of precipitation extremes substantially affect different stages of crop growth development and yield [1]. The rising trends in regional and global temperature are the key factors that affect crop productivity during the crop plantation, maturity and harvesting stages. A temperature rise can suppress the filling duration and sterility of grain and ultimately reduce the net crop yield  [2]. From a climate change perspective, it is important to understand the impact of climate extremes on the crop to optimize and secure the yield.
The South Asian countries are facing numerous challenges due to increasing trends in floods, droughts and heatwaves, and their likely impacts on the agricultural sector [3,4]. The fifth assessment report of the Inter-Governmental Panel on Climate Change (IPCC) reveals that rainfall events are characterized by strong variability, with both increasing and decreasing trends in different parts of Asia [5]. The agricultural-based economy of the South Asian countries mostly depends on the Indian summer monsoon rainfall. Climate change and related extremes are a severe threat to crop farming in South Asia, including Pakistan, India, Bangladesh and Nepal [6]. It is worth mentioning that the projected changes in the monsoon rainfall patterns could potentially impact the agricultural productivity of the whole region.
Due to its diverse climatology and complex topographical feature, Pakistan is more vulnerable to extreme hydrometeorological events such as heatwaves, droughts and floods [7]. The agriculture sector is the backbone of Pakistan's economy and is facing great challenges due to climate change impacts and extreme events [8][9][10]. Variations in rainfall patterns and changes in the diurnal temperature, timing of sowing and harvesting, availability of water and evapotranspiration losses could potentially impact Pakistan's crop productivity [11]. Pakistan has been ranked 12th by the Global Climate Risk Index as a highly affected country due to climate change and extreme events [12].
Punjab is an industrial and agricultural keystone of Pakistan. The agriculture sector contributes approximately 21.4% to the GDP of the agro-based economy of Pakistan and provides 62% of the livelihood of the population in the rural areas [13,14]. It produces 56.1-61.5% of the total agricultural products and approximately 49% of the labor is engaged in the agriculture sector [15]. Though the agricultural sector of Punjab has a substantial contribution to Pakistan's overall economy, it is under the stress of major climate change-induced impacts in the form of floods and droughts [13]. Wheat and cotton are two major crops grown in Punjab and have a vital contribution to Pakistan's economy. Cotton, also called white gold, is a cash crop in Pakistan and contributes the largest share in export revenues [16]. Wheat is also one of the major crops and accounts for 2.2% of the GDP and approximately 10.3% of the agrarian economy of Pakistan [17].
In the past decade, various studies have been conducted to assess and model the extreme events in Pakistan. For example, Zahid and Rasul [18] studied the potential impact of heatwaves in different provinces of Pakistan. They showed that Punjab, Sindh and Balochistan are the most vulnerable to heatwave-related climate extremes. Ikram et al. [19] studied past and future rainfall extreme events using the global climate model (GCM) data from the Coupled Model Inter-Comparison Phase 5 (CMIP5) and inferred that the lower Punjab, Khyber Pakhtunkhwa and the upper Balochistan areas will face increasing precipitation trends under different representative conservative pathways (RCPs). Abbas et al. [20] concluded that the frequency of extreme precipitation events has already been increased in Sindh, Northern Areas, Azad Jammu and Kashmir and Balochistan. They inferred that the southern parts of Pakistan have experienced more wet spells than dry spells in the past few decades.
Raza and Ahmad [21] provided a detailed analysis of climate change's impact on cotton and wheat crops in the Sindh and Punjab provinces. They analyzed the effect of mean variability in temperature and precipitation on cotton without considering the impact of extreme events. Abid et al. [9] collected household farm data from three districts of Punjab to examine how farmers adapt to and perceive farming under climate change scenarios. Rehman et al. [16] studied the economic perspectives of Pakistan for major field crops and determined that the output of rice and wheat has a significant positive relationship with the agricultural GDP.
Ahmed et al. [22] adopted a modeling approach to simulate different phenological stages such as leaf area index, maturity and flowering days of wheat in Pakistan under rainfed conditions. Aslam [23] found that Pakistan's agricultural productivity (cotton, wheat, rice, sugarcane and maize) is below the potential yields due to environmental, technological and socio-economic constraints. To mitigate the potential impact of climate change, Ahmad et al. [24] assessed the relationship between climate change and phenological changes in the rice-wheat system in Punjab, Pakistan. However, the authors did not investigate the quantitative effect of climate variability, particularly the extreme events, on crop yield. Khan et al. [10] examined the future losses in the agricultural productivity of Pakistan by using the global economic modeling approach. Akbar and Gheewala [25] analyzed the effect of climate change on cotton and sugarcane yield by employing two regional climate models and by considering seasonal variability in mean temperature and precipitation. Ullah et al. [26] evaluated the drought characteristics in Pakistan by using PMD station and reanalysis datasets.
Some of the above-mentioned studies used only station data that do not cover the entire study area due to the sparsity of the data in the region. Moreover, in some other studies, the authors used mean climate variability from the global/regional model instead of reliable extreme indices from the reanalysis of observation. The impact of the climate extremes on crop production is a certain gap in the previous studies on which we focused. The previous studies do not include the effect of extreme climate events on crop yield over the study area. Therefore, to address these gaps, the current study aims to discuss the effect of extreme indices on two major crops in Southern Punjab, which will provide preliminary information for researchers to perform more detailed studies.
It is crucial to assess the effect of extreme events on crop yield for food security and the economic sustainability of resilience agriculture [27]. In our study, we present a comprehensive analysis of climate extremes' impact on wheat and cotton crops in six selected districts of Southern Punjab by using in-situ data from the Pakistan Meteorological Department (PMD) and ERA5 reanalysis for data-sparse regions. The objectives of the study include the analysis of trends in the indices of climate extreme events and the evaluation of their effects on the productivity of wheat and cotton from 1986 to 2015. The results from the present study will assist the farmers, investors and policymakers in the agricultural sector to mitigate and adapt to future climate extreme events in Southern Punjab, Pakistan.

Study Area
In the present research, we focused on six districts of Southern Punjab, namely Dera Ghazi Khan, Layyah, Muzaffargarh, Bahawalpur, Rahimyar Khan and Multan (Figure 1), to analyze the impact of climate extremes on the two major crops, i.e., wheat and cotton. Dera Ghazi Khan is located in the strip between the Indus River and the range of Koh-Suleman mountains. The Indus River flows along the eastern side of Muzaffargarh and Layyah, whereas Dera Ghazi Khan lies on the west. Multan and Bahawalpur are at the eastern bank of the Chenab River and Rahimyar Khan is on the eastern edge of the Indus River. The Dera Ghazi Khan and Muzaffargarh districts are most vulnerable to flash floods caused by torrential rainfall. The study area comprises a diverse landscape. For example, Layyah and Multan have deserts with a hot climate and lie in Sindh Sagar Doab, a region between the Indus and Jhelum Rivers. Southern Punjab is the most vulnerable region to climate extreme events and has experienced drought hazards and flood episodes in the past few decades [28].

Datasets
In this study, we used daily observed PMD and ERA reanalysis datasets for temperature and precipitation. Daily observed data from 1985 to 2015 were acquired for three stations, i.e., Multan, Rahimyar Khan and Bahawalpur, from the PMD. We considered two reanalysis datasets (i.e., ERA5 and MERRA-2) for three data-sparse districts, i.e., Muzaffargarh, Layyah and Dera Ghazi Khan. Two major crop seasons were selected in the six selected districts of Southern Punjab. The first season was Kharif, in which cotton sowing begins in May and is harvested in October. The second one was Rabi, in which wheat sowing starts in November and is harvested in April.
The reanalysis datasets are important tools to assess the past global as well as regional climate conditions in various research fields [29]. The European Centre for Medium-Range Weather Forecasts (ECMWF) fifth-generation reanalysis dataset (hereafter ERA5) provides hourly estimations for different land, ocean and atmospheric climate variables. ERA5 covers 137 atmospheric vertical levels from the surface to 80 km height with 30 km horizontal resolution. Based on advanced data assimilation and modeling techniques, ERA5 aggregates enormous historic observations into global approximations.
The second dataset that was used for data-sparse districts was the Modern-Era Retrospective analysis for Research and Applications (hereafter MERRA-2). MERRA-2 is a National Aeronautics and Space Administration (NASA) space-based long-term global reanalysis dataset started in 1980. MERRA-2 has an upgraded data assimilation system of the Goddard Earth Observing System Model, version 5 (GEOS-5) [30]. Two major crops (cotton and wheat) were considered for this study, owing to their importance in Central and Southern Punjab. Agriculture data contained the crop area and crop yield as measured by the Bureau of Statistics, Punjab for the six selected districts, from 1985 to 2015 (Table 1).

Methodology
Meteorological and agricultural data were acquired to assess the likely impacts of climate extremes on the crops in Southern Punjab. The proximity assessment was carried out using the linear regression model for the reanalysis and PMD datasets. The purpose was to determine the best-suited reanalysis of observation from ERA5 and MERRA-2 for the three data-sparse districts. The linear regression model was used to develop a relationship for two variables and to fit the observed data in a linear equation.
To assess extreme climatic changes at regional and global scales, the World Meteorological Organization (WMO) and World Climate Research Program (WCRP) have jointly established the Expert Team on Climate Change Detection and Indices (ETCCDI) and defined 28 representative climate indices [31]. After the post-processing, we calculated the extreme indices for daily precipitation and temperature (maximum and minimum) datasets from PMD and ERA-5 through R-ClimDex for the six districts of Southern Punjab. The RClimDex library provides a Graphical User Interface (GUI) in R. It was developed by the Climate Research Branch of the Canada Meteorological Services. RClimDex provides data quality assessment and calculation of all 28 core climate extreme indices for precipitation and temperature [32,33].
Before the calculation of core indices, the homogeneity tests and quality control procedures were adopted in the R-ClimDex software for station and reanalysis datasets. Data Quality Control (QC) is a prerequisite for climate index calculations. The unreasonable values include a) daily precipitation amounts less than zero and b) daily maximum temperatures less than the daily minimum temperature. In addition, outliers in daily minimum and maximum temperature were identified [31]. The outliers are daily values outside a defined region as the mean minus or plus n times the standard deviation of the daily value, i.e., mean-n * std, mean + n * std. Here, std represents the standard deviation for the day, n is input by the user and the mean is computed from the climatology of the day [31]. Once the QC was completed, the data were checked for unreasonable values and outliers in daily temperature and precipitation. Moreover, a bootstrap procedure based on the Monte Carlo simulation experiment was adopted to avoid inhomogeneity in percentile-based indices of temperature and precipitation [31]. The bootstrap procedure produces much more homogeneous estimates for time series of exceedance rate.
Index calculation involves threshold estimation and base period extreme index calculation. The procedure is based on empirical quantile estimation [34,35]. The quantile of a distribution is defined as The empirical quantile is set to the smallest or largest value in the sample when j < 1 or j > n, respectively. In other words, quantile estimates corresponding to p < 1/(n + 1) are set to the smallest value in the sample, and those corresponding to p > n/(n + 1) are set to the largest value in the sample [31].
After the calculation of extreme indices, we additionally used a non-parametric Mann-Kendall test for trend analysis of each district for data that were less sensitive to outliers and did not follow a normal distribution [36]. The Z value was used to find a trend's direction. A Z value less than zero indicates a downward trend and Z values greater than zero indicate an upward trend. A higher significance in the trend indicates a higher Z value. Absolute values of Z that are greater than 1.96 show that the trend is significant at the α = 0.05 level. Based on the Mann-Kendall test, a Z value greater than 3.30 indicates significance at α = 0.001 and Z values greater than 2.58 show significance at α = 0.01 [37]. Moreover, a linear regression was used to analyze the specific relationships between two or more variables by fitting a linear equation to the observed data.
To assess the impact of extreme climate events on two seasonal crops, we considered three precipitation and four temperature indices ( Table 2). The extreme temperature indices, such as warm days (TX90p), cool days (TX10p), warm nights (TN90p) and cool nights (TN10p), described the number of days when temperatures exceeded the threshold values. In the case of precipitation indices, the annual total wet-day precipitation (PRCPTOT), consecutive wet days (CWD) and consecutive dry days (CDD) were used to analyze the impact on crop yield. The comparison of reanalysis datasets with observed data using a simple linear regression model is given in Table 3.
Threshold criteria were carefully chosen to calculate required indices from daily precipitation and temperature (max and min) data for both seasons and all districts in Southern Punjab ( Table 4). The selected threshold values are explained in Table 4. The whole methodology is described schematically in Figure 2. The variations in seasonal crop yield were assessed through the seasonal mean variation in temperature, precipitation and area of production. Finally, we analyzed the likely effects of seven extreme climate indices on the variation in crop yield by considering the constant area of production in all six districts. To analyze the particular relationships between two or more variables, a linear regression analysis was used by fitting a linear equation to the observed data. Pearson correlation was used to investigate the link between extreme indices and crop yield. The values of the correlation coefficient were used to assess the level of significance in each district.
The Inverse Distance Weightage (IDW) method was used for interpolation and distribution of trends in extremes. The IDW method is a highly adaptable interpolation method used for algorithm optimization, spatial data and image interpolation [38]. Furthermore, IDW was used for the spatial distribution of the link between extreme indices and crop yield.

Comparison of Observed and Reanalysis Datasets
The reanalysis datasets (ERA-5 and MERRA-2) were compared with observed PMD station data for temperature (maximum and minimum) and precipitation using a simple linear regression model for three districts: Bahawalpur, Multan and Rahimyar Khan. The purpose was to determine the best-suited reanalysis of observation for the other three data-sparse districts, i.e., Muzaffargarh, Layyah and Dera Ghazi Khan. For each district, a simple linear regression model was used to explore the regression between reanalysis datasets and the independent variable of observed data. The comparison of reanalysis with observed data is given in Table 3. Based on our regression analysis (Table 3), we found ERA-5 to be the best-suited reanalysis dataset relative to MERRA-2 for the data-sparse areas. Though ERA-5 and MERRA-2 show good agreement for temperature, the R 2 values for precipitation are higher for ERA-5 as compared to MERRA-2, as shown in Table 3. In a few previous studies, ERA-5 has shown enormous progress and performance as compared to other available reanalysis datasets for data-sparse regions [29,39,40], particularly over western and southwestern parts of Pakistan [24]

Temperature-and Precipitation-Related Extreme Climate Indices
Three precipitation and four temperature indices (Table 3) from ETCCDI were used to analyze the climatic extreme patterns over six selected districts of the study area ( Figure  1). An increase in warm days (TX > 90th percentile) and warm nights (TN > 90th percentile) was observed in most of the districts of the study area in the Rabi season. On the other hand, a decreasing trend in the cool days (TX < 10th percentile) and cool nights (TN < 10th percentile) was observed in both seasons and for all selected districts, except Muzaffargarh and Multan. Choi et al. [41] also revealed an increase in the number of warm days/nights and decrease in the cool days/nights in Southeast Asia and the South Pacific. In the Kharif season, all districts presented a decreasing trend in warm days as well as in warm nights, as shown in Figures 3-5.    Over the Multan district, a significant increasing trend in warm nights (TN90P) was observed in the Kharif and Rabi seasons (Figures 3, 6a and 7a). On the other hand, a highly significant increasing trend in CWD was observed at a significance level of 0.05 in Rabi (Figure 4) over the Bahawalpur district (Figures 6f and 7f). Nevertheless, a significant decreasing trend in cool nights (TN10P) was observed over the Multan district in the Rabi season, as shown in Figure 7e. Our analysis revealed a decreasing trend in TN10P with a significance level of 0.05 over the Bahawalpur district in the Kharif season (Figure 6d). Over the Rahimyar Khan district, TX90P was decreasing significantly in the Kharif season ( Figure 6b). However, we observed a significant increasing trend in TN90P and significant decreasing trends in TN10P for both Kharif and Rabi seasons over the Rahimyar Khan district (Figures 3, 4, 6 and 7).

Seasonal Variations in Temperature, Precipitation and Crop Yield
Seasonal averages of daily precipitation and temperature (max and min) in the six districts were used to analyze the inter-annual variability of seasonal crop production. For brevity, we present only the results of three districts, i.e., Multan, Bahawalpur, Muzaffargarh.
From 1985 to 2000, the area, as well as the net Kharif yield, showed an upward trend in the Bahawalpur district, as shown in Figure 8. It is interesting to note that, after 2000, prominent variations in the area were observed in the seasonal crop yield under the constant area of production. Similar effects were observed for the cotton crop due to seasonal variation in precipitation and temperature in 2011 and 2014 in the Bahawalpur district. Significant upward and downward trends in the area of production were observed in the Muzaffargarh district ( Figure 9). Initially, the area increased from 1985 to 2005 and then decreased until 2015. For the latter, we conclude that ups and downs in the area of production are directly proportional to the crop yield. Hence, our analysis revealed that a decrease in the seasonal maximum temperature and an increase in seasonal precipitation were likely associated with a decrease in the cotton yield in Multan. Similarly, the area of production was the same from 2003 to 2007, with considerable changes in the Kharif season yield in the Multan district. We noted that, in 2003, the net production was 133.4 tons/ha, with daily average precipitation of 4.9 mm/day. In 2004, the net cotton production increased to 184 tons/ha, with a decrease in precipitation of up to 1.2 mm/day. An increase in yield from 143.3 to 177.8 tons /ha in 2005 and 2006 was probably due to an increase in the minimum temperature ( Figure 10).
The area of production was slightly increased in the first five years of the study period. From 1996 to 1998, the area of production decreased and very few changes were observed for the remaining years. Linear fitted analysis for maximum and minimum temperatures showed a constant trend line with very little inter-annual variability. From 1988 to 1990, the area of production remained the same, but considerable fluctuation was observed in net production in the Multan district. The production decreased from 284 tons/ha in 1988 to 198.8 tons/ha in 1989. It is clear from Figure 10 that the seasonal average of daily average precipitation increased from 4.6 to 7.2 mm/day in 1988 and 1999. The linear fitted trend analysis showed that the seasonal average of daily precipitation decreased from 1985 to 2015 in Multan for the Kharif season (see Figure 10). However, the inter-annual variations were prominent throughout the study period.

Impact of Climate Extremes on Crop Yield
In addition to mean climate variability, we also assessed the potential impact of extreme indices on crop yield for both seasons. The results of the Z value from the Mann-Kendall test and correlation coefficient from Pearson's correlation were used to assess the spatial distribution of the correlation coefficient between extreme indices and crop yield. We found that four temperature indices had a significant correlation with crop yield over the study area. Therefore, for brevity, we only present the spatial map between temperature extremes and crop yield in the following sections.

Kharif Season
In the Kharif season, the crop yield showed numerous fluctuations in all districts with changes in the area of production. As mentioned earlier, the area of production is the primary cause of variability in crop yield. However, in some cases, the variation in crop yield was observed with an unchanged area of production. This means that secondary factors may have affected the net crop production. Here, our aim was to associate the impact of climate extremes with the crop yield variation for the Kharif season by considering those years in which the area of production was constant.
We noted a significant negative correlation of warm days (TX90p) and warm nights (TN90P) with crop yield over the districts of Layyah, Dera Ghazi Khan and Bahawalpur (Figure 11a,b). For instance, a decrease in TX90p and TN90p frequency was associated with an increase in cotton production from 25.2 to 28 tons/ha in the Layyah district (Supplementary Figure S1). On the other hand, a significant positive correlation was observed between cotton yield and TX90P over Rahimyar Khan (Figure 11b). This means that the net crop yield as well as the TX90P increased during the 30-year period. A higher number of warmer days could potentially impact the cotton yield.
However, Silvertooth [42] revealed that, regardless of whether the maximum daytime temperature increases or decreases, the yield depends on the soil moisture availability and the stage of crop development. Other potential factors, such as the use of advanced technology, fertilizer, seeds, etc., may lead to improved crop production. Regarding the correlation of cool days (TX10P) and cool nights (TN10P) with cotton yield, we found a significant positive correlation over Multan and Bhawalpur. On the other hand, a significant negative correlation was found over Layyah. We observed an overall increase in cotton production with the decreased numbers of TX10p and TN10p, while the area of production remained constant. Remarkably, we found a net increase in cotton production with a positive correlation with TX10P over Multan and TN10P over Bahawalpur ( Figure  11c,d).

Rabi Season
In the Rabi season, the wheat production either increased or decreased with a constant area of production in several districts of the study area, owing to the likely impact of extreme temperature indices. As mentioned earlier, there is a range of factors related to crop production. However, climate extremes are most influential as compared to other factors that could affect crop yield. Moreover, in our analysis, we only considered those years in which the area production remained the same. For example, the area of production was constant in the Bahawalpur district, and an upward trend was found for the frequency of TX90p and TN90p (Supplementary Figure S2). This upward trend was associated with a decline in wheat yield from 161.5 to 136.7 tons/ha (Supplementary Figure S2). This confirms a strong positive correlation between wheat yield and TN90P (Figure 11d). Similar results were found over the Multan district for TX90P (Supplementary Figure S3). In contrast to TX90P, a significant positive correlation was observed for TN90P over the Multan and Dera Ghazi Khan districts (Figure 12a,b). In the Rahimyar Khan district, the wheat yield declined, with a strong positive correlation with TN10P (Figure 12c). Very similar behavior of TN10P was observed for the Layyah district. As with the Kharif season, the Rabi season crop yield was found to be positively associated with an increase in TX10P in the Bahawalpur district. In this study, we considered only two factors that could affect the crop yield, i.e., the area of production and climate extremes. Other important factors such as soil fertility, availability of water and diseases or pests must also be considered in future studies to provide a complete evaluation of the causes of yield losses and resilience of agricultural production.

Discussion
The climate extreme trends were analyzed by using R-ClimDex for 1985-2015 using ERA-5 and PMD data. In general, we found that the intensity and frequency of precipitation extremes have increased in the past 30 years. A positive trend was observed for CWD and PRCPTOT in all districts (Figures 3-5). Klein Tank and Können [43] reported an increase in heavy precipitation events from 1946 to 1999 over the European states. We observed that the consecutive dry days (CDD) either remained stable or slightly increased in a few districts from 1985 to 2015. The changes in CDD length provide vital information about drought conditions [44]. The number of CDD is an effective measure of seasonal drought and extreme precipitation. The increasing trend in CDD may indicate greater risk and higher dryness for drought in all seasons [45]. However, the precipitation-related extremes had a less significant impact on crop yield since the precipitation data for ERA5 were least correlated with observations as compared to temperature.
We also analyzed the inter-annual variability in daily mean precipitation and temperature and their likely impact on crop yield over the six districts of the study area. Our focus was to consider the variability in crop yield for those years having a constant area of crop production. We found a close association of cotton yield with the seasonal maximum and minimum temperature. Noteworthy increasing and decreasing trends in the area of production were observed in a few districts that were directly proportional to crop yield. For example, it is clear from Figure 9 that the area of production was significantly reduced after 2010 and the net crop production was also decreased for the Muzaffargarh district. We concluded that the overall net crop production followed the trend line of fluctuation in the area. Therefore, the seasonal changes in temperature or precipitation cannot be considered as a single entity that could affect the yield of cotton or wheat.
There is a range of factors related to crop production. The crop yield mainly depends on the area of production. However, some other important factors could affect the crop yield, e.g., the availability of water, soil fertility, diseases or pests and climate extremes. Irrespective of the behavior of other factors, we only focused on climate extremes and analyzed their likely effects on wheat production for the past thirty years in Southern Punjab. The most important overlooked factor associated with crop yield is climate extremes. Our results suggest that the precipitation-related extreme indices may have very little importance across both crops. Therefore, we only focused on the impact of temperature-related extreme indices for Kharif and Rabi season crops. According to Reddy et al. [46], the temperature is a basic controlling factor for cotton development and growth rate. The frequent occurrence of optimal temperatures in many cotton belt areas can alter the growth and development of the cotton crop. Vogel et al. [47] also revealed that crop yield was negatively affected due to an increase in the number of unusually warm days.
Extreme temperature indices-more importantly, a higher number of warm dayscan cause a serious reduction in cotton yield. The growing number of warm days can affect the physiology of cotton and potentially impact the net yield [42]. Similarly, the freezing temperatures and increased number of cool nights may also affect the seedling establishment and germination percentage [42]. The cotton crop grows more slowly when TX10p and TN10p reach high values for the whole season [48]. This means that the colder than usual temperatures may reduce the crop yields because the growth of the plant depends on exposure to accumulated temperatures.
We also analyzed the association of the frequency of TX90p and TN90p with wheat crop yield in all six districts. The fact of the matter is that an increase in the frequency of TX90p and TN90p may suppress the grain size, reduce its quality and ultimately decrease the yield of the wheat crop. Heat stress can strongly affect different phases of crop development and, consequently, reduce the net yield [49]. Giménez et al. [50] revealed that warmer nights may reduce the grain yield and size in a critical period due to a weakened grain number per unit area, which is related to a lower spike number. Cotton uses more water during dry and hot days than humid or cold days. The availability of abundant water gives the cotton a capacity to function and live under very high temperatures. Therefore, although the temperature is an important factor, the relative humidity, soil moisture and air movement are also vital for crop yield. For crop yield improvement, there is a need for a further understanding of the biochemical and molecular basis of heat stress under a warmer environment [51].

Limitations
There is a range of important factors that could affect crop yield-for example, the availability of water, soil fertility, diseases, pests and climate extremes. In the present study, we only considered the impact of extreme climate events on two major crop yields in Southern Punjab. The green revolution in Pakistan started in 1960 and has had some implications for crop yields. However, there is a lack of consistent literature on the green revolution of Pakistan. In a few previous studies, it was reported that the adoption of new high yielding varieties (HYV) increased the food grain yield from 1960 to 1990 [52]. Nevertheless, we noted that, in a few districts, production decreased between 1985 and 2015 (Figures 8-10).
We used ERA-5 for the data-sparse districts of the study area from 1985 to 2015. The ensemble of data assimilation (EDA) techniques addresses a few uncertainties in the ERA-5 modeling system, e.g., in the observation and sea surface temperature. However, some uncertainties in ERA-5, such as systematic errors and radiative forcing of the model, are not accounted for. Moreover, the current study is limited to present-day climatology and does not include future extreme events and their potential impact on crop yield. Additionally, we assumed that the technological development and agronomic practices remained the same during the study period. Regardless of the limitations, the present study can provide a basis for better planning, policy-making and adaptation strategies under the threat of near-future extreme climate events and their potential impact on major field crops. This study also provides preliminary information for researchers to perform more detailed studies by considering the aforementioned limitations.

Conclusions
It is crucial to improve our understanding of the impacts of extreme events on crop yields to enhance the resilience of agriculture systems. The current study was intended to assess the effects of climate extreme indices on two seasonal crops (i.e., wheat and cotton) in Southern Punjab for the period of 1985 to 2015. We evaluated reanalysis datasets with in-situ data from the PMD by using a linear regression model in three districts (i.e., Multan, Bahawalpur and Rahimyar Khan) and found ERA-5 to be a more suitable dataset than MERRA-2 for this region. Our results suggest that the precipitation-related indices have very little importance across both crops. Therefore, we only focused on the impact of temperature-related extreme indices on Kharif and Rabi season crops over the six districts in Southern Punjab. It was observed that the frequency of warm days (warm nights) was increased and the frequency of cool days (cool nights) was decreased in the past 30 years in most of the districts. Our results showed an increasing trend in annual total precipitation and consecutive wet days in all districts.
An increase in maximum temperature showed a negative effect on both crops. In contrast, an increase in minimum temperature was probably associated with an increase in the net crop yield. An increasing trend in both crops' production was found from 1985 to 2015 in a few districts in the study area. However, seasonal changes in temperature and precipitation could not be considered the only factor that affected the production. For instance, there are several other factors, such as technological advancements and an increase in fertilizer usage, that may have significantly boosted the net yield. The potential impact of selected extreme indices was observed on the crops by considering a constant area of production. We concluded that although the net crop production increased in some districts, a considerable decline in cotton and wheat yield was significantly correlated with temperature extremes. The substantial increase in TX90p and TN90p may have suppressed the gain in the grain size and reduced the crop yield. This study was limited to a 30-year reference period and will help to identify the most relevant climate extreme factors for Rabi and Kharif season crops, as well as aiding in the selection of proper adaptation strategies to prevent substantial yield losses in all districts of Southern Punjab. Results may be extended in future studies by considering the projection of climate extremes from the latest climate model scenarios. Funding: The authors did not receive grants or other financial support from any organization for the submitted work.