Climate resilience of the top ten wheat producers in the Mediterranean and the Middle East

Wheat is the main staple crop and an important commodity in the Mediterranean and the Middle East. These are among the few areas in the world where the climate is suitable for growing durum wheat but also are among the most rapidly warming ones, according to the available scenarios of climate projections. How much food security and market stability in the Mediterranean and the Middle East, both depending on wheat production and its interannual variability, are going to be compromised by global warming is an overarching question. To contribute in addressing it, we use a recently established indicator to quantify crop production climate resilience. We present a methodological framework allowing to compute the annual production resilience indicator from nonstationary time series. We apply this approach on the wheat production of the 10 most important producers in the Mediterranean and the Middle East. Our findings shows that if no adaptation will take place, wheat production reliability in the Mediterranean and the Middle East will be threatened by climate change already at 1.5 °C global warming. Average climate-related wheat production losses will exceed the worst past event even if the 2 °C mitigation target is met. These results call for urgent action on adaptation to climate change and support further efforts for mitigation, fully consistently with the Paris Agreement recommendations.


Introduction
Wheat originated approximately 10,000 years ago in the Fertile Crescent in Western Asia, and it was the main staple crop that accompanied the birth of Western civilization there (Preece et al. 2017) and in the Mediterranean, where the structural importance of wheat in building social resilience has been demonstrated (White 2011;Xoplaki et al. 2018). Durum wheat is a species of wheat mainly grown in the Mediterranean region nowadays (Rharrabti et al. 2001;Nazco et al. 2014;Royo et al. 2014;Guzmán et al. 2016;Tidiane Sall et al. 2019).
Durum wheat flour is characterized by higher protein content than soft wheat, stronger gluten, lower glycemic index, and longer durability, which are all essential properties to make pasta ). As such, durum wheat provides important cultural and commercial benefits to the Mediterranean region, but it may be seriously threatened by climate change. Climate model simulations project that the ancient "Fertile Crescent," where durum wheat was originated, will disappear in this century (Kitoh et al. 2008). Furthermore, evidences of recent climate events with serious social implications have been already observed in the Eastern Mediterranean (Kelley et al. 2015). Whether this trend will affect the entire Mediterranean regions is an alarming open question.
The Mediterranean and Middle East wheat production is tightly coupled with climate variability and is largely affected by the occurrence of heat waves and droughts Fontana et al. 2015;Guzmán et al. 2016;Zampieri et al. 2017). Heat waves and droughts are projected to increase in future climate (Dosio et al. 2018;Naumann et al. 2018), threatening the general food production and security (Cramer et al. 2018;Tebaldi and Lobell 2018a, b).

Communicated by Luis Lassaletta
This study addresses wheat production in countries in the Mediterranean region and the Middle East. These regions are usually characterized by sparse natural vegetation (Zampieri and Lionello 2010) and dry summers that are typical of the socalled Mediterranean climate (i.e., Köppen classifications Csa and Csb, warm temperate climate with dry summer; Kottek et al. 2006). We here focus on the top ten wheat producers characterized by Mediterranean climate including Egypt, where irrigation is used to alleviate the effects of a dry climate along the Nile, allowing the wheat cultivation. Despite the scarcity of detailed data on durum wheat production for all the countries in the region, it is well established that durum wheat represents a large proportion of the total wheat grown in these countries Tidiane Sall et al. 2019). Contrarily, soft wheat is the most common wheat species predominantly cultivated in cooler climates.
We quantify the wheat production losses due to drought and heat stress using a statistical method linking the observed climate anomalies to the official production data recorded by the Food and Agriculture Organization (FAO) for the top ten wheat producing countries in the Mediterranean and the Middle East. For these countries, we estimate the climaterelated wheat production losses at different global warming levels using high-resolution climate model simulations. We interpret the results in terms of the recently defined crop production resilience indicator (Zampieri et al. 2019c, d), which we further develop in order to be applicable to the nonstationary time series characterizing observed and projected crop production data.

