Crop diversity effects on temporal agricultural production stability across European regions

Stabilizing agricultural production is fundamental to food security. At the national level, increasing the effective diversity of cultivated crops has been found to increase temporal production stability, i.e., the year-to-year stability of total caloric production of all crops combined. Here, we specifically investigated these effects at the regional level for the European Union and tested the effect of crop diversity in relation to agricultural inputs, soil properties, climate instability, and time on caloric, protein, and fat stability, as we hypothesized that the effect of diversity is context dependent. We further investigated these relationships for specific countries. We found that greater crop diversity was consistently associated with an increase in production stability, particularly in regions with large areas equipped for irrigation and low soil type diversity. For instance, in Spain and Italy, crop diversity showed the strongest positive effect among all predictors, while on the European level, the stabilizing effect of nitrogen use was substantially higher. In Germany, the crop diversity-stability relationship was weak, suggesting that crops react similarly to climatic, economic, and political factors or are grown in the same periods. With this study, we substantiate previous findings that crop diversity stabilizes agricultural caloric production and extend these with regard to protein and fat. The results elucidate the key drivers that enhance production stability for different European countries and regions, which is of key importance for a comparably productive agricultural region like Europe.


Introduction
Increasing the temporal stability of agricultural production is key to foster the resilience and reliability of agricultural systems and is thus fundamental to food security (Knapp and van der Heijden 2018), in particular in the light of climate change and rising demands for food (Challinor et al. 2014;Valin et al. 2014). While for a long time agricultural assessments have focused on agricultural production and productivity, understanding the temporal stability of agricultural production and its underlying drivers has only recently received more attention (Gaudin et al. 2015;Ray et al. 2015;KC et al. 2016; Raseduzzaman and Jensen 2017; Knapp and van der Heijden 2018; Mehrabi and Ramankutty 2019; Renard and Tilman 2019).
Resilience theory suggests that diversity reduces risks and increases resilience in agricultural systems or the ability to maintain functioning in the face of disturbance and change (Biggs et al. 2012;Mijatović et al. 2013;Kremen and Merenlender 2018). For example, if one crop fails due to climate instability and extreme events, other crops can compensate the losses (Yachi and Loreau 1999;Lin 2011;Raseduzzaman and Jensen 2017). Diversity encompasses variety (how many different elements), balance (how many of each element), and disparity (how different the elements are from one another) (Stirling 2007;Biggs et al. 2012). At the national level, it has been found that a greater effective diversity of cultivated crops based on the Shannon index increases temporal production stability, i.e., the year-toyear stability of the total caloric production of all crops combined (Renard and Tilman 2019;Egli et al. 2020).
Agricultural inputs such as fertilizer and irrigation further stabilize agricultural production (Rist et al. 2014;Ray et al. 2015;Knapp and van der Heijden 2018;Renard and Tilman 2019), potentially instead of or in addition to crop diversity. During the last decades, dependency on fertilizer and irrigation and the spatial coverage of their application have heavily increased, and hence yielded large production gains (Tilman et al. 2001;Foley et al. 2011). Favorable conditions to grow crops related to soils, topography, and climate might be further important to foster stable agricultural production (Pan et al. 2009). Increased soil type diversity is likely to increase temporal production stability, in particular if crop diversity is low because it fosters the asynchrony between different plots of the same crop (Coomes et al. 2019;Mehrabi and Ramankutty 2019). Finally, climate instability reduces temporal yield and production stability at the national level (Ray et al. 2015;Renard and Tilman 2019;Egli et al. 2020).
To better understand the effect of crop diversity and other factors on the temporal production stability of agricultural systems, these processes need to be investigated at different scales. To the best of our knowledge, regional assessments at a large spatial extent such as Europe are missing. Moreover, better understanding of these processes within political or economic unions or individual countries is particularly important because opportunities to implement policies to foster resilient food system are potentially higher. For example, within the European Union, the Common Agricultural Policy provides a framework to direct agricultural management, also with regard to sustainability targets (Louhichi et al. 2017;Früh-Müller et al. 2019). Subnational assessments can inform such policies in a more nuanced way. For example, if crop diversity would be found to stabilize agricultural production more effectively than agricultural inputs, this would provide a possibility for more sustainable management.
Here, we compiled long-term data at the regional level data in Europe and used linear and linear mixedeffects models to (i) study the effects of crop diversity on the year-to-year stability of total caloric, protein, and fat production (hereafter 'production stability,' 'caloric stability,' 'protein stability,' or 'fat stability') in relation to agricultural inputs, agricultural suitability, soil productivity, soil type diversity, climate instability, and time, and (ii) to investigate these relationships in specific countries. We hypothesize that besides crop diversity, soil productivity, soil type diversity, and agricultural inputs are positively associated with production stability, while climate instability and time are negatively related (Table 1). We further assume that the positive relationship of crop diversity with production stability is stronger in regions with lower soil productivity and soil type diversity, and with higher climate instability.

Materials and methods
We collected agricultural, soil, and climate data across the NUTS 2 regions (Nomenclature of Territorial Units for Statistics) as defined in the year 2013 for the European Union (including the UK) for up to four 10-year time intervals (1978-1987, 1988-1997, 1998-2007, 2008-2017). NUTS 2 regions typically reflect the second subdivision of a nation into subnational units. For rather small countries including Cyprus, Estonia, Latvia, Lithuania, Luxemburg, and Malta, the NUTS 2 level reflected the entire nation. For Germany, sufficient data was only available at the state level (NUTS 1). To exclude regions where crop cultivation is of very minor relevance, we sorted them by the average area of arable land between 1978 and 2017 in descending order and included regions up to a cumulated area of 99% of the total area of all regions. We further excluded regions in Ireland and the Netherlands because they use much of their fertilizer on pastures (Renard and Tilman 2019). We only included arable crops, for which nutrient data could be clearly assigned (FAO 2001). To exclude crops of minor relevance, we sorted the remaining crops by the average harvested areas between 1978 and 2017 in descending order and only included crops up to a cumulated area of 99% of the total harvested area of all crops. Furthermore, we removed data points where harvested Time interval Four ten-year time intervals (1978-1987, 1988-1997, 1998-2007, 2008-2017) − Temporal production stability decreases over time due to climate change and other global change processes Interaction effects with crop diversity Agricultural inputs (nitrogen use, irrigation) + / − Crop diversity enhances the effect of agricultural inputs, or the effect of crop diversity is highest in regions with low agricultural inputs Soil productivity − The effect of crop diversity is highest in regions with marginal agricultural lands Soil type diversity − The effect of crop diversity is highest in regions with homogenous agricultural landscapes Coomes et al. 2019 Climate instability (temperature, precipitation) + The effect of crop diversity is highest in regions with high climate instability Yachi and Loreau 1999;Biggs et al. 2012 area was above zero, but production was not reported or zero, or vice versa. For each time interval and region, we only included crops where harvested areas and production were reported for at last 8 years. In a last step, we removed areas where only one crop remained. Our final datasets consisted of 409 data points representing 161 regions and 18 crops (Table S1). Each data point represents one region in one time interval. Data for multiple time intervals was available for 81.99% of these regions. On average, the crops considered represent 44.78% of the arable land of these regions and 41.94% of the arable land reported in the Eurostat database. Details on the underlying data sources can be found in Table 2.

Response variables
To calculate production stability, we used agricultural production data of arable cops from the EUROSTAT database (Table 2). We converted crop-specific production from tons to calories, protein, and fat using standardized nutritive factors (Table 2). We note that other aspects could be considered to aggregate crops ). To account for crop stability independent of long-term trends, we time-detrended annual production data by regressing annual total caloric production of all or individual crops on year squared for each time interval and region (Renard and Tilman 2019). We calculated caloric, protein, or fat stability as where R is the mean of the non-time-detrended total regional caloric, protein, or fat production for each time interval and R is the standard deviation of the respective time-detrended production (Mehrabi and Ramankutty 2019; Renard and Tilman 2019).