Data and methods
Wheat production in the Mediterranean and the Middle East is dominated by Turkey with more than 20 million tons (FAOSTAT average 2008(FAOSTAT average -2017Fig. 1). Italythe country that is most famous for pastais ranked at the fourth position, with 7.5 million tons of wheat produced on average in the period 2008-2017. We note that this selection excludes France, which is the largest European wheat producer. This choice is motivated by the fact that French production is mostly located north of the country, where the climate is no longer Mediterranean, and mainly attributable to soft wheat varieties (Ceglar et al. 2016). In Northern France, in fact, wheat production is more sensitive to overwet conditions (Zampieri et al. 2016;Ben-Ari et al. 2018) rather than droughts as typically occurs in the Mediterranean countries growing durum wheat Fontana et al. 2015;Dettori et al. 2017;Dixit et al. 2018).
It is worth noting that the production statistics reported by FAOSTAT do not distinguish among the different wheat species. However, the countries here analyzed are among the main durum wheat producers in the world Tidiane Sall et al. 2019). Therefore, we can assume that the links we find between the reported wheat production and the climate anomalies represented by heat stress and drought are the best current estimates achievable for durum wheat at the country level using a statistical approach.
In order to quantify the links between climate and wheat production variability, we adopt a statistical model defined as a calibrated combination of heat and water stress indicators, i.e., the Combined Stress Index (CSI; Zampieri et al. 2017). The CSI has been widely used to estimate the effects of climate anomalies on crop yields anomalies in current climate (Zampieri et al. 2017(Zampieri et al. , 2019b; it was applied in the context of seasonal forecasting (Ceglar et al. 2018) and in economical modeling (Chatzopoulos et al. 2019). The CSI procedure was further developed to quantify crop optimal production and the climate-related losses from the recorded crop production data available at the country level and to estimate the future production losses according to climate projections (Zampieri et al. 2019a).
Here we apply the most recent version of the CSI procedure computing the heat and water stress indicators over wheat fields (Fig. 2) using climate data derived from observations (Ruane et al. 2015). The CSI is computed during the period of the year when wheat is more sensitive to climate anomalies (Zampieri et al. 2017) and then aggregated at the country level and calibrated with production data from the official FAO statistics for the top ten wheat producers in the Mediterranean and the Middle East in the period 1980-2010 (FAOSTAT, http://www.fao.org/faostat/).
It is worth noting that previous versions of the CSI model (e.g., Zampieri et al. 2017) were based on the links between climate and yield anomalies. However, focusing on yield neglects significant impacts of climate on the harvested area variability (Cohn et al. 2016). Therefore, the version of the CSI model developed for maize (Zampieri et al. 2019a) and adopted here is modified in order to account for the correlations between climate and crop production (i.e., yield multiplied by harvested area) instead of yield. This version of the CSI model explains a large portion of the recorded crop production interannual variability, raising the global r 2 of almost the 10% compared with the overall level that could be achieved considering the yield instead of production (Zampieri et al. 2019a). Moreover, production from different regions can be upscaled on a larger area by simply summing them, which is relevant for the specific purposes of the present paper.
We use a high-resolutions ensemble of climate models simulation (HELIX, high-end climate impacts and extremes) (Dosio et al. 2018;Naumann et al. 2018). The HELIX dataset was obtained by dynamical downscaling of the coarse GCM simulations listed in Tables 1. More specifically, sea surface temperature and sea ice resulting from the GCMs were used as boundary conditions for the atmospheric component of the EC-EARTH model, running at 0.35 degrees resolution. This procedure provides a global meteorological dataset with enough spatial resolution to resolve complex orography and land-sea contrasts that are found in regions such as the Mediterranean area.
The HELIX ensemble considers only a high greenhouse gases emission scenario (RCP8.5). In order to circumvent this limitation and to provide results that are independent from the emission scenario, we organize the analysis in terms of the global warming levels instead of specific future time periods. This strategy allows isolating the impacts of global warming from the climate models' sensitivities to different emission scenarios (Dosio et al. 2018;Naumann et al. 2018). As a result, the specific periods used in this study are based on the model-dependent timings when the different global warming levels are reached in the future (see Table 1). In principle, these periods will occur later for more moderate emission scenarios than the RCP8.5. However, the climate impacts under RCP8.5 diverge significantly from those under the more moderate RCP4.5 scenario only around 2050 (Tebaldi and Friedlingstein 2013). Moreover, the timings when the 1.5 and 2°C global warming levels are reached do not depend significantly on the particular emission scenario, even considering the stronger mitigation levels of the RCP2.6 scenario (Bärring and Strandberg 2018). Therefore, we can assume that our results for the near-future period are realistic albeit based on a single scenario. The method to identify time windows follows the guidelines of the HELIX project . The time windows are centered on the years when the 20-year running mean of global average temperature exceeds 1.5, 2, and 3°C. The timings of the 1.5, 2, and 3°C warming levels are listed in Table 1 and account for a warming of 0.81°C in 2005, compared with period 1880-1900.
The CSI method implemented here estimates the production losses due to drought and heat stress (Zampieri et al. 2019a), which are the predominant yield limitation factors for durum wheat in the Mediterranean regions (Fontana et al. 2015;Dettori et al. 2017;Dixit et al. 2018), and it is also able to estimate the optimal production in the calibration period (Zampieri et al. 2019a). However, future production trends driven by the increase of potential yields due to technological improvement and adaptation are unknown. Furthermore, the physiological effects of elevated CO 2 on crops is largely uncertain as it depends on adaptation as well (Kimball 2016). Therefore, we assume two extreme scenarios of optimal production trends: one where the current trend continues linearly and one where the optimal production (i.e., the production that would be observed with perfect climate conditions) remains constant. Then, we remove the production departures from the optimal yield due to the concurrent climate anomalies estimated by the CSI applied to the climate model simulations. Shifts in the growing season, that may occur under changing climate conditions, are not accounted for. This limitation can be justified by considering that shifts will not be very large (with respect to the temporal framework of the method) and that changes in cultivars would limit the anticipated advance of wheat phenology (Dixit et al. 2018;Rezaei et al. 2018). Other potential forms of adaptations such as shift in growing areas (Ceglar et al. 2019) are not considered in this study. In order to reduce the impacts of these assumptions, we limit the analysis for the time period in the near future when the climate normals are in the range of the observed variability in the calibration period.
In order to quantify the production stability, we use a new indicator combining the mean production with the production variance, i.e., the crop production resilience indicator (Zampieri et al. 2019d). For a stationary time series, the crop production resilience indicator (R C ) is defined as: where μ and σ represent the mean and standard deviation of the production annual values, respectively. This indicator carries the following advantages: 1. It is directly related to the ecological definition of resilience, thus, theoretically more grounded than similar indicators based on different functions of the μ over σ ratio. 2. It is inversely/directly proportional to the frequency/ return period of the extreme events leading to large production losses. 3. It takes into account spatial heterogeneities and diversity in a simple and intuitive manner, i.e., R C computed on the sum of n uncorrelated time series with same μ and σ is exactly n-times R C of the individual time series.
Observed production data often display significant trends and low-frequency variability that are related to technological improvements, changes in cropping areas, and other exogenous factors. Since the crop production resilience indicator definition (Eq. 1) only depends on the ratio between mean production and fluctuation, it is possible to normalize the production time series by the baseline trend or moving window means. This allows for computing the resilience indicator using nonstationary data and implies that there is no change of system resilience if the standard deviation of the crop production anomalies varies proportionally to the baseline mean values.
In order to account for the resilience related to the interannual variability with respect to varying baseline production values, we first compute the smoothed production time series using LOESS procedure (Cleveland and Devlin 1988); then, we compute the normalized production anomalies and the standard deviation: where the p i represents the production values of the time series under evaluation; π i are the normalized anomalies with respect to the baseline values, i.e., the smoothed time series P i ; and σ' is the standard deviation of the normalized anomalies. The nonstationary crop production resilience indicator is simply given by the inverse squared standard deviation of the normalized anomalies: In case the production time series is stationary, R' C is exactly equal to R C . This can be demonstrated mathematically, and it is also proven heuristically by computing the resilience indexes of artificial crop production time series with prescribed statistical properties as the examples depicted in Fig. 3.
Typical situations that can be found with real production data include, among other cases, stationary time series (Fig. 3, blue line), increasing time series (Fig. 3, green line), and increasing and then stagnating time series (Fig. 3, yellow line).
The blue line is a random realization of a stationary production time series with mean equal to 10 and standard deviation equal to 1. The theoretical crop production resilience of such a time series is R C = 100, while the value computed for this particular realization is 89.5. This difference is within the acceptable accuracy range for the R C estimation which is about 15% for 100 years (Zampieri et al. 2019d). The green Fig. 3 Three randomly generated time series representing stationary (blue), linearly increasing (green), and nonlinear (yellow) production data (generic units). Dashed lines represent the smoothed baseline values after the filtering procedure (LOESS) line represents a time series obtained by multiplying data of the blue one for a factor that increases linearly, so the final production is doubled. The size of the anomalies with respect to the smoothed time series, i.e., the dashed line, representing the baseline values, increases of the same factor. The nonlinear crop production resilience value for this time series is R' C = 92.5. Data associated with the yellow line are derived similarly to the ones represented by the green line, just multiplied by a factor that grows faster and then stagnates. The nonlinear crop production resilience value for this time series is R' C = 89.0. Thus, this procedure guarantees that the estimated resilience depends only on the μ over σ running ratio despite the presence of long-term trends and low-frequency variability.

Results
The climate-limited production estimated by the CSI for the top ten wheat producers in the Mediterranean and the Middle East matches the observed values reasonably well (Fig. 4). The nonstationary crop resilience indicator computed on these time series is 216 for the observations and 278 for the CSI estimate. The discrepancy is within the range of uncertainty of the resilience estimation for time series with limited duration, which is about the 30% for 30 years (Zampieri et al. 2019d). This analysis highlights the increasing climate influence on production over the period 1980-2010, represented by the increasing gap between the blue and the red lines shown in Fig. 4. This trend is mainly related to increasing amplitude of heat waves hitting the region (see e.g. Zampieri et al. 2016). The worst observed loss of 21.7 million tons occurred in 2008, and it is well captured by the CSI. The optimal production estimated by the CSI grows approximately linearly until 2005, followed by a slight decrease in the rate. The observed production for the top ten producers displays a stagnation lasting also for the latest years in the records (2011)(2012)(2013)(2014)(2015)(2016)(2017).
A linear trend analysis for the individual countries (Fig. 5) shows positive values for most of the top ten wheat producing countries in the Mediterranean and the Middle East (Fig. 1). However, the trend is stagnating in Iran and Tunisia, while Italy, Syria, and Greece display a recent negative production trend.
1. The analysis of the optimal production for the past (Fig. 4), as well as the recent observed trend of the individual countries (Fig. 5), suggests considering two main overall scenarios for the near-future optimal production: optimal production, linked to technological improvement, is constantly increasing (by filling the so-called yield gap, not increasing cultivated areas). 2. Optimal production is stagnating.
The ensemble average departure of the wheat production levels with respect to the optimal values simulated for the past ( Fig. 6a and b, blue lines) is on average 11 million tons and displays an increasing trend that is consistent with the observed values (Fig. 4). In both scenarios, differing only in the last 5 years, the ensemble mean production averaged between 1980 and 2010 is 58.3 million tons for scenario 1 and 57.6 for scenario 2. The ensemble averaged resilience indicators are R C = 365 for scenario 1 and R C = 336 for scenario 2. The latter one is only slightly overestimating the estimate obtained by using the observed climate data (R C = 278).
The departures of the estimated future wheat production levels with respect to the optimal values continue to increase in the present and near-future climate projections in both scenarios ( Fig. 6a and b, light blue lines). When the optimal production increases linearly (scenario 1), however, the climate-related losses are smaller than the gains, and the balance is positive compared to the earlier period. The interannual variability increases as well, and it is similar between the scenarios because they are based on the same climate data ( Fig. 6c and d). After 2010 but before reaching the 1.5°C global warming level, the ensemble mean production is estimated to be 76 million tons for scenario 1 (R C = 336) and 62.5 million tons for scenario 2 (R C = 223) with a departure from the optimal values of 17.5 million tons in both cases due to heat stress and droughts and a 40% difference between the average resilience indicators (Fig. 6e and f). The ensemble mean computed from the seven climate models over a variable period integrates information of more than 100 years altogether (see Table 1), ensuring a reliable estimation of resilience, with a relative error lower than 15% (Zampieri et al. 2019d).
In the time window between 1.5 and 2°C global warming levels, the departures from the optimal productions surpass the Fig. 4 Sum of top ten wheat production time series in the Mediterranean and the Middle East: estimated optimal production without climate effects in the period 1980-2010 (i.e., heat stress and drought, blue line), estimated production including the effects of climate (red line), and the official records from FAOSTAT in the period 1980-2015 worst event recorded in 2008. The ensemble average of the climate-related losses, consisting of about 22.4 million tons, leads to an average production of 86 million tons for scenario 1 (R C = 387) and 57.6 million tons for scenario 2 (R C = 175). The difference of resilience between the scenarios is 75% during this period, meaning that the frequency of extreme events leading to severe yield loss for the sum of the top ten wheat producers in the Mediterranean and the Middle East is 75% larger in the scenario of increasing optimal production (scenario 1) compared with the scenario of constant optimal production (scenario 2). Results of such a statistical model can be considered valid up to the period in the future when the 2°C global warming level is reached (Zampieri et al. 2019a), which is projected to happen in the 2040s according to the RCP8.5 scenario. After that period, the CSI results are to be considered qualitative (gray lines in Fig. 6). Fig. 6 Future wheat production estimated considering two scenarios of optimal production: scenario 1 with constantly increasing optimal production (panels a, c, and e) and scenario 2 with constant optimal production (panels b, d, and f).   East (1998East ( -2017 according to the FAOSTAT data set (FAOSTAT, 2019). Positive and negative trends are statistically significant at the 90% level (p < 0.1)

Discussion
Wheat production in the Mediterranean and the Middle East is currently threatened by increasing heat waves and droughts (Fontana et al. 2015;Zampieri et al. 2017Zampieri et al. , 2019aTebaldi and Lobell 2018a). Given the socioeconomic importance of wheat in the Mediterranean (Asseng et al. 2018), as well as the global consequences of regional yield anomalies and losses (Chatzopoulos et al. 2019), we have here addressed the need of reliable indicators to estimate agricultural resilience with respect to climate extremes.
This study presents an application of a new indicator to measure crop production resilience (Zampieri et al. 2019d). The crop production resilience indicator is defined as the reciprocal of the squared coefficient of variance, which is a special case of the generalized entropy index (Shorrocks 1980). The crop production resilience indicator applied to crop production time series is directly/inversely proportional to the return period/frequency of extreme events causing severe production losses (Zampieri et al. 2019d) and has been used to estimate also general vegetation resilience and the reliability of the related ecosystem services (Zampieri et al. 2019c). Compared with other indicators of stability, a concept close to resilience, such as the reciprocal of the coefficient of variance (Kahiluoto et al. 2019a;Renard and Tilman 2019), the crop production resilience indicator is consistent with the original ecological definition of resilience (Holling 1973(Holling , 1996. In this paper, we have shown how to remove the trend from the time series in order to compute the resilience indicator from nonstationary crop production data. This method has been illustrated assuming two contrasting agronomic scenarios. The resulting time series, where also the climate limitation factor has been added, are representative of the range of projections produced by crop models in the Mediterranean area (Dettori et al. 2017;Dixit et al. 2018;Tebaldi and Lobell 2018a, b;Webber et al. 2018). It is worth noting that our methods do not consider change of seasonality, which is consistent with the implementation of varieties characterized by slower phenology, adapted to warmer climates (Dixit et al. 2018;Rezaei et al. 2018). The estimated decline of wheat resilience due to climate change is even more alerting when considering also the observed effects of varieties' diversity loss estimated for the current climate (Kahiluoto et al. 2019a). In this respect, while the limited data available in the observations can limit the reliability of the resilience estimation (Kahiluoto et al. 2019b;Snowdon et al. 2019;Zampieri et al. 2019d), we highlight the advantage of multi-model-based results allowing to reach robust conclusions.

Conclusions
In this study, a framework for agricultural resilience assessment has been introduced and applied to analyzing wheat production projections for the top ten producing countries in the Mediterranean and the Middle East (of primary importance for the durum wheat cultivation). Further evaluation of the annual crop production resilience, here preliminarily estimated under strong assumptions, is needed especially based on a wider ensemble of crop production simulations such as the ones done within the AgMIP-coordinated exercise (Rosenzweig et al. 2013).
Nevertheless, this simplified analysis already provides serious reasons to act. If no adaptation to climate change will take place, the Mediterranean wheat production is expected to experience severe losses associated with higher interannual variability already at 1.5°C global warming. These findings support the need of intensifying climate change mitigation efforts consistently with the message of the Paris Climate Agreement and act immediately to develop effective adaptation strategies.
Adaptation measures, related to technological and agromanagement improvements, might balance out the effects of the increasing variability when increasing constantly the production levels. However, wheat production is already stagnating in some countries due to maximum potential yield having already been reached (Brisson et al. 2010) and recent social crises (Kelley et al. 2015).
With less ambitious global climate mitigation targets, and with no effective adaptation strategies (needed to be consistently sustainable considering the decreasing water resources in the Mediterranean region), we anticipate that societies will have to face significantly less reliable wheat production.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.