Explanatory variables
We used crop diversity, fertilizer, irrigation, soil productivity, soil type diversity, temperature and precipitation instability, and time as explanatory variables of production stability. Regarding crop diversity, we calculated effective diversity as the exponential of the Shannon diversity of the harvested areas of the different crops for each year and region (as reported in the EUROSTAT database) and calculated mean values for each time interval (Hill 1973).
To account for land use intensity, we included fertilizer use and irrigation. We used national data from the FAOSTAT database, as we wanted to cover a large temporal extent regarding agricultural production data (40 years), for which regional data is largely incomplete (Table 2). Therefore, we calculated mean nitrogen use for each time (1) P R = R ∕ R interval and nation relative to the respective mean cropland area. For irrigation, we used the mean area equipped for irrigation relative to agricultural land for each time interval and nation.
To describe soil productivity, we used a soil productivity index within croplands based on soil, topographic, and climate criteria and in relation to potential productivity increase due to fertilization (Table 2). For each region we extracted the mean productivity index. We used the digital soil map of the world to calculate the soil type diversity as the exponential of the Shannon diversity of the different soil types within the agricultural areas of each region based on the CORINE land cover dataset of the year 2000 (Table 2).
To describe climate instability, we used gridded monthly temperature (°C) and precipitation (cm) data on 0.5-degree resolution ( Table 2). For each region, we extracted mean temperature and precipitation during the growing season of major crops and within the agricultural areas as defined above (Table 2). We then calculated temperature and precipitation instability as the negative of the mean temperature and precipitation, respectively, over its standard deviation for each time interval and region (Renard and Tilman 2019). Finally, we included the four time intervals for which the other temporal variables were aggregated as an explanatory variable.

Statistical analyses
We used linear models to test the dependence of production stability on crop diversity, as well as on fertilizer, irrigation, soil productivity, soil type diversity, temperature and precipitation instability, and time. Moreover, we included pairwise interactions of crop diversity with agricultural management, soil, and climate variables. Production stability was log-normally distributed (AIC = 2429.65, 2407.69, and 2183.87 compared to AIC = 2663.72, 2689.97, and 2445 for a normal distribution for caloric, protein, and fat stability, respectively). To equalize spread and reduce leverage, we square root-transformed fertilizer use (explanatory variable). We scaled explanatory variables (mean-centered and divided by their standard deviations) to reduce remaining collinearity between main effects (Gelman and Hill 2006) and to make regression coefficients comparable. We tested both a linear model and a linear mixed-effects model fitting a random intercept for each nation. We based our results on the linear mixed-effects model (AIC = 609.15, 618.07, and 622.55 for caloric, protein, and fat stability, respectively) as it exceeded the performance of the linear model (AIC = 677.00, 680.77, and 624.63). The maximum variance inflation factor was 2.73, indicating that multicollinearity was successfully removed (Marquardt 1970).
To assess whether these relationships are similar for individual countries, we used separate linear regression models for countries with more than 40 data points (Germany, Spain, France, and Italy). As fertilizer and irrigation data was only available on the national level, we removed these variables. Moreover, we did not include the interaction terms to reduce the number of regression parameters in relation to the sample size (Heinze et al. 2018).
We refitted all selected models with 1000 bootstrap replications to assess the robustness of our conclusions (Heinze et al. 2018). We used the statistical software package R 3.5.1 (R Core Team 2019) run via RStudio (RStudio Team 2015) for data analysis. To assess distributional assumptions in the response variables, we used the 'fitdistrplus' package in R (Delignette- Muller and Dutang 2015). For linear mixedeffect models, we used the 'nlme' library in R based on restricted maximized log-likelihood ('REML') (Pinheiro et al. 2019). All data that support the findings of this study, as well as related codes for data preparation and analyses, are openly available here: https:// doi. org/ 10. 5281/ zenodo. 54580 18.

Results
Production stability was highest in regions of Denmark, Germany, Greece, Italy, Poland, Portugal, southern Sweden, and the UK (Fig. 1a-c), as well as southern Finland for fat stability. Crop diversity was particularly high in Portugal and the Baltic states (Fig. 1d). Nitrogen use was highest in Belgium, Germany, and the UK, and irrigation was highest in Greece, Italy, and Romania ( Fig. 1e-f). Soil productivity was generally higher in Central and Northern Europe (Fig. 1g). Regions with high soil type diversity were found in all countries except in northern Europe and the Baltic states (Fig. 1h). Temperature instability was highest in Eastern and Northern Europe, while precipitation instability was particularly high on the Iberian Peninsula and in Southeastern Europe (Fig. 1i-j).
Regarding the main effects, nitrogen use was positively and significantly associated with caloric stability (p < 0.05), while temperature and precipitation instability showed a significant negative relationship (p < 0.05; Fig. 2; Table S2). Crop diversity (p = 0.14), irrigation (p = 0.62), soil productivity (p = 0.09), and time (p = 0.75) showed a positive, and soil type diversity (p = 0.12) a negative relationship, but all of them were insignificant. However, in regions with large areas equipped for irrigation, crop diversity was positively and significantly associated with caloric stability (p < 0.05; Fig. 3; Table S2). Furthermore, the positive effect of crop diversity was larger in regions with low nitrogen use (p = 0.29), soil productivity (p = 0.52), soil type diversity (p = 0.06), and precipitation instability (p = 0.43), and in regions with high temperature instability (p = 0.20), yet these interactions were insignificant. The variance in caloric stability explained by the fixed effects (R 2 marginal) was 0.30 and the total variance explained (R 2 conditional) was 0.58. Except for time, the observed direction of main effect and interactions were identical in at least 78.8% of the bootstrap replicates (Table S2). Crop diversity, for example, was positively associated with caloric stability in 91.7% of the replications, and its positive interaction with irrigation in 99.4% of the replications.
These relationships were nearly identical for protein stability (Fig. 2, Table S2). Regarding fat stability, the overall patterns were also very similar, but crop diversity and time showed a significant and stronger positive relationship (p < 0.05). Moreover, the interaction with irrigation was not significant (p = 0.16) and the effect of crop diversity was larger in areas with high precipitation instability (p = 0.25), but not in areas with high temperature instability (p = 0.74). The variance in protein and fat stability explained by the fixed effects (R 2 marginal) were 0.30 and 0.22, and the total variances explained (R 2 conditional) were 0.57 and 0.29, respectively.
In Germany, none of the variables was significant (p > 0.05) ( Fig. 4; Table S3) and the variance explained was relatively low (R 2 = 0.22). In Spain (R 2 = 0.63), France (R 2 = 0.55), and Italy (R 2 = 0.28), crop diversity showed the strongest positive association with caloric stability among all variables, though the effect was insignificant in France (p = 0.13). Soil productivity was negatively associated with caloric stability in France (p < 0.05), while soil type diversity showed a positive but insignificant effect in Italy (p = 0.44). Climate stability generally showed a negative relationship across all four countries, except regarding temperature instability in Spain. In Spain and France, the effect of time was significant and negative (p < 0.05), indicating that caloric stability decreased during the last decades.  Determinants of caloric (green), protein (blue), and fat (orange) stability. Regression coefficients (± SE) for all main effects included in the linear mixed-effects regression models (n = 409). Caloric, protein, and fat stability (mean divided by standard deviation) were logtransformed, and nitrogen use was square root-transformed. Each predictor variable was standardized to zero mean and one standard deviation across all regions and time intervals. *p < 0.05;**p < 0.01; ***p < 0.001; NS = not significant Fig. 3 Effect of crop diversity in combination with the other explanatory variables on caloric stability (n = 409). Interaction of crop diversity with a nitrogen use, b area equipped for irrigation, c soil productivity, d soil type diversity, e temperature instability, and f precipitation instability and associated p values. Predicted values for caloric stability were back-transformed from log-transformation. Predictions were calculated using the observed range of crop diversity, while keeping all the other variables at their mean values or the second variable of the interaction at mean − 1 SD (low) or mean + 1 SD (high)

Discussion
Our results show that crop diversity is in general positively associated with production stability at the regional level in the European Union, particularly in regions with large areas equipped for irrigation and with low soil type diversity. We found evidence that the positive effect of crop diversity is higher in regions with instable temperatures, supporting resilience theory according to which diversity provides an insurance to maintain functioning in the face of disturbance and change (Yachi and Loreau 1999;Biggs et al. 2012). While our results generally confirm the stabilizing effect of crop diversity, its importance is generally less pronounced compared to the national level globally (Renard and Tilman 2019;Egli et al. 2020). Regarding fat stability, however, the effect of crop diversity was stronger. This highlights the importance of evaluating nutritional aspects in the food security debate that go beyond the assessment of calories produced (Hertel and Baldos 2016;Berners-Lee et al. 2018).
Across all countries and variables considered, nitrogen use had the largest positive effect on production stability. Accordingly, agricultural inputs are not only relevant to close yields gaps (Mueller et al. 2012), but also to increase temporal production stability (Knapp and van der Heijden 2018; Renard and Tilman 2019). The effect of nitrogen fertilizer in the European Union was more important than globally (Renard and Tilman 2019), which is not surprising given that a large share of European croplands is characterized by medium to high intensity (Václavík et al. 2013;Levers et al. 2018). However, this dependency on agricultural inputs has serious implications for biodiversity and ultimately human well-being (Tscharntke et al. 2012;Loos et al. 2014;West et al. 2014;Beckmann et al. 2019), and associated specialization may also increase climatic sensitivity of agriculture (Ortiz-Bobea et al. 2018). Furthermore, nitrogen fertilizer already reached a peak year in 1983 (Seppelt et al. 2014). Our results indicate that Fig. 4 Determinants of caloric stability in Germany (a, n = 43), Spain (b, n = 56), France (c, n = 64), and Italy (d, n = 50). Regression coefficients (± SE) for all effects included in the linear regression models. Caloric stability (mean divided by standard deviation) was log-transformed. Each predictor variable was standardized to zero mean and one standard deviation across all regions and time intervals. *p < 0.05;**p < 0.01; ***p < 0.001; NS = not significant crop diversity could partially compensate the stabilizing effect of nitrogen use (negative yet insignificant interaction). Accordingly, more integrated approaches to tackle the multiple challenges agricultural landscapes face today such as crop diversification, renewable inputs, and other agroecological principles are needed (Gurr et al. 2016;Landis 2017;Kremen and Merenlender 2018;Lee et al. 2019;Seppelt et al. 2020). For example, more diverse cropping systems, beyond just organic agriculture, could potentially achieve higher production with lower environmental impacts than today (Seufert et al. 2012;Iverson et al. 2014;Tscharntke et al. 2015;Muller et al. 2017;Knapp and van der Heijden 2018;Sirami et al. 2019).
Regarding soil characteristics, our results were mixed. As expected, soil productivity showed a positive association with production stability that was comparable to the effect of crop diversity (Pan et al. 2009). Yet insignificant, but largely stable, the negative interaction between crop diversity and soil productivity indicates that crop diversity is a promising mechanism to stabilize agricultural production in marginal areas (Shahid and Al-Shankiti 2013). While the main effect of soil type diversity was insignificant and negative, we found evidence that in regions with low soil type diversity, crop diversity is needed to increase field diversity and thus the stability of agricultural production (Coomes et al. 2019). While this also indicates that crop diversification could decrease production stability in areas with high soil type diversity, such trade-offs might be avoided with diversification at the farm level .
We found that climate instability substantially decreased temporal production stability in Europe. Thus, expected increases in climate instability in the light of climate change may not only lead to production losses (Challinor et al. 2014;Liang et al. 2017;Zhao et al. 2017), but also to lower temporal production stability. However, climate seems to play a minor role compared to the effect size of nitrogen use. This indicates that at the smaller scale, management factors are more relevant. Accordingly, until recent, climate instability in Europe can be buffered by agricultural management, which was also observed for agricultural productivity in the USA (Liang et al. 2017). However, this might change due to increasing risks of climate change-induced soil erosion, droughts, and other aspects that affect agricultural productivity (Mezösi et al. 2013;Deutsch et al. 2018).
The individual country models revealed important differences regarding the determinants of production stability. For example, in Spain and Italy, crop diversity showed the strongest positive effect among all predictors. In contrast, crop diversity in Germany had no clear positive effect, indicating the cultivation of crops that react similarly to climatic, economic, and political factors or are grown in the same periods (Egli et al. 2020, which poses risks to production stability, in particular in the light of climate change effects such as the extreme droughts in 2018-2020, as well as rising perturbations in global markets (Suweis et al. 2015;Ortiz-Bobea et al. 2018). In contrast to the patterns observed at the European level, caloric stability decreased over time in Spain and France. Given their importance in European agriculture, this might have implications for food security. These country-specific findings should be considered in strategies to stabilize agricultural production.
The study presented here faces two major limitations. First, data availability and quality of the underlying datasets are limited. After processing, several regions and crops were excluded due to data gaps and our dataset represented only around 45% of the regional arable land. Nevertheless, this should not systematically bias our results as production stability and crop diversity refer to the same data pool. This is also indicated by the robustness of our main conclusions. Furthermore, long-term information on nitrogen use and irrigation was only available at the national level, for all crops combined and with regard to the total agricultural area including pastures for irrigation. The effects of management on production stability might also differ across different regions and climates. The aggregation of climate and soil data to the regional level may introduce additional biases (Hoffmann et al. 2016). Second, the explanatory power of the fixed effects (R 2 marginal = 0.22-0.30) indicates that other important determinants of production stability are missing, for example, related to farmers' behavior and knowledge, as well as other aspects determining total factor productivity, including pesticide use, capital stock in agriculture (sum of investment and physical assets, including machinery), storage facilities, and modern technologies, but also ecosystem services that are important for maintaining productivity on the long term (Coomes et al. 2019). However, temporal and spatial data coverage of these potential determinants is yet limited.
This study opens new avenues for future research. Future assessments could include more crop types, as well as different crop varieties and functional traits, as they most likely play a fundamental role for resilience (Ficiciyan et al. 2018). Further studies could also target other aspects of resilience, disturbances, and agricultural diversity (Letourneau et al. 2011;Lin 2011;Cottrell et al. 2019;Egli et al. 2019). Moreover, the combined effects of crop diversity and other factors on agricultural production and temporal production stability and in particular the identification of synergies, e.g., through intercropping, and strategies to reduce fertilizer dependency merit further attention (Raseduzzaman and Jensen 2017; Martin-Guay et al. 2018).

Conclusion
Our study suggests that besides nitrogen use, crop diversity increases temporal production stability at the regional level in Europe, particularly in regions with large areas equipped for irrigation and with low soil type diversity. Our results also show that climate instability reduces temporal production stability and that the stabilizing effect of crop diversity seems to be higher in regions with high temperature instability. These findings emphasize the need for integrated land-use planning including diversification and a sustainable use of agricultural inputs.
Funding Open Access funding enabled and organized by Projekt DEAL. L.E. acknowledges funding from the Helmholtz Association (Research School ESCALATE, VH-KO-613). We thank Volker Grimm and Zia Mehrabi for constructive discussions and Mick Wu for statistical support.
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:// creat iveco mmons. org/ licen ses/ by/4. 0/.