Assessing the capacity of agricultural research and development to increase the stability of global crop yields under climate change

Assessing the capacity of agricultural research and development to increase the stability of global crop yields under climate change | Research Square window.SnipcartSettings = { analytics: { enabled: false } }; (function() { var accessVector = localStorage.getItem('access_vector') || ''; window.dataLayer = window.dataLayer || []; if (accessVector) { window.dataLayer.push({ user: { profile: { profileInfo: { snid: accessVector } } } }); } })(); (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0],j=d.createElement(s),dl=l!='dataLayer'?'&l='+l:'';j.async=true;j.src='https://www.googletagmanager.com/gtm.js?id='+i+dl;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-K279D39R'); Browse Preprints In Review Journals COVID-19 Preprints AJE Video Bytes Research Tools Research Promotion AJE Professional Editing AJE Rubriq About Preprint Platform In Review Editorial Policies Our Team Advisory Board Help Center Sign In Submit a Preprint Cite Share Download PDF Research Article Assessing the capacity of agricultural research and development to increase the stability of global crop yields under climate change Toshichika Iizumi, Toru Sakai, Yoshimitsu Masaki, Kei Oyoshi, and 4 more This is a preprint; it has not been peer reviewed by a journal. https://doi.org/ 10.21203/rs.3.rs-3895355/v2 This work is licensed under a CC BY 4.0 License Status: Published Journal Publication published 21 Mar, 2025 Read the published version in PNAS Nexus → Version 2 posted You are reading this latest preprint version Show more versions Abstract Agricultural research and development (R&D) has increased crop yields, but little is known about its ability to increase yield stability in the context of increasingly frequent extreme weather events. Using a grid yield dataset, we show that from 2000 to 2019, the standard deviation (SD) of yield anomalies for maize, rice, wheat and soybean, increased in 20% of the global harvested area. Based on random forest models relating yield anomaly to climate, soil, management and public R&D expenditure, we show that cumulative agricultural R&D expenditure, proportion of growing season exposed to optimal hourly temperatures, and dry and very wet days are key factors explaining crop yield variability. An attribution analysis based on large ensemble climate simulations with and without human influence on the global climate shows that unfavorable agro-climatic conditions due to climate change has increased SD, while higher R&D expenditure has led to more contrasting trends in SD over 2000–2019. Although R&D has continued steadily in most countries, this study indicates that the progress made in R&D since 2000 may have lagged behind the unfavorable effect of climate change on yield variability. Agroecology Climate Analysis and Modeling Figures Figure 1 Figure 2 Figure 3 Figure 4 Figure 5 Figure 6 Figure 7 Figure 8 Significance Statement Improved technologies and management practices - fostered by agricultural research and development (R&D) - have helped to increase yield stability, but it is not known whether R&D investments have been sufficient to fully offset the increase in yield variability induced by climate change. Based on a global dataset, our analysis quantifies the relationships between R&D expenditure, climate change, and yield stability for major agricultural commodities worldwide. Our analysis reveals that R&D investments have probably not been high enough since 2000 to fully mitigate the increased variability of yields resulting from climate change. Future research is needed to establish a causal relationship between R&D and yield stability underlying the empirical results found in this study. Introduction Among the four pillars of food security (availability, access, utilization and stability) ( 1 ), ensuring the production stability is crucial for preventing any deleterious consequences cascading to agri-food systems downstream. The disruptions to food supply chains caused by recent conflicts and pandemic reaffirm the importance of stability in the interconnected world ( 2 – 4 ). In the past, production stability has been improved by optimizing growing areas, harvestable area fraction, annual number of harvests and cultivar tolerance to abiotic and biotic stresses ( 5 , 6 ). Increasing yield stability remains an important avenue to make crop production more climate-resilient, especially because yield growth has been the predominant contributor to production increase in many parts of the world in recent decades, while the increase in cultivated areas have been a much smaller contributor since 2000 ( 5 , 7 ). From 1981 to 2010, yield stability for major crops (maize, rice, wheat and soybean) decreased in 9–22% of the global harvested area in response to increased crop exposure to heat, yet 19–33% of the area saw stabilizing yields ( 8 ). After the study by Iizumi and Ramankutty ( 8 ) published in 2016, many production shocks associated with extreme climate events have been reported ( 9 – 12 ). Research and development (R&D) on climate-resilient agri-food systems is expanding, including breeding targeting emerging agro-climatic conditions ( 13 ) and climate-smart agriculture. However, although R&D innovations have been extensively tested at experimental stations, there is limited evidence showing that advances in R&D and producers’ adoption of improved seeds, technologies and management practices have reduced climate risks at large spatial scales ( 14 ). Despite improvements in breeding and weather forecasting in the last decades, some studies suggest that the yield sensitivity to extreme climate events has not been substantially reduced in high-yielding regions ( 15 , 16 ). It is therefore crucial to update our understanding on yield stability and its drivers. In this paper, we contribute to this effort by studying how climate change and agricultural R&D have affected yield stability since 2000. Recent advances in climate attribution science and resulting climate model simulations ( 17 – 19 ) allow researchers to estimate the contributions of human-induced climate change to average yields ( 7 , 20 ) and crop failures ( 21 , 22 ). However, contrasting to the accumulated evidence showing that climate change has slowed down the growth of average yields ( 7 ), little is known about climate-induced changes in yield stability. In addition, it remains difficult to assess the capacity of R&D on reducing climate risks, although this information is often requested by research project funders ( 23 ). To fill these gaps, we present here a global analysis of the changes in yield stability due to climate change and R&D from 2000 to 2019, using a large ensemble climate model simulation ( 17 ) and machine learning models that relate interannual variation in yields of maize, rice, wheat and soybean at the 0.5° grid scale (~ 55 km at the equator) to climate, soil, management and agricultural R&D expenditure. The four crops considered here are of high importance as they provide nearly two-thirds of the world agricultural calories ( 24 ). Results R&D capacity to increase and stabilize crop yields Yield anomaly—the departure of annual yield from trend of average yield—is a key metric measuring year-to-year yield fluctuation. The analysis of country-level yield anomalies over 2000–2019 showed that increases in public R&D expenditure are proportional to increases in average yields. Considering the four crops taken together, countries with large areas under cultivation (≥ 1 Mha) had invested considerably more in R&D than countries with smaller cropping areas (Fig. 1 A). Average yield and yield anomaly SD are positively correlated (Fig. 1 B) but, for a given average yield level, countries with higher levels of R&D expenditure often achieved higher stability (lower SD and/or lower CV = SD/average) than countries with lower levels of R&D expenditure (Fig. 1 C). [Fig. 1 here] Recent patterns of agricultural R&D expenditure The level of cumulative agricultural R&D expenditure varied considerably between countries, although increasing trends are quite common (Fig. 2 A). Compared to annual R&D expenditure, their cumulative values, which take into account the time lag between research and adoption of new technologies by producers and technology obsolescence, are relevant to explain the increase in average yield and yield stability. Since 2000, annual R&D expenditure was the highest in China, followed by the United States of America (the US) (Fig. 2 B). China’s annual spending increased rapidly from 3.6 billion USD ( $ B) in 2000 to 22.9 $ B in 2019, while the US annual spending increased from 4.2 $ B to 5.2 $ B for the period. The level of annual R&D expenditure in other countries were much less than these two countries. Although cumulative R&D expenditure increased in most countries, for some countries, cumulative R&D expenditure stagnated (e.g., in France) (Fig. 2 A, C). When the level of annual R&D expenditure was not sufficient to maintain the level of cumulative R&D expenditure, which is subject to annual rates of technology obsolescence, cumulative R&D expenditure declined with time (e.g., in Russia and Sudan). [Fig. 2 here] Recent changes in yield stability To analyze the evolution of yield stability, SD was selected as the main metric because SD was more sensitive than CV in detecting climate-induced changes in yield anomalies. Under climate change, yield anomaly SD for a certain period may increase as a result of both occurrences of poor harvests (due to severe agro-climatic conditions) and good harvests (due to favorable agro-climatic conditions, or to elevated carbon dioxide concentrations, which enhance photosynthesis). However, with CV, such a change in SD can be masked by an increase in average yield because CV decreases when average yield increases faster than yield anomaly SD over time. Our grid-wise analysis showed that 20% of the global harvested area (129.2 Mha, out of 646.2 Mha in 2010) had increased yield anomaly SD in 2000–2019 (Fig. 3 B), while areas with increased CV were only 58% of the areas with increased SD (75.8 Mha; Fig. 3 A). The overall increasing trends in average yields over the same period explained this difference in results between SD and CV (Fig. 3 C). [Fig. 3 here] Modelling the effects of climate and R&D on yield anomaly Machine learning random forest (RF) models were trained with current climate conditions (represented by the retrospective meteorological forcing data; Materials and methods ) to attribute the changes in modeled yield anomaly SD to the R&D and climate-related explanatory variables (Table 1). The models explained 19% (rice) to 43% (soybean) of the variance of the grid yield anomaly between the locations, seasons and years when evaluated using out-of-bag data. Intermediate percentage of explained variance was obtained for maize and winter wheat (26% for both) and spring wheat (39%). Concerning the model ability to reproduce interannual yield variations, the leave-one-year-out cross-validation led to average correlations ranging from 0.639 (major-season rice) to 0.793 (soybean), and to average root mean squared errors ranging from 0.66 (spring wheat) to 0.81 (major-season rice) in SD units (figs. S4–S5). Importantly, the RF model projections obtained with d4PDF factual climate simulations (representation of the actual climate simulated by the climate model, including human-induced climate change) were able to reproduce the geographic pattern of the changes in yield anomaly SD for the crops, with moderate agreement with the actual pattern derived from the grid yield dataset [Cohen’s kappa (κ) of 0.55; Fig. 4 A, C]. The agreement with the actual pattern was as high as κ = 0.71 when the RF models were forced by the retrospective meteorological forcing data (Fig. 4 A, B). We also checked that d4PDF factual and counterfactual climate simulations fell within the ranges covered by the training dataset as the RF model performance often diminish outside the training data (only 0.18% of cases across all climate variables, locations, seasons, years and ensemble members were outside the ranges covered by the training dataset). [Table 1 here] [Fig. 4 here] Variable importance The RF models were used to rank the relative importance of R&D, climate, soil and management variables in explaining the yield anomaly. Two methods were considered to determine the importance of variables ( Materials and methods ). Cumulative R&D expenditure was identified as the primary factor, irrespective of the crops and methods (Fig. 5 ). N application rate was often detected as the next important factor after R&D expenditure for the crops except soybean (fig. S6). Temperatures ranging from 15°C to 35°C, which encompass the exposures to temperature fluctuations around optimal thermal conditions for crops, were identified as important, albeit with some crop-specific variation. For maize, rice and soybean, the relatively higher temperatures from the 5th bin (15–20°C) to the 9th bin (≥ 35°C) were more important. For winter and spring wheat, the relatively lower temperatures from the 4th bin (10–15°C) to the 5th bin (15–20°C) showed stronger effects. Among the precipitation variables, the 1st bin (< 1 mm d –1 ) and the 5th bin (≥ 30 mm d –1 ) tended to the most important, with variations between the crops and methods. On average, SOC was one of the most important soil variables. Irrigation intensity was less influential. [Fig. 5 here] Climate effect Effects of climate and R&D on the changes in yield stability at SD level were estimated by comparing the yield anomalies simulated by the RF models for different scenarios. Two yield anomaly time series were generated using d4PDF factual and counterfactual climate simulations, representing the historical climate (including human-induced climate change) and the preindustrial climate (without human-induced climate change), respectively ( fc and ct runs; Fig. 6 ). R&D expenditures were fixed to those observed in both scenarios. We classified the change in agro-climatic condition due to climate change as “unfavorable” (“favorable”) when the simulated yield anomaly SD was higher (lower) with the factual climate than with the counterfactual climate. The global grid-wise attribution revealed that unfavorable agro-climatic conditions were relatively more prevalent in the cropland area where an increase in yield anomaly SD occurred than in the area with decreasing yield anomaly SD. The area with unfavorable agro-climatic conditions was estimated to account for the half (0.5) of the area with an increase in yield anomaly SD; this proportion (0.5) was derived by comparing the area with both unfavorable agro-climatic conditions and increasing SD (10% of the global harvested area; Fig. 7 A, D) and the area with increasing yield anomaly SD all agro-climatic conditions taken together (20% of the global harvested area; Fig. 7 C, D). This proportion (0.5) is twice as high as the proportion (0.25) of area with unfavorable agro-climatic conditions and decreasing SD (4%) in the area with decreasing yield anomaly SD all agro-climatic conditions taken together (16%) (Fig. 7 D). In other words, the frequency of unfavorable agro-climatic conditions is twice as high in areas where the yield anomaly SD is increasing than in areas where the yield anomaly SD is decreasing. Consistent with the above result, favorable agro-climatic conditions were relatively more common in the area where a decrease in yield anomaly SD occurred than in the area with increasing yield anomaly SD. Favorable agro-climatic conditions were found in a quarter (0.25) of the area with decreasing yield anomaly SD; this proportion (0.25) compared the area with both favorable agro-climatic conditions and decreasing SD (4%) and the area with decreasing yield anomaly SD all agro-climatic conditions taken together (16%) (Fig. 7 D). This proportion (0.25) is 1.7 times higher than the proportion (0.15) between the area with both favorable agro-climatic conditions and increasing yield anomaly SD (3%) and the area with increasing yield anomaly SD all agro-climatic conditions taken together (20%). R&D effect The R&D effect was estimated by comparing two series of yield anomaly simulations, i.e., a series simulated with the actual R&D expenditures in 2000–2019 and another series simulated without R&D investment ( fc and fc.em runs; Fig. 6 ). The results showed that a significant yield-stabilizing effect of R&D was found in a large part (0.44) of the cropland area where a decrease in yield anomaly SD occurred. This proportion (0.44) was derived by comparing the area with a significant R&D effect of 7% and the area with decreasing yield anomaly SD of 16%, both relative to the global harvested area (Fig. 7 B, C, E). This proportion (0.44) is 4.4 times greater than the proportion (0.10) calculated between the area with a significant R&D effect (2%) and the area with increasing yield anomaly SD (20%). On the other hand, the yield-stabilizing effect of R&D was not found to be significant in a dominant proportion (0.85) of the area where an increase in yield anomaly SD occurred. This value (0.85) was derived from the area with a non-significant R&D effect of 17% and the area with increasing yield anomaly SD of 20% (Fig. 7 E), which is 1.7 times greater than the value (0.50) calculated from the area with a non-significant R&D effect (8%) and the area with decreasing yield anomaly SD (16%). Importantly, the area with a significant R&D effect became smaller when estimated over a more recent, shorter period than over the entire study period (Fig. 7 F). [Fig. 6 here] [Fig. 7 here] Sensitivity to different R&D assumptions The estimated area shares for the climate and R&D effect categories were found to be robust to the different assumptions used in the R&D variable calculation. We took into account whether the cumulative R&D expenditure was for agriculture or for the four crops, whether the time lag between research and adoption of technologies by producers was relatively shorter (6 years) or longer (12 years), and whether the technology obsolescence was relatively slower (10% per year) or faster (20% per year) ( Materials and methods ). The main results described in the former subsections (specifically, the fact that R&D is the primary explanatory variable explaining yield anomaly, and unfavorable agro-climatic conditions and/or non-significant R&D effect are proportionately more prevalent in the cropland area experiencing a decrease in yield stability (higher SD)) were robust to the different R&D assumptions (figs. S9–S10). Region and crop-specific characteristics We identified the regions and crop species for which R&D were able to compensate for the negative impacts of climate change on yield stability after 2000. For the four crops taken together (479.8 Mha), R&D had no significant effect in 75% of this area. However, R&D was found to have a significant stabilizing effect on large shares of croplands in several specific crop species and regions, especially for soybean [42% (significant R&D effect) vs. 58% (non-significant)], in Sub-Saharan Africa (41% vs. 59%), East Asia (35% vs. 65%) and North America (29% vs. 71%) (Fig. 8 A). Climate change had no significant effect on yield stability in almost half of this area (43%), but led to unfavorable agro-climatic conditions (higher SD) and favorable agro-climatic conditions (lower SD) in 36% and 21% of the area, respectively. The negative effect of climate change on yield stability was particularly pronounced for maize [45% (unfavorable) vs. 22% (favorable)] and soybean (57% vs. 23%), and in North America (81% vs. 3%), followed by East Asia (33% vs. 18%). On the contrary, the shift to favorable agro-climatic conditions was relatively prominent in South Asia (10% vs. 36%) (Fig. 8 B). [Fig. 8 here] Discussion Yield stability has increased in 16% of the global cropland. Although not mutually exclusive, nearly half of this area is associated with a significant stabilizing effect of R&D, and a quarter of this area is associated with a shift to favorable agro-climatic conditions. In contrast, yield stability has decreased in 20% of the global cropland. Half of this area is associated with a shift to unfavorable agro-climatic conditions, while about nine-tenths of this area is associated with a lack of significant R&D effect, although these areas may overlap. These results provide insights into the association between investment in agricultural R&D, climate change, and resulting changes in yield stability. In addition, based on these results, we could formulate the hypothesis that, in the global cropland area where yield stability has decreased, the progress made in R&D since 2000 has lagged behind the negative impacts of climate change. This question needs to be addressed in future research. Although R&D effect is relatively more prevalent in the area with increasing yield stability, significant R&D effect is also detected in some part of the global cropland with unchanged and decreasing yield stability (8% and 3%, respectively). In this study, the R&D effect is considered significant when R&D expenditures have contributed to a decrease of yield anomaly SD over the study period. Even if the R&D effect is found to be significant, decreasing or unchanged yield stability may occur if the R&D effect is not strong enough to reverse the trend in yield anomaly SD. However, without any R&D effect, the yield stability status in these areas would probably deteriorate further. It should be noted that the R&D effect on yield stability identified here is empirical rather than causal, as some confounding factors potentially correlated with R&D expenditure were not considered in this study. Process-resolving analyses would be useful to determine the causal effects of individual technologies and practices, and confirm our empirical evidence of the yield-stabilizing effect of R&D, which could ultimately support the implementation of effective interventions. Although applied to study average yield and not yield stability, a relevant example is presented in ( 25 ), which decomposes the historical yield increases into individual drivers, including increases in the use of land, water, and other inputs, and increases in the efficiency of input use, with an analysis of policy implications. Important factors of yield anomaly The yield-stabilizing effect of R&D is not totally surprising. For instance, new varieties bred through R&D generally have greater tolerance to sub-optimal growing conditions and higher yield potential than earlier varieties ( 26 – 28 ). In addition to breeding, R&D informs producers about practices to reduce climate risks, such as sowing date, soil, water and fertilizer management, and pests, diseases and weed control. The temperatures in the growing season, as well as occurrences of dry and heavy rainfall days are identified to be important climatic factors explaining yield anomaly. This result is consistent with previous findings ( 29 – 31 ). Under warming, a decrease in the frequency of near-optimal temperatures is accompanied by an increase in the frequency of supra-optimal temperatures, which decreases yield stability ( 32 ). The effects on yield stability from heavily wet conditions identified here are more uncertain than that the effects of dry conditions, although there is a growing body of literature reporting the negative impacts of waterlogging and resulting excessive soil water ( 33 ). The importance of N application in explaining yield anomaly can be explained by previous findings that abiotic stresses due to growing season climate become increasingly important stressors as N inputs increase ( 34 , 35 ). The importance of SOC found in this study confirms previous findings that drought damage is mitigated by SOC build-up ( 36 ). Despite the initial expectation that irrigation would increase yield stability, the contribution of irrigation appears to be limited in this study. The fact that the existence of irrigation equipment does not necessarily ensure water availability during droughts ( 37 ), and that yield sensitivity at regional levels to heat no longer decreases as irrigated area saturates ( 38 ), partly explain the counterintuitive outcome. Limitations and uncertainties There are a few limitations to this study. First, we only consider public expenditures and do not include expenditures by private sector and international organizations. Therefore, the R&D effect estimated in this study may be underestimated to some extent, given that the private sector's share of total agricultural R&D expenditure can be substantial, particularly in high-income countries (30%) compared with low- and middle-income countries (5%) ( 39 ). Moreover, the role of international organizations for agricultural R&D is crucial in low-income countries, and investments in agricultural R&D for one country may generate spillover benefits for neighboring countries and regions with similar agro-climatic conditions. Data on private R&D expenditure are difficult to obtain. Even for public expenditure, the information on the technology on which R&D expenditure was spent is rarely available ( 40 ). Future research should fill these data gaps to distinguish the relative contributions of individual technologies and management practices. Second, a direct comparison with previous estimates of climate-induced changes in yield stability is not possible due to a lack of relevant literature. There is literature showing changes in average yield ( 20 , 23 ). However, none of these studies assessed changes in yield stability and the R&D effect. The remaining limitations and uncertainties are mostly related to the underlying data. Third, this study is based on a single grid yield dataset. Different yield datasets often produce different results. Subnational census statistics are expected to be more reliable than satellite- or model-based yield datasets, including one adopted here. However, the spatial resolution of census statistics is coarser than the resolution of d4PDF climate simulation, which compromises the ability to capture the yield impacts from extreme agro-climatic conditions. In addition, some recent grid yield datasets ( 41 ) use temperature and precipitation as inputs in their yield estimation procedure, positing the circularity issues because we relate yield anomaly to temperature and precipitation. Fourth, attribution outcomes may be influenced by the choice of meteorological forcing datasets and climate models ( 42 , 43 ). Fifth, our analyses are based on the harvested area in 2010, although the area for the four crops in 2022 increased 12% since 2010. This increase in area would include not only a response to increasing demand but also a shift of areas suitable for cultivation in response to the observed warming ( 44 , 45 ). However, we used the grid harvested area map for 2010, the midpoint of the study period, as no annual crop-specific grid harvested area map is available at the global scale. Last, future research should consider more management variables than those considered here, such as pesticide application rates, which are thought to increase yield stability ( 46 ). Materials and methods Experimental design This section describes the data on crop yields, climate (actual, factual and counterfactual) and agricultural R&D expenditure, as well as the method used to detect changes in yield stability. Also included are the machine learning models and the yield impact attribution analyses. Crop yields We obtained country annual yield data from the FAO statistical database for maize, rice, wheat and soybean ( 47 ). In addition, we used the 0.5° grid global dataset of historical yields (GDHY) ( 48 ). In the GHDY dataset, yield data are available for two seasons for maize and rice (major and secondary) and for wheat (winter and spring). For soybean, data are available for a single main season. The GDHY dataset is a hybrid of satellite vegetation index and national census statistics reported by FAO. The information on accumulated net primary production over the crop-specific growing period, calculated using remotely-sensed leaf area index and fraction of photosynthetically-active radiation as well as reanalysis downward shortwave radiation and crop-specific radiation-use efficiency, is combinedly used to spatially disaggregate FAO country average yields to the grid level. Grid yield values are available in 76–92% of the global harvested area with variations by crop and lacking when crop calendar information ( 49 ) is missing. Climate We used 0.5° grid daily data on maximum and minimum air temperatures and precipitation, covering three distinct climate conditions in 2000–2019: actual , factual and counterfactual . The actual climate refers to observed conditions. The factual climate represents a modelled approximation of the actual climate, influenced by both human activities and natural forces (e.g., volcanic eruptions and solar activity). In contrast, the counterfactual climate is a modeled preindustrial climate without any appreciable human impacts on the global climate. The actual climate data were obtained from the global retrospective meteorological forcing data. We used the JRA55-CDFDM-S14FD forcing data (called JCS for simplicity), which is a bias-corrected Japanese 55-year Reanalysis (JRA55) using the cumulative distribution function-based downscaling and bias-correction method (CDFDM) and the reference climatology from S14FD forcing dataset ( 50 ). For the factual and counterfactual climates, we used a large ensemble of long-term climate simulations with and without historical trends in external forcing to the global climate, called d4PDF (the Database for Policy Decision Making for Future Climate Change) ( 17 , 51 ). The climate simulations were performed using the Meteorological Research Institute Atmospheric General Circulation Model, version 3.2 (MRI-AGCM3.2) with a grid interval of 60 km ( 52 ). Each of the factual and counterfactual climate simulations has 100 ensemble members associated with slightly different initial conditions and small perturbations in sea surface temperatures that represent observational uncertainties. Comparisons between the factual and counterfactual climate simulations with large ensemble members enable researchers to conduct a robust attribution analysis of the effects of human activities on average climate conditions, extreme climate events ( 53 – 56 ) and crop yields ( 57 , 58 ). The daily outputs from the climate model were interpolated to the 0.5° resolution and bias-corrected using the CDFDM method along with the reference climatology from S14FD forcing dataset. Agricultural R&D expenditure Public agricultural R&D includes government, higher education and nonprofit entities, but excludes the private for-profit sector ( 39 ). It includes salaries, operating and program costs and capital investments for all entities, excluding the private for-profit sector, involved in agricultural R&D expenditure in a country. For the Organisation for Economic Co-operation and Development (OECD) member economies, we sourced annual data on gross domestic expenditure on R&D expressed as a percentage of gross domestic production (GDP) from the Main Science and Technology Indicators (MSTI) database ( 59 ). Agricultural R&D expenditure for OECD member economies was calculated by multiplying the GDP share of domestic R&D expenditure with the agricultural share of GDP and the total GDP, both obtained from World Bank ( 60 ). For 88 low and middle-income countries, data were obtained from the International Food Policy Research Institute (IFPRI)’s Agricultural Science and Technology Indicators (ASTI) database ( 39 ). The country annual R&D expenditure data collected were adjusted for inflation to be expressed in 2015 USD. Detection of changes in yield stability Changes in yield stability were detected based on the method of Iizumi and Ramankutty ( 8 ). For each grid cell, an annual yield time series (1981–2019) was detrended with a cubic smoothing spline representing time-evolving average yield (fig. S12 A). Although the present study analyzed changes in yield stability between 2000 and 2019, a longer yield time series was used for a more robust identification of yield trends. Yield anomalies (in tonnes per hectare) were computed as the differences between the original yield time series and the fitted trends (fig. S12 B). By applying a 9-year centered moving window, an annual time series of the standard deviation (SD) of yield anomalies was derived. A regression line was fitted to the recent portion (2004–2015) of the yield anomaly SD series, and the estimated slope was tested for a monotonic trend using the two-sided Mann-Kendall trend test ( 61 ) (fig. S12 C). Note that the yield stability trends analyzed here were based on the yield anomalies over the 2000–2019 period (fig. S12 C). A significant positive (negative) slope indicated a decreasing (increasing) yield stability. Monotonic trends in average yield and the coefficient of variation (CV) time series of yield anomaly in the same period were tested in the same manner as the SD time series. To aid visualizing multi-season and multi-crop results as a single map, the detected changes in yield stability for individual seasons and crops were combined according to the following rules. The major season result was selected as the representative when both major (winter) and secondary (spring) seasons were operated. When only one season was operated, the result for that season was considered even when it was the secondary (spring) season. Then the number of crops with specific yield change (either of a significant increase in yield anomaly CV, significant increase in yield anomaly SD and significant increase in average yield), out of four (maize, wheat, rice and soybean), was counted and divided by the number of crops grown for a given location (Fig. 3 D). Random forest yield anomaly models Machine learning models were developed to relate yield anomaly to climate, soil, management and R&D expenditure. We used random forest regression models ( 62 ) and did not explore alternative machine learning approaches as random forests are straightforward to implement, and comparable to other approaches in many applications ( 63 , 64 ). Random forests are generally insensitive to correlations between explanatory variables in terms of predictive accuracy and overfitting ( 63 ), which further supports their use in this study. A random forest algorithm deals with collinearity by randomly sampling a small subset of input variables to build each splitting node of each tree. As a result, the learning algorithm can only choose among a small number of inputs at each splitting node, thus reducing the risk of collinearity. We calculated the z-score values of yield anomaly by dividing annual yield anomaly by the SD of yield anomaly for the same period as that used to estimate yield trends (fig. S12 D). The z-scored yield anomalies were related to climate, soil, management and R&D variables (Table 1). Details of the explanatory variables are provided in the subsequent sections. A specific model was constructed for each crop. Data from major and secondary maize (rice) seasons were combined and used to develop a single maize (rice) model. However, the winter and spring wheat models were developed separately due to potential differences in the yield responses to environmental conditions because winter and spring wheat has different growing seasons and physiological characteristics. Model fitting was conducted using the statistical software R ( 65 ), using the randomForest package ( 66 ) with the following settings: number of trees (ntree = 500), number of predictors sampled at each split (mtry = 3), and minimum size of terminal nodes (nodesize = 5). These hyperparameter settings were selected based on the package’s default values. For model validation purposes, we used the leave-one-year-out cross validation technique and assessed the Pearson’s correlation coefficients and the root mean squared errors calculated between the actual and simulated yield anomalies. Furthermore, the model evaluation was done to see whether the RF models with d4PDF factual climate simulation could reproduce the actual global pattern of the changes in yield anomaly SD derived from the grid yield dataset. We also checked whether values of d4PDF factual and counterfactual climate simulations did not exceed the ranges covered by the training data (i.e., the JCS retrospective meteorological forcing) to ensure the reliability of the RF models. Last, we assessed the importance of the explanatory variables. We used two methods to rank model input variables according to their importance; the permutation method and the node purity method that rely on the residual sum of squares and the Gini index as importance metrics, respectively. Climate variables Temperature and precipitation were considered to be the primary climatic drivers of yield anomalies. To account the effects of short-term extreme temperatures, relative counts of hourly temperature exposure for a growing season were calculated using 9 bins with constant 5°C intervals (T < 0, 0 ≤ T < 5, 5 ≤ T < 10, …, 25 ≤ T < 30, 30 ≤ T < 35, and 35 ≤ T). Hourly temperature values were obtained by fitting a sine curve to daily maximum and minimum temperatures. Relative counts of daily precipitation exposure during the growing season were also calculated using 5 bins with unequal intervals (P < 1, 1 ≤ P < 5, 5 ≤ P < 20, 20 ≤ P < 30, and 30 ≤ P). The growth periods were fixed based on the reported planting and harvesting dates ( 49 ). For winter wheat only, the period before the completion of vernalization was excluded from the calculation by identifying the first date at which the fraction of growing-degree days to the crop total thermal requirements exceeds 0.1; this was estimated using the crop phenology model with the crop-specific base (0°C) and maximum temperature (26°C) for winter wheat. These climate variables were calculated for each of the actual, factual and counterfactual climate conditions. Soil variable We focused on soil organic carbon (SOC) content as a main soil variable, given the fact that in arid and semiarid regions of the world, high levels of SOC improve soil water holding capacity and moderate drought damage ( 36 ). The topsoil (0–30 cm) SOC data were obtained from the Regridded Harmonized World Soil Database v1.2 ( 67 ) and aggregated to the 0.5° resolution. In addition to SOC, we also used fractions of clay, silt and sand, bulk density, cation exchange capacity of the clay fraction, pH, and available water storage capacity, all at the topsoil, to consider possible effects of soil characteristics on yield anomaly. Management variables We considered two variables related to producer’s crop management, i.e., irrigation and N fertilizer. Irrigation mitigates drought damages and heat stress through evaporative cooling of the crop canopy ( 68 ). The irrigation intensity used in this study considered changes due to the expansion of the irrigation-equipped area based on the global historical irrigation dataset for the period of 1900–2005 (HID) ( 69 ) and crop-specific variations based on the global monthly irrigated and rainfed crop areas around the year 2000 (MIRCA2000) ( 70 ). The data were aggregated to the 0.5° resolution. Due to the lack of data, irrigation intensity values for the period 2006–2019 were extrapolated using a linear regression using the data from 1998 to 2005, although this assumption may not be totally accurate for some regions experiencing decreases in irrigation area, such as Russia ( 37 ). In addition to irrigation, we used the grid annual N application rate for each crop ( 71 ) as it is known that abiotic stresses become important as N inputs increase ( 34 , 35 ). R&D variable We incorporated the effect of agriculture R&D and resulting technologies and management practices that help stabilize yields into the RF models. To do so, an economic indicator was used. The indicator, R i , t , is the accumulated agricultural R&D expenditure for country i from 1995 to year t and expressed in constant 2015 USD. The accumulation took into account the duration of research and the obsolescence of technology: $$\:{R}_{i,\:t}={E}_{i,t-6}+\left(1-\delta\:\right){R}_{i,t-1}$$ 1 where E is the country annual agricultural R&D expenditure (constant 2015 USD); and δ is the annual rate of technology obsolescence (= 0.1 or 10% per year) ( 71 ). The lag time between research and technology adoption by producers was set at six years, although lag time may vary by country and over time ( 72 ). As the R&D variable is of key importance in this study, we performed the sensitivity analysis on this variable, as described in the latter subsection. Attribution analyses Simulation experiments were conducted using the RF models and d4PDF factual and counterfactual climate simulations. As shown in Fig. 6 , the fc and ct runs compared the factual and counterfactual climates and determine the climate effect to yield stability under actual R&D conditions. The fc and fc.em runs compared yield stability with and without R&D after 2000 under climate change. Climate effect The changes in yield stability associated with human-induced climate change was detected through a comparison of fc and ct runs (Fig. 6 ). Using 100 samples of the slope value of yield anomaly SD derived from the RF models forced by d4PDF climate simulation, a histogram was derived for each of the factual and counterfactual yields (fig. S13. When the average slope value of the factual yield anomaly SD is higher (lower) than that of the counterfactual yield anomaly SD, it indicates that climate change contributed to decreasing (increasing) yield stability in 2000–2019. The two-sided Wilcoxon rank sum test ( 73 ) was used to test the null hypothesis that the two types of slopes were equal. R&D effect The R&D effect after 2000 to yield stability under climate change was identified through a comparison of the fc and fc.em runs (Fig. 6 ). We assumed that the fc run, which included the time-dependent R&D and management variables, reflect advances in R&D and associated practices. Then, we performed a counterfactual run in which advances in R&D and management practices are lacking to varying degrees. To do this, we fixed the R&D and management variables to constant values. In the 1st model specification (M1), values of the R&D and management variables were set at the 2000 levels to calculate yield anomalies from 2000 to 2019, while in the 5th model specification (M5), they were set equal to their yearly values until 2004 and were held constant thereafter (fig. S14). This model setting resulted in 20 model specifications per crop (M1 to M20). For each year in which advances in R&D and practices were eliminated in the RF models, a histogram of the slope value for the corresponding model specification was derived and compared with a histogram of the slope value for actual R&D and practices, each consisting of 100 members from the RF models forced by d4PDF factual climate simulation (fig. S15). When the average slope value for actual R&D and practice was significantly smaller than that without advances in R&D and practice, it indicates that changes in R&D and practices after 2000 increased yield stability. The absence of a significant difference indicates that changes in R&D and practices since 2000 did not have a significant effect on increasing yield stability under climate change. Sensitivity analysis We examined the sensitivity of variable importance and attribution results to the different assumptions used in the R&D variable calculation. The country annual R&D expenditure data collected were for agriculture and not limited to the four crops considered here. Therefore, we first examined if the use of R&D expenditure for these crops, not for agriculture, affected results. To that end, we assumed that annual R&D expenditure designated for the four crops was proportional to the share of value derived from production of the four crops to total value of agricultural production. We collected country annual data on value of production for the four crops and entire agriculture from FAO ( 47 ) and derived them for each country and year. Second, we considered shorter (6 years) and longer (12 years) time lags between research and adoption of technologies by producers. Finally, we considered relatively slower (10% per year) and faster (20% per year) technology obsolescence by setting δ value in Eq. ( 1 ) to 0.1 and 0.2, respectively. We performed five experiments in the sensitivity analysis: Ex. 1—R&D expenditure for agriculture (c0), shorter time lag between research and adoption (r6) and slower technology obsolescence (o1) (denoted “c0r6o1”); Ex. 2—R&D expenditure for agriculture, longer time lag (r12) and slower obsolescence (“c0r12o1”); Ex. 3—R&D expenditure for the four crops (c1), shorter time lag and slower obsolescence (“c1r6o1”); Ex. 4—R&D expenditure for the four crops, longer time lag and slower obsolescence (“c1r12o1”); and Ex. 5—R&D expenditure for the four crops, shorter time lag and faster obsolescence (“c1r6o2”). For each experiment, the RF models were built and then the variable importance and the area shares for the climate and R&D effect categories were estimated. In the main text, we presented the results from Ex. 1 (c0r6o1). We found that the assumptions used in the R&D variable calculation affected the level of but not the pattern of cumulative R&D expenditure (fig. S16). The difference in the level of cumulative R&D expenditure was absorbed in the RF models, as ordinary regression models would do, and did not affect the variable importance and attribution results. Declarations Funding TI was supported by the Environment Research and Technology Development Fund (Grant Number: JPMEERF23S21120) of the Environmental Restoration and Conservation Agency provided by the Ministry of the Environment of Japan; and Grants-in-Aid for Scientific Research were received from the Japan Society for the Promotion of Science (Grant Numbers: 22H00577 and 23H00351). YI and HS were supported by the Advanced Studies of Climate Change Projection (SENTAN, Grant Number: JPMXD0722680395) of the Ministry of Education, Culture, Sports, Science and Technology of Japan. DM and TI were supported by INRAE-NARO Linkage Call 2023. Author contributions TI conceptualized, developed methodology, and conducted formal analysis; TS, YM, KO, and TT curated crop data and assisted with the investigation; HS and YI curated climate data and assisted with the investigation; DM assisted with methodology development; all authors contributed to writing the original draft. Data availability The country yield and area harvested data from the FAO statistical database can be obtained from https://www.fao.org/faostat/en/#home. Gross domestic R&D spending data for some developed countries can be obtained from the OECD’s MSTI database: https://data-explorer.oecd.org/. Agricultural R&D expenditure data for some developing countries are accessible from the IFPRI’s ASTI database: https://doi.org/10.7910/DVN/9OXBIB. The grid yield dataset is accessible at the Data Integration and Analysis System (DIAS): https://doi.org/10.20783/DIAS.564. The actual climate data and the bias-corrected d4PDF factual and counterfactual climate simulations can be respectively obtained from https://doi.org/10.20783/DIAS.523 and https://doi.org/10.20783/DIAS.544. References FAO, IFAD, UNICEF, WFP, WHO. The State of Food Security and Nutrition in the World 2020 (FAO, 2020). https://doi.org/10.4060/ca9692en Nóia Júnior, R. de S., Ewert, F., Webber, H., Martre, P., Hertel, T.W., van Ittersum, M.K., Asseng, S. Needed global wheat stock and crop management in response to the war in Ukraine. Glob. Food Secur. 35 , 100662 (2022). https://doi.org/10.1016/j.gfs.2022.100662. Aday, S., Aday, M.S. Impact of COVID-19 on the food supply chain. Food Qual. Saf. 4 , 2020, 167–180. https://doi.org/10.1093/fqsafe/fyaa024 Gaupp, F. Extreme Events in a Globalized Food System. One Earth 2 , 518–521 (2020). https://doi.org/10.1016/j.oneear.2020.06.001 Zhu, P., Burney, J., Chang, J. et al. Warming reduces global agricultural production by decreasing cropping frequency and yields. Nat. Clim. Chang 12 , 1016–1023 (2022). https://doi.org/10.1038/s41558-022-01492-5 Wei, D. et al. Key role of planted and harvested area fluctuations in US crop production shocks. Nat. Sustain. 6 , 1177–1185 (2023). https://doi.org/10.1038/s41893-023-01152-2 Bezner Kerr, R. et al. “Food, Fibre, and Other Ecosystem Products” in Climate Change 2022: Impacts, Adaptation and Vulnerability (Cambridge Univ. Press, 2022). doi:10.1017/9781009325844.007 Iizumi, T., Ramankutty, N. Changes in yield variability of major crops for 1981–2010 explained by climate change. Environ. Res. Lett. 11 , 034003 (2016). https://dx.doi.org/10.1088/1748-9326/11/3/034003 Ben-Ari, T., Boé, J., Ciais, P. et al. Causes and implications of the unforeseen 2016 extreme yield loss in the breadbasket of France. Nat Commun. 9 , 1627 (2018). https://doi.org/10.1038/s41467-018-04087-x Thirumalai, K. et al. Extreme temperatures in Southeast Asia caused by El Niño and worsened by global warming. Nat. Commun. 8 , 15531 (2017). https://doi.org/10.1038/ncomms15531 Beillouin, D., Schauberger, B., Bastos, A., Ciais, P., Makowski, D. Impact of extreme weather conditions on European crop production in 2018. Phil. Trans. R. Soc. B 375 , 20190510 (2020). https://doi.org/10.1098/rstb.2019.0510 White, R. H. et al. The unprecedented Pacific Northwest heatwave of June 2021. Nat. Commun. 14 , 727 (2023). https://doi.org/10.1038/s41467-023-36289-3 Cooper, M., Messina, C. D. Breeding crops for drought-affected environments and improved climate resilience. The Plant Cell 35 , 162–186 (2023). https://doi.org/10.1093/plcell/koac321 Berrang-Ford, L. et al. A systematic global stocktake of evidence on human adaptation to climate change. Nat. Clim. Chang. 11 , 989–1000 (2021). https://doi.org/10.1038/s41558-021-01170-y Lobell, D. B. et al. Greater Sensitivity to Drought Accompanies Maize Yield Increase in the U.S. Midwest. Science 344 , 516–519 (2014). doi:10.1126/science.1251423 Miller, S. J., Clarke, E., Mathews, S. L. Little evidence of avoided yield loss in US corn when short-term forecasts correctly predict extreme heat. Environ. Res. Lett. 18 , 124041 (2023). https://dx.doi.org/10.1088/1748-9326/ad0bd5 Mizuta, R. et al. Over 5,000 years of ensemble future climate simulations by 60-km global and 20-km regional atmospheric models. Bull. Amer. Meteor. Soc. 98 , 1383–1398 (2017). https://doi.org/10.1175/BAMS-D-16-0099.1 Otto, F. E. L. Attribution of weather and climate events. Annu. Rev. Environ. Resour. 42 , 627–646 (2017). https://doi.org/10.1146/annurev-environ-102016-060847 Gillett, N. P. et al. Constraining human contributions to observed warming since the pre-industrial period. Nat. Clim. Chang 11 , 207–212 (2021). https://doi.org/10.1038/s41558-020-00965-9 Romanovska, P., Undorf, S., Schauberger, B., Duisenbekova, A., Gornott, C. 2024. Human-induced climate change has decreased wheat production in northern Kazakhstan. Environmental Research: Climate 3 , 031005. https://dx.doi.org/10.1088/2752-5295/ad53f7 Funk, C. et al. Anthropogenic enhancement of moderate-to-strong El Niño events likely contributed to drought and poor harvests in southern Africa during 2016. Bull. Amer. Meteor. Soc. 99 , S91–S96 (2018). https://doi.org/10.1175/BAMS-D-17-0112.1 Verschuur, J. et al. Climate change as a driver of food insecurity in the 2007 Lesotho-South Africa drought. Sci. Rep. 11 , 3852 (2021). https://doi.org/10.1038/s41598-021-83375-x FAO. The Impact of Disasters on Agriculture and Food Security 2023 – Avoiding and reducing losses through investment in resilience (FAO, 2023). https://doi.org/10.4060/cc7900en Ray, D.K., Mueller, N.D., West, P.C., Foley, J.A. (2013). Yield Trends Are Insufficient to Double Global Crop Production by 2050. PLoS ONE 8 , e66428. https://doi.org/10.1371/journal.pone.0066428 Fuglie, K.O., Morgan, S., Jelliffe, J. (2024). World agricultural production, resource use, and productivity, 1961–2020 (Report No. EIB-268). U.S. Department of Agriculture, Economic Research Service. DeBruin, J.L., Schussler, J.R., Mo, H. and Cooper, M. (2017), Grain Yield and Nitrogen Accumulation in Maize Hybrids Released during 1934 to 2013 in the US Midwest. Crop Science, 57: 1431-1446. https://doi.org/10.2135/cropsci2016.08.0704 H.E. Ahrends et al. 2018 Genetic yield gains of winter wheat in Germany over more than 100 years (1895–2007) under contrasting fertilizer applications. Environmental Research Letters 13, 104003. https://dx.doi.org/10.1088/1748-9326/aade12 Shew, A.M., Tack, J.B., Nalley, L.L. et al. Yield reduction under climate warming varies among wheat cultivars in South Africa. Nat Commun 11, 4408 (2020). https://doi.org/10.1038/s41467-020-18317-8 Schmidt, M., Felsche, E. The effect of climate change on crop yield anomaly in Europe. Clim. Resilience Sustainability 3 , e61 (2024). https://doi.org/10.1002/cli2.61 Kim, W., Iizumi, T., Nishimori, M. Global patterns of crop production losses associated with droughts from 1983 to 2009. J. Appl. Meteorol. Climatol. 58 , 1233-1244 (2019). https://doi.org/10.1175/JAMC-D-18-0174.1 Vogel, E. et al. The effects of climate extremes on global agricultural yields. Environ. Res. Lett. 14 , 054010 (2019). https://dx.doi.org/10.1088/1748-9326/ab154b Tigchelaar, M., Battisti, D.S., Naylor, R.L., Ray, D.K. Future warming increases probability of globally synchronized maize production shocks. Proc. Natl. Acad. Sci. USA 115 , 6644–6649. https://doi.org/10.1073/pnas.1718031115 Kim, Y.-U., Webber, H., Adiku, S.G.K., Nóia Júnior, R.d.S., Deswarte, J.-C., Asseng, S., Ewert, F., 2024. Mechanisms and modelling approaches for excessive rainfall stress on cereals: Waterlogging, submergence, lodging, pests and diseases. Agric. For. Meteorol. 344, 109819. https://doi.org/10.1016/j.agrformet.2023.109819 Müller, C., et al. Global patterns of crop yield stability under additional nutrient and water inputs. PLoS ONE 13 , e0198748 (2018). https://doi.org/10.1371/journal.pone.0198748 Ahrends, H. E. et al. Nutrient supply affects the yield stability of major European crops—a 50 year study. Environ. Res. Lett. 16 , 014003 (2021). https://dx.doi.org/10.1088/1748-9326/abc849 Qiao, L. et al. Soil quality both increases crop production and improves resilience to climate change. Nat. Clim. Chang. 12 , 574–580 (2022). https://doi.org/10.1038/s41558-022-01376-8 Mehta, P. et al. Half of twenty-first century global irrigation expansion has been in water-stressed regions. Nat. Water 2 , 254–261 (2024). https://doi.org/10.1038/s44221-024-00206-9 Zaveri, E., Lobell, D. B. The role of irrigation in changing wheat yields and heat sensitivity in India. Nat. Commun. 10 , 4144 (2019). https://doi.org/10.1038/s41467-019-12183-9 IFPRI. Agricultural Science and Technology Indicators: 2019 Global Food Policy Report Annex Table 1 (IFPRI, 2019). https://doi.org/10.7910/DVN/9OXBIB Nalley L, Dixon B, Chaminuka P, Naledzani Z, Coale MJ (2018) The role of public wheat breeding in reducing food insecurity in South Africa. PLoS ONE 13 , e0209598. https://doi.org/10.1371/journal.pone.0209598 Wu, H. et al. AsiaRiceYield4km: seasonal rice yield in Asia from 1995 to 2015. Earth Syst. Sci. Data 15 , 791–808 (2023). https://doi.org/10.5194/essd-15-791-2023 Ruane, A.C. et al. 2021. Strong regional influence of climatic forcing datasets on global crop model ensembles. Agricultural and Forest Meteorology, 300, 108313. https://doi.org/10.1016/j.agrformet.2020.108313 Mengel, M., Treu, S., Lange, S., and Frieler, K.: ATTRICI v1.1 – counterfactual climate for impact attribution, Geosci. Model Dev., 14, 5269–5284, https://doi.org/10.5194/gmd-14-5269-2021, 2021. Sloat, L.L., Davis, S.J., Gerber, J.S. et al. Climate adaptation by crop migration. Nat Commun 11, 1243 (2020). https://doi.org/10.1038/s41467-020-15076-4 G. Zhang et al. 2017 Spatiotemporal patterns of paddy rice croplands in China and India from 2000 to 2015. Science of The Total Environment, 579, 82-92. https://doi.org/10.1016/j.scitotenv.2016.10.223. Maggi, F. et al. PEST-CHEMGRIDS, global gridded maps of the top 20 crop-specific pesticide application rates from 2015 to 2025. Sci. Data 6 , 170 (2019). https://doi.org/10.1038/s41597-019-0169-4 FAO. FAOSTAT (FAO, 2023). https://www.fao.org/faostat/en/#home Iizumi, T., Sakai, T. The global dataset of historical yields for major crops 1981–2016. Sci. Data 7 , 97 (2020). https://doi.org/10.1038/s41597-020-0433-7 Sacks, W. J., Deryng, D., Foley, J. A. Ramankutty, N. Crop planting dates: an analysis of global patterns. Glob. Ecol. Biogeogr. 19 , 607-620 (2010). https://doi.org/10.1111/j.1466-8238.2010.00551.x Iizumi, T. et al. Increasing heavy rainfall events and associated excessive soil water threaten a protein-source legume in dry environments of West Africa. Agric. For. Meteorol. 344 , 109783 (2024). https://doi.org/10.1016/j.agrformet.2023.109783 Ishii, M., Mori, N. d4PDF: large-ensemble and high-resolution climate simulations for global warming risk assessment. Prog. Earth Planet. Sci. 7 , 58 (2020). https://doi.org/10.1186/s40645-020-00367-7 Mizuta, R. et al. Climate simulations using MRI-AGCM3.2 with 20-km grid. J. Meteorol. Soc. Jpn 90A , 233-258 (2012). https://doi.org/10.2151/jmsj.2012-A12 Shiogama, H. et al. Attributing historical changes in probabilities of record-breaking daily temperature and precipitation extreme events. SOLA 12 , 225-231 (2016). https://doi.org/10.2151/sola.2016-045 Imada, Y., Watanabe, M., Kawase, H., Shiogama, H., Arai, M. The July 2018 high temperature event in Japan could not have happened without human-induced global warming. SOLA 15A , 8-12 (2019). https://doi.org/10.2151/sola.15A-002 Kawase, H. et al. Contribution of historical global warming to local-scale heavy precipitation in western Japan estimated by large ensemble high-resolution simulations. J. Geophys. Res. Atmos. 124 , 6093–6103 (2019). https://doi.org/10.1029/2018JD030155 Imada, Y. et al. Advanced risk-based event attribution for heavy regional rainfall events. npj Clim. Atmos. Sci. 3 , 37 (2020). https://doi.org/10.1038/s41612-020-00141-y Iizumi, T. et al. Crop production losses associated with anthropogenic climate change for 1981–2010 compared with preindustrial levels. Int. J. Climatol. 38 , 5405–5417 (2018). https://doi.org/10.1002/joc.5818 Sultan, B., Defrance, D., Iizumi, T. Evidence of crop production losses in West Africa due to historical global warming in two crop models. Sci. Rep. 9 , 12834 (2019). https://doi.org/10.1038/s41598-019-49167-0 OECD. MSTI database (OECD, 2023). https://stats.oecd.org/Index.aspx?DataSetCode=MSTI_PUB World Bank. World Development Indicators (World Bank, 2023) https://databank.worldbank.org/reports.aspx?source=world-development-indicators# Hipel, K. W., McLeod, A. I. Time Series Modelling of Water Resources and Environmental Systems (Elsevier, 1994). Breiman, L. Random forests. Mach Learn. 45 , 5–32 (2001). https://doi.org/ 10.1023/A: 10109 33404 324 Folberth, C. et al. Spatio-temporal downscaling of gridded crop model yield estimates based on machine learning. Agric. For. Meteorol. 264 , 1–15 (2019). https://doi.org/10.1016/j.agrformet.2018.09.021 Chen, M., Guilpart, N., Makowski, D. Comparison of methods to aggregate climate data to predict crop yield: an application to soybean. Environmental Research Letters 19, 054049 (2024). https://dx.doi.org/10.1088/1748-9326/ad42b5 R Core Team. R: a language and environment for statistical computing (R Foundation for Statistical Computing, 2023) Liaw, A., Wiener, M. Classification and regression by randomForest. R News 2 , 18–22 (2002). https://journal.r-project.org/articles/RN-2002-022/ Wieder, W. R., Boehnert, J., Bonan, G. B. Evaluating soil biogeochemistry parameterizations in Earth system models with observations. Glob. Biogeochem. Cycles 28 , 211–222 (2014). https://doi.org/10.1002/2013GB004665 McDermid, S. et al. Irrigation in the Earth system. Nat. Rev. Earth Environ. 4 , 435–453 (2023). https://doi.org/10.1038/s43017-023-00438-5 Siebert, S. et al. A global data set of the extent of irrigated land from 1900 to 2005. Hydrol. Earth Syst. Sci. 19 , 1521–1545 (2015). https://doi.org/10.5194/hess-19-1521-2015 Portmann, F. T., Siebert, S., Döll, P. MIRCA2000—Global monthly irrigated and rainfed crop areas around the year 2000: A new high‐resolution data set for agricultural and hydrological modeling. Glob. Biogeochem. Cycles 24 , GB1011 (2010). https://doi.org/10.1029/2008GB003435 Iizumi, T. et al. Responses of crop yield growth to global temperature and socioeconomic changes. Sci. Rep. 7 , 7800 (2017). https://doi.org/10.1038/s41598-017-08214-4 Baldos, U.L.C., Fuglie, K.O., Hertel, T.W. The research cost of adapting agriculture to climate change: A global analysis to 2050. Agric. Econ. 51 , 207–220 (2020). https://doi.org/10.1111/agec.12550 Hollander, M., Wolfe, D. A., Chicken, E. Nonparametric Statistical Methods, Third Edition (John Wiley & Sons, 2015). doi:10.1002/9781119196037 Yu, Q. et al. A cultivated planet in 2010 – Part 2: The global gridded agricultural-production maps. Earth Syst. Sci. Data 12 , 3545–3572 (2020). https://doi.org/10.5194/essd-12-3545-2020 Table Table 1 is available in the Supplementary Files section Additional Declarations The authors declare no competing interests. Supplementary Files sm20250304.docx Supplementary material Figs. S1 to S16 table11.docx Cite Share Download PDF Status: Published Journal Publication published 21 Mar, 2025 Read the published version in PNAS Nexus → Version 2 posted You are reading this latest preprint version Show more versions Research Square lets you share your work early, gain feedback from the community, and start making changes to your manuscript prior to peer review in a journal. As a division of Research Square Company, we’re committed to making research communication faster, fairer, and more useful. We do this by developing innovative software and high quality services for the global research community. Our growing team is made up of researchers and industry professionals working together to solve the most critical problems facing scientific publishing. Also discoverable on Platform About Our Team In Review Editorial Policies Advisory Board Help Center Resources Author Services Accessibility API Access RSS feed Manage Cookie Preferences © Research Square 2026 | ISSN 2693-5015 (online) Privacy Policy Terms of Service Do Not Sell My Personal Information {"props":{"pageProps":{"initialData":{"identity":"rs-3895355","acceptedTermsAndConditions":true,"allowDirectSubmit":true,"archivedVersions":[],"articleType":"Research Article","associatedPublications":[],"authors":[{"id":434945936,"identity":"9f359042-612e-434e-86dd-c9df6c1cb8b8","order_by":0,"name":"Toshichika Iizumi","email":"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAZAAAAAyAQMAAABI0h/eAAAABlBMVEX///8AAABVwtN+AAAACXBIWXMAAA7EAAAOxAGVKw4bAAAAz0lEQVRIiWNgGAWjYBACNmYQacAgx8B88AFU7ABRWgyMGdiSDYjTAgUGiQ0ILQQAHzv7M6kbBX/S+9uY2T78YKiVY2A8i98aoMPSpHMMDHJnHGNmntnDcNyYgeFcAiEtx8BaGu73Hwb661hiA8MZ/C5kY2ZsA2lJlwfaQqwWZjaQlgQDiJYaYrSwMVvnGBgbbgRqYewxOGDMRsgv8v3HH97O+SMnLwfUwvCjok6OX4JAiKEBg8MMbBJnSNHBwFDHwMDfQ5qWUTAKRsEoGPYAAKbUOYN1O216AAAAAElFTkSuQmCC","orcid":"https://orcid.org/0000-0002-0611-4637","institution":"Institute for Agro-Environmental Sciences, National Agriculture and Food Research Organization","correspondingAuthor":true,"prefix":"","firstName":"Toshichika","middleName":"","lastName":"Iizumi","suffix":""},{"id":434945937,"identity":"5c32e7d8-82c2-441f-92b8-0814206deb09","order_by":1,"name":"Toru Sakai","email":"","orcid":"","institution":"Japan International Research Center for Agricultural Sciences","correspondingAuthor":false,"prefix":"","firstName":"Toru","middleName":"","lastName":"Sakai","suffix":""},{"id":434945938,"identity":"a5693823-ff95-4070-86f6-777dc47ffbdf","order_by":2,"name":"Yoshimitsu Masaki","email":"","orcid":"","institution":"Institute for Agro-Environmental Sciences, National Agriculture and Food Research Organization","correspondingAuthor":false,"prefix":"","firstName":"Yoshimitsu","middleName":"","lastName":"Masaki","suffix":""},{"id":434945939,"identity":"c2534cab-e620-4f22-b737-668684cdcf1c","order_by":3,"name":"Kei Oyoshi","email":"","orcid":"","institution":"Earth Observation Research Center, Japan Aerospace Exploration Agency","correspondingAuthor":false,"prefix":"","firstName":"Kei","middleName":"","lastName":"Oyoshi","suffix":""},{"id":434945940,"identity":"8b9ac755-0773-4283-b53e-a754450b5762","order_by":4,"name":"Takahiro Takimoto","email":"","orcid":"","institution":"Institute for Agro-Environmental Sciences, National Agriculture and Food Research Organization","correspondingAuthor":false,"prefix":"","firstName":"Takahiro","middleName":"","lastName":"Takimoto","suffix":""},{"id":434945941,"identity":"e54d0f31-2031-425d-9cf3-df6a89eb476d","order_by":5,"name":"Hideo Shiogama","email":"","orcid":"","institution":"National Institute for Environmental Studies","correspondingAuthor":false,"prefix":"","firstName":"Hideo","middleName":"","lastName":"Shiogama","suffix":""},{"id":434945942,"identity":"939eee23-921b-4d21-bf6d-31119d8a4a91","order_by":6,"name":"Yukiko Imada","email":"","orcid":"","institution":"Atmosphere and Ocean Research Institute, The University of Tokyo","correspondingAuthor":false,"prefix":"","firstName":"Yukiko","middleName":"","lastName":"Imada","suffix":""},{"id":434945943,"identity":"e6508960-7786-4d2a-b108-c360b22277cc","order_by":7,"name":"David Makowski","email":"","orcid":"","institution":"Université Paris-Saclay, INRAE, AgroParisTech","correspondingAuthor":false,"prefix":"","firstName":"David","middleName":"","lastName":"Makowski","suffix":""}],"badges":[],"createdAt":"2024-01-24 22:45:41","currentVersionCode":2,"declarations":"","doi":"10.21203/rs.3.rs-3895355/v2","doiUrl":"https://doi.org/10.21203/rs.3.rs-3895355/v2","draftVersion":[],"editorialEvents":[{"content":"https://doi.org/10.1093/pnasnexus/pgaf099","type":"published","date":"2025-03-22T00:00:00+00:00"}],"editorialNote":"","failedWorkflow":false,"files":[{"id":79380285,"identity":"cfbea64b-91e5-4d39-a5fa-6336b35f32e0","added_by":"auto","created_at":"2025-03-18 09:18:32","extension":"png","order_by":1,"title":"Figure 1","display":"","copyAsset":false,"role":"figure","size":114650,"visible":true,"origin":"","legend":"\u003cp\u003eRelationship between R\u0026amp;D, yield increase and stability. (A) Natural log of the cumulative annual public R\u0026amp;D expenditure for agriculture since 1995 versus average yield. (B) Average yield versus the yield anomaly SD. (C) Average yield versus the yield anomaly CV (CV=SD/average). These results were obtained over 2000–2019 for 115 crop-producing countries. The selected countries have an annual production of 0.1 Mt or more for at least one of the crops considered here (maize, rice, wheat and soybean). The cumulative R\u0026amp;D expenditure is in units of constant 2015 million USD ($M, log-transformed). The size and color of the symbols indicate the total area harvested and accumulated agricultural R\u0026amp;D expenditure, respectively. The solid lines show the best-fitted nonparametric quantile regression for the median that divide data into two equal portions. For A to C, all slope values are significant (p\u0026lt;0.01; F test).\u003c/p\u003e","description":"","filename":"1.png","url":"https://assets-eu.researchsquare.com/files/rs-6174385/v1/1406ad6d455bf45ab7aa60aa.png"},{"id":79380292,"identity":"cfcea686-8a6e-493a-a1f8-6c6a31836380","added_by":"auto","created_at":"2025-03-18 09:18:32","extension":"png","order_by":2,"title":"Figure 2","display":"","copyAsset":false,"role":"figure","size":150682,"visible":true,"origin":"","legend":"\u003cp\u003ePatterns of agricultural R\u0026amp;D expenditure. (A) Global map of trends in cumulative R\u0026amp;D expenditure in 2000–2019. (B) Country annual R\u0026amp;D expenditure in constant 2015 USD. (C) Country cumulative R\u0026amp;D expenditure since 1995. Classes of increasing (decreasing) cumulative R\u0026amp;D expenditure are classified by growth rate, i.e., “rapidly” and “slowly”, respectively, if a given annual rate is faster and equal to or slower than the median value. The median annual rate is 4.8% for the increasing trend and –1.5% for the decreasing trend, both are relative to the 2015 level. The colors in the time series (B, C) correspond to those in the map (A).\u003c/p\u003e","description":"","filename":"2.png","url":"https://assets-eu.researchsquare.com/files/rs-6174385/v1/023028d3e6628abb5a7a59b7.png"},{"id":79380294,"identity":"faf51ecf-17d5-429a-8735-ec82fcfedfdb","added_by":"auto","created_at":"2025-03-18 09:10:32","extension":"png","order_by":3,"title":"Figure 3","display":"","copyAsset":false,"role":"figure","size":445222,"visible":true,"origin":"","legend":"\u003cp\u003eRelative counts of crops with specific pattern of yield change. (A) Increase in yield anomaly CV. (B) Increase in yield anomaly SD. (C) Increase in average yield. (D) Number of crops considered here that are harvested in 2010 (the midpoint of the study period). Statistical significances of monotonic trends in 2000–2019 are examined using the two-sided Mann-Kendall test. Maize, rice, wheat and soybean are considered together. The results for each crop are available in figs. S1–S3.\u003c/p\u003e","description":"","filename":"3.png","url":"https://assets-eu.researchsquare.com/files/rs-6174385/v1/06691e58a7e49106f33944ce.png"},{"id":79380296,"identity":"ed39f729-a33b-45a3-94bf-1bbf414549d9","added_by":"auto","created_at":"2025-03-18 09:10:32","extension":"png","order_by":4,"title":"Figure 4","display":"","copyAsset":false,"role":"figure","size":373527,"visible":true,"origin":"","legend":"\u003cp\u003eActual and reproduced changes of yield anomaly SD. (A) Actual pattern (same as in Fig. 3 B). (B) Reproduction obtained with the trained random forest models using the retrospective meteorological forcing data. (C) Reproduction using d4PDF factual climate simulation. In C, a cell is classified to as “increasing” if at least one-third or more of the 100 climate-simulation members have a significant increasing trend in the factual yield anomaly SD series. Cohen’s kappa statistics (κ) calculated for each of B and C against A are shown as an indication of the similarity of the geographic pattern between a given pair.\u003c/p\u003e","description":"","filename":"4.png","url":"https://assets-eu.researchsquare.com/files/rs-6174385/v1/a5a4a5cd67dd6895a468db86.png"},{"id":79380298,"identity":"e7319f5c-3918-400d-b62d-1d42b0e897ca","added_by":"auto","created_at":"2025-03-18 09:10:33","extension":"png","order_by":5,"title":"Figure 5","display":"","copyAsset":false,"role":"figure","size":126502,"visible":true,"origin":"","legend":"\u003cp\u003eRanking of variable importance according to the random forest models. The bars show the average rank calculated over the crops (the closer the rank is to one, the more important the variable is in explaining yield anomaly). The thin lines indicate the standard deviation of the rank between the crops. Two variable importance methods (permutation and node purity) are used to account for uncertainty. The results of individual crops are available in fig. S6.\u003c/p\u003e","description":"","filename":"5.png","url":"https://assets-eu.researchsquare.com/files/rs-6174385/v1/71efbc19bccbe4fdfc41a6ac.png"},{"id":79380299,"identity":"9f223d1d-4e55-46dd-9a29-25df35e73fce","added_by":"auto","created_at":"2025-03-18 09:18:32","extension":"png","order_by":6,"title":"Figure 6","display":"","copyAsset":false,"role":"figure","size":143544,"visible":true,"origin":"","legend":"\u003cp\u003eAttribution of changes in yield stability. Model simulations were conducted to estimate the effects of climate and R\u0026amp;D on the changes in yield stability in 2000–2019.\u003c/p\u003e","description":"","filename":"6.png","url":"https://assets-eu.researchsquare.com/files/rs-6174385/v1/178e0b27594c8ed3b924ff31.png"},{"id":79380301,"identity":"0b9908d9-ab9e-4683-9ad0-bcc92ead9b31","added_by":"auto","created_at":"2025-03-18 09:10:33","extension":"png","order_by":7,"title":"Figure 7","display":"","copyAsset":false,"role":"figure","size":440135,"visible":true,"origin":"","legend":"\u003cp\u003eClimate and R\u0026amp;D effects on yield anomaly SD. (A) Relative counts of crops showing that climate change in 2000–2019 increased the yield anomaly SD, compared to what would occur under counterfactual climate. (B) Relative counts of crops showing that actual R\u0026amp;D had a significant effect in decreasing yield anomaly SD, compared to what would occur in the absence of investments to R\u0026amp;D after 2000. The differences in annual rate of changing in yield anomaly SD in A and B are respectively tested using the two-sided Wilcoxon rank sum test (p\u0026lt;0.05). (C) Relative counts of crops showing an increase in yield anomaly SD derived using the grid yield dataset (same as in Fig. 3 B). (D, E) The pie charts show the percentage of the global harvested area in 2010 (\u003cem\u003e74\u003c/em\u003e) with color-coded categories for changes in yield anomaly SD as well as the climate and R\u0026amp;D effects. (F) Area share for the R\u0026amp;D effect categories estimated from year to year. In F, levels of R\u0026amp;D expenditure, N application and irrigation were kept unchanged after the year indicated in the x-axis to predict yields with the random forest models. The results obtained for individual crops are available in figs. S7–S8.\u003c/p\u003e","description":"","filename":"7.png","url":"https://assets-eu.researchsquare.com/files/rs-6174385/v1/5b133e4b6f5514461b30353e.png"},{"id":79380303,"identity":"7db6ed91-ac1b-4774-8c06-da74ff54c801","added_by":"auto","created_at":"2025-03-18 09:18:33","extension":"png","order_by":8,"title":"Figure 8","display":"","copyAsset":false,"role":"figure","size":269704,"visible":true,"origin":"","legend":"\u003cp\u003eAreas shares of R\u0026amp;D and climate effects categories by crop and region. (A) The R\u0026amp;D effect is classified into two categories, i.e., R\u0026amp;D after 2000 has a significant, or non-significant effect, in increasing yield stability in 2000–2019 under climate change. (B) The climate effect is classified into three categories, i.e., climate change, relative to preindustrial levels, leads to a favorable, unfavorable, or no significant change in agro-climatic condition in terms of yield stability. The region codes are as follows: EAs, East Asia; EEuCAs, Eastern Europe and Central Asia; LAmCa, Latin America and Caribbean; MENAf, Middle East and North Africa; NAm, North America; SAs, South Asia; SEAsO, Southeast Asia and Oceania; SSAf, Sub-Saharan Africa; and WEu, Western Europe (see also fig. S11).\u003c/p\u003e","description":"","filename":"8.png","url":"https://assets-eu.researchsquare.com/files/rs-6174385/v1/aa13f5a8afc58fef8c87e91a.png"},{"id":79380282,"identity":"07c94174-537f-4f49-9eef-c5c6df7eb6b0","added_by":"auto","created_at":"2025-03-18 09:34:34","extension":"pdf","order_by":0,"title":"","display":"","copyAsset":false,"role":"manuscript-pdf","size":2681796,"visible":true,"origin":"","legend":"","description":"","filename":"manuscript.pdf","url":"https://assets-eu.researchsquare.com/files/rs-6174385/v1/58900288-5ac1-4e73-aa68-0e39109fc501.pdf"},{"id":79380284,"identity":"0f5d9d74-9e83-49d7-8861-a6f6a7e6b657","added_by":"auto","created_at":"2025-03-18 09:26:33","extension":"docx","order_by":1,"title":"","display":"","copyAsset":false,"role":"supplement","size":2459266,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eSupplementary material\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eFigs. S1 to S16\u003c/p\u003e","description":"","filename":"sm20250304.docx","url":"https://assets-eu.researchsquare.com/files/rs-6174385/v1/8ece5f6a8eb87694a6e54c68.docx"},{"id":79380290,"identity":"16b60e3d-fffc-4d2d-8a31-ad97cfe03152","added_by":"auto","created_at":"2025-03-18 09:26:32","extension":"docx","order_by":2,"title":"","display":"","copyAsset":false,"role":"supplement","size":42657,"visible":true,"origin":"","legend":"","description":"","filename":"table11.docx","url":"https://assets-eu.researchsquare.com/files/rs-6174385/v1/85067b84eba3043e84bd244f.docx"}],"financialInterests":"The authors declare no competing interests.","formattedTitle":"Assessing the capacity of agricultural research and development to increase the stability of global crop yields under climate change","fulltext":[{"header":"Significance Statement ","content":"\u003cp\u003eImproved technologies and management practices - fostered by agricultural research and development (R\u0026amp;D) - have helped to increase yield stability, but it is not known whether R\u0026amp;D investments have been sufficient to fully offset the increase in yield variability induced by climate change. Based on a global dataset, our analysis quantifies the relationships between R\u0026amp;D expenditure, climate change, and yield stability for major agricultural commodities worldwide.\u0026nbsp;Our analysis reveals that R\u0026amp;D investments have probably not been high enough since 2000 to fully mitigate the increased variability of yields resulting from climate change. Future research is needed to establish a causal relationship between R\u0026amp;D and yield stability underlying the empirical results found in this study.\u003c/p\u003e"},{"header":"Introduction","content":"\u003cp\u003eAmong the four pillars of food security (availability, access, utilization and stability) (\u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e), ensuring the production stability is crucial for preventing any deleterious consequences cascading to agri-food systems downstream. The disruptions to food supply chains caused by recent conflicts and pandemic reaffirm the importance of stability in the interconnected world (\u003cspan additionalcitationids=\"CR3\" citationid=\"CR2\" class=\"CitationRef\"\u003e2\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR4\" class=\"CitationRef\"\u003e4\u003c/span\u003e). In the past, production stability has been improved by optimizing growing areas, harvestable area fraction, annual number of harvests and cultivar tolerance to abiotic and biotic stresses (\u003cspan citationid=\"CR5\" class=\"CitationRef\"\u003e5\u003c/span\u003e, \u003cspan citationid=\"CR6\" class=\"CitationRef\"\u003e6\u003c/span\u003e). Increasing yield stability remains an important avenue to make crop production more climate-resilient, especially because yield growth has been the predominant contributor to production increase in many parts of the world in recent decades, while the increase in cultivated areas have been a much smaller contributor since 2000 (\u003cspan citationid=\"CR5\" class=\"CitationRef\"\u003e5\u003c/span\u003e, \u003cspan citationid=\"CR7\" class=\"CitationRef\"\u003e7\u003c/span\u003e).\u003c/p\u003e \u003cp\u003eFrom 1981 to 2010, yield stability for major crops (maize, rice, wheat and soybean) decreased in 9\u0026ndash;22% of the global harvested area in response to increased crop exposure to heat, yet 19\u0026ndash;33% of the area saw stabilizing yields (\u003cspan citationid=\"CR8\" class=\"CitationRef\"\u003e8\u003c/span\u003e). After the study by Iizumi and Ramankutty (\u003cspan citationid=\"CR8\" class=\"CitationRef\"\u003e8\u003c/span\u003e) published in 2016, many production shocks associated with extreme climate events have been reported (\u003cspan additionalcitationids=\"CR10 CR11\" citationid=\"CR9\" class=\"CitationRef\"\u003e9\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR12\" class=\"CitationRef\"\u003e12\u003c/span\u003e). Research and development (R\u0026amp;D) on climate-resilient agri-food systems is expanding, including breeding targeting emerging agro-climatic conditions (\u003cspan citationid=\"CR13\" class=\"CitationRef\"\u003e13\u003c/span\u003e) and climate-smart agriculture. However, although R\u0026amp;D innovations have been extensively tested at experimental stations, there is limited evidence showing that advances in R\u0026amp;D and producers\u0026rsquo; adoption of improved seeds, technologies and management practices have reduced climate risks at large spatial scales (\u003cspan citationid=\"CR14\" class=\"CitationRef\"\u003e14\u003c/span\u003e). Despite improvements in breeding and weather forecasting in the last decades, some studies suggest that the yield sensitivity to extreme climate events has not been substantially reduced in high-yielding regions (\u003cspan citationid=\"CR15\" class=\"CitationRef\"\u003e15\u003c/span\u003e, \u003cspan citationid=\"CR16\" class=\"CitationRef\"\u003e16\u003c/span\u003e). It is therefore crucial to update our understanding on yield stability and its drivers. In this paper, we contribute to this effort by studying how climate change and agricultural R\u0026amp;D have affected yield stability since 2000.\u003c/p\u003e \u003cp\u003eRecent advances in climate attribution science and resulting climate model simulations (\u003cspan additionalcitationids=\"CR18\" citationid=\"CR17\" class=\"CitationRef\"\u003e17\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR19\" class=\"CitationRef\"\u003e19\u003c/span\u003e) allow researchers to estimate the contributions of human-induced climate change to average yields (\u003cspan citationid=\"CR7\" class=\"CitationRef\"\u003e7\u003c/span\u003e, \u003cspan citationid=\"CR20\" class=\"CitationRef\"\u003e20\u003c/span\u003e) and crop failures (\u003cspan citationid=\"CR21\" class=\"CitationRef\"\u003e21\u003c/span\u003e, \u003cspan citationid=\"CR22\" class=\"CitationRef\"\u003e22\u003c/span\u003e). However, contrasting to the accumulated evidence showing that climate change has slowed down the growth of average yields (\u003cspan citationid=\"CR7\" class=\"CitationRef\"\u003e7\u003c/span\u003e), little is known about climate-induced changes in yield stability. In addition, it remains difficult to assess the capacity of R\u0026amp;D on reducing climate risks, although this information is often requested by research project funders (\u003cspan citationid=\"CR23\" class=\"CitationRef\"\u003e23\u003c/span\u003e). To fill these gaps, we present here a global analysis of the changes in yield stability due to climate change and R\u0026amp;D from 2000 to 2019, using a large ensemble climate model simulation (\u003cspan citationid=\"CR17\" class=\"CitationRef\"\u003e17\u003c/span\u003e) and machine learning models that relate interannual variation in yields of maize, rice, wheat and soybean at the 0.5\u0026deg; grid scale (~\u0026thinsp;55 km at the equator) to climate, soil, management and agricultural R\u0026amp;D expenditure. The four crops considered here are of high importance as they provide nearly two-thirds of the world agricultural calories (\u003cspan citationid=\"CR24\" class=\"CitationRef\"\u003e24\u003c/span\u003e).\u003c/p\u003e"},{"header":"Results","content":"\u003cdiv id=\"Sec3\" class=\"Section2\"\u003e \u003ch2\u003eR\u0026amp;D capacity to increase and stabilize crop yields\u003c/h2\u003e \u003cp\u003eYield anomaly\u0026mdash;the departure of annual yield from trend of average yield\u0026mdash;is a key metric measuring year-to-year yield fluctuation. The analysis of country-level yield anomalies over 2000\u0026ndash;2019 showed that increases in public R\u0026amp;D expenditure are proportional to increases in average yields. Considering the four crops taken together, countries with large areas under cultivation (\u0026ge;\u0026thinsp;1 Mha) had invested considerably more in R\u0026amp;D than countries with smaller cropping areas (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003eA). Average yield and yield anomaly SD are positively correlated (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003eB) but, for a given average yield level, countries with higher levels of R\u0026amp;D expenditure often achieved higher stability (lower SD and/or lower CV\u0026thinsp;=\u0026thinsp;SD/average) than countries with lower levels of R\u0026amp;D expenditure (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003eC).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003e[Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003e here]\u003c/p\u003e \u003c/div\u003e\n\u003ch3\u003eRecent patterns of agricultural R\u0026D expenditure\u003c/h3\u003e\n\u003cp\u003eThe level of cumulative agricultural R\u0026amp;D expenditure varied considerably between countries, although increasing trends are quite common (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003eA). Compared to annual R\u0026amp;D expenditure, their cumulative values, which take into account the time lag between research and adoption of new technologies by producers and technology obsolescence, are relevant to explain the increase in average yield and yield stability. Since 2000, annual R\u0026amp;D expenditure was the highest in China, followed by the United States of America (the US) (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003eB). China\u0026rsquo;s annual spending increased rapidly from 3.6\u0026nbsp;billion USD (\u003cspan\u003e$\u003c/span\u003eB) in 2000 to 22.9 \u003cspan\u003e$\u003c/span\u003eB in 2019, while the US annual spending increased from 4.2 \u003cspan\u003e$\u003c/span\u003eB to 5.2 \u003cspan\u003e$\u003c/span\u003eB for the period. The level of annual R\u0026amp;D expenditure in other countries were much less than these two countries. Although cumulative R\u0026amp;D expenditure increased in most countries, for some countries, cumulative R\u0026amp;D expenditure stagnated (e.g., in France) (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003eA, C). When the level of annual R\u0026amp;D expenditure was not sufficient to maintain the level of cumulative R\u0026amp;D expenditure, which is subject to annual rates of technology obsolescence, cumulative R\u0026amp;D expenditure declined with time (e.g., in Russia and Sudan).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003e[Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003e here]\u003c/p\u003e\n\u003ch3\u003eRecent changes in yield stability\u003c/h3\u003e\n\u003cp\u003eTo analyze the evolution of yield stability, SD was selected as the main metric because SD was more sensitive than CV in detecting climate-induced changes in yield anomalies. Under climate change, yield anomaly SD for a certain period may increase as a result of both occurrences of poor harvests (due to severe agro-climatic conditions) and good harvests (due to favorable agro-climatic conditions, or to elevated carbon dioxide concentrations, which enhance photosynthesis). However, with CV, such a change in SD can be masked by an increase in average yield because CV decreases when average yield increases faster than yield anomaly SD over time.\u003c/p\u003e \u003cp\u003eOur grid-wise analysis showed that 20% of the global harvested area (129.2 Mha, out of 646.2 Mha in 2010) had increased yield anomaly SD in 2000\u0026ndash;2019 (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003eB), while areas with increased CV were only 58% of the areas with increased SD (75.8 Mha; Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003eA). The overall increasing trends in average yields over the same period explained this difference in results between SD and CV (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003eC).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003e[Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003e here]\u003c/p\u003e\n\u003ch3\u003eModelling the effects of climate and R\u0026D on yield anomaly\u003c/h3\u003e\n\u003cp\u003eMachine learning random forest (RF) models were trained with current climate conditions (represented by the retrospective meteorological forcing data; \u003cb\u003eMaterials and methods\u003c/b\u003e) to attribute the changes in modeled yield anomaly SD to the R\u0026amp;D and climate-related explanatory variables (Table\u0026nbsp;1). The models explained 19% (rice) to 43% (soybean) of the variance of the grid yield anomaly between the locations, seasons and years when evaluated using out-of-bag data. Intermediate percentage of explained variance was obtained for maize and winter wheat (26% for both) and spring wheat (39%). Concerning the model ability to reproduce interannual yield variations, the leave-one-year-out cross-validation led to average correlations ranging from 0.639 (major-season rice) to 0.793 (soybean), and to average root mean squared errors ranging from 0.66 (spring wheat) to 0.81 (major-season rice) in SD units (figs. S4\u0026ndash;S5).\u003c/p\u003e \u003cp\u003eImportantly, the RF model projections obtained with d4PDF factual climate simulations (representation of the actual climate simulated by the climate model, including human-induced climate change) were able to reproduce the geographic pattern of the changes in yield anomaly SD for the crops, with moderate agreement with the actual pattern derived from the grid yield dataset [Cohen\u0026rsquo;s kappa (κ) of 0.55; Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003eA, C]. The agreement with the actual pattern was as high as κ\u0026thinsp;=\u0026thinsp;0.71 when the RF models were forced by the retrospective meteorological forcing data (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003eA, B). We also checked that d4PDF factual and counterfactual climate simulations fell within the ranges covered by the training dataset as the RF model performance often diminish outside the training data (only 0.18% of cases across all climate variables, locations, seasons, years and ensemble members were outside the ranges covered by the training dataset).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003e[Table\u0026nbsp;1 here]\u003c/p\u003e \u003cp\u003e[Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003e here]\u003c/p\u003e\n\u003ch3\u003eVariable importance\u003c/h3\u003e\n\u003cp\u003eThe RF models were used to rank the relative importance of R\u0026amp;D, climate, soil and management variables in explaining the yield anomaly. Two methods were considered to determine the importance of variables (\u003cb\u003eMaterials and methods\u003c/b\u003e). Cumulative R\u0026amp;D expenditure was identified as the primary factor, irrespective of the crops and methods (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003e). N application rate was often detected as the next important factor after R\u0026amp;D expenditure for the crops except soybean (fig. S6). Temperatures ranging from 15\u0026deg;C to 35\u0026deg;C, which encompass the exposures to temperature fluctuations around optimal thermal conditions for crops, were identified as important, albeit with some crop-specific variation. For maize, rice and soybean, the relatively higher temperatures from the 5th bin (15\u0026ndash;20\u0026deg;C) to the 9th bin (\u0026ge;\u0026thinsp;35\u0026deg;C) were more important. For winter and spring wheat, the relatively lower temperatures from the 4th bin (10\u0026ndash;15\u0026deg;C) to the 5th bin (15\u0026ndash;20\u0026deg;C) showed stronger effects. Among the precipitation variables, the 1st bin (\u0026lt;\u0026thinsp;1 mm d\u003csup\u003e\u0026ndash;1\u003c/sup\u003e) and the 5th bin (\u0026ge;\u0026thinsp;30 mm d\u003csup\u003e\u0026ndash;1\u003c/sup\u003e) tended to the most important, with variations between the crops and methods. On average, SOC was one of the most important soil variables. Irrigation intensity was less influential.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003e[Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003e here]\u003c/p\u003e \u003cdiv id=\"Sec8\" class=\"Section2\"\u003e \u003ch2\u003eClimate effect\u003c/h2\u003e \u003cp\u003eEffects of climate and R\u0026amp;D on the changes in yield stability at SD level were estimated by comparing the yield anomalies simulated by the RF models for different scenarios. Two yield anomaly time series were generated using d4PDF factual and counterfactual climate simulations, representing the historical climate (including human-induced climate change) and the preindustrial climate (without human-induced climate change), respectively (\u003cem\u003efc\u003c/em\u003e and \u003cem\u003ect\u003c/em\u003e runs; Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003e). R\u0026amp;D expenditures were fixed to those observed in both scenarios. We classified the change in agro-climatic condition due to climate change as \u0026ldquo;unfavorable\u0026rdquo; (\u0026ldquo;favorable\u0026rdquo;) when the simulated yield anomaly SD was higher (lower) with the factual climate than with the counterfactual climate.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eThe global grid-wise attribution revealed that unfavorable agro-climatic conditions were relatively more prevalent in the cropland area where an increase in yield anomaly SD occurred than in the area with decreasing yield anomaly SD. The area with unfavorable agro-climatic conditions was estimated to account for the half (0.5) of the area with an increase in yield anomaly SD; this proportion (0.5) was derived by comparing the area with both unfavorable agro-climatic conditions and increasing SD (10% of the global harvested area; Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003eA, D) and the area with increasing yield anomaly SD all agro-climatic conditions taken together (20% of the global harvested area; Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003eC, D). This proportion (0.5) is twice as high as the proportion (0.25) of area with unfavorable agro-climatic conditions and decreasing SD (4%) in the area with decreasing yield anomaly SD all agro-climatic conditions taken together (16%) (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003eD). In other words, the frequency of unfavorable agro-climatic conditions is twice as high in areas where the yield anomaly SD is increasing than in areas where the yield anomaly SD is decreasing.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eConsistent with the above result, favorable agro-climatic conditions were relatively more common in the area where a decrease in yield anomaly SD occurred than in the area with increasing yield anomaly SD. Favorable agro-climatic conditions were found in a quarter (0.25) of the area with decreasing yield anomaly SD; this proportion (0.25) compared the area with both favorable agro-climatic conditions and decreasing SD (4%) and the area with decreasing yield anomaly SD all agro-climatic conditions taken together (16%) (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003eD). This proportion (0.25) is 1.7 times higher than the proportion (0.15) between the area with both favorable agro-climatic conditions and increasing yield anomaly SD (3%) and the area with increasing yield anomaly SD all agro-climatic conditions taken together (20%).\u003c/p\u003e \u003c/div\u003e\n\u003ch3\u003eR\u0026D effect\u003c/h3\u003e\n\u003cp\u003eThe R\u0026amp;D effect was estimated by comparing two series of yield anomaly simulations, i.e., a series simulated with the actual R\u0026amp;D expenditures in 2000\u0026ndash;2019 and another series simulated without R\u0026amp;D investment (\u003cem\u003efc\u003c/em\u003e and \u003cem\u003efc.em\u003c/em\u003e runs; Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003e). The results showed that a significant yield-stabilizing effect of R\u0026amp;D was found in a large part (0.44) of the cropland area where a decrease in yield anomaly SD occurred. This proportion (0.44) was derived by comparing the area with a significant R\u0026amp;D effect of 7% and the area with decreasing yield anomaly SD of 16%, both relative to the global harvested area (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003eB, C, E). This proportion (0.44) is 4.4 times greater than the proportion (0.10) calculated between the area with a significant R\u0026amp;D effect (2%) and the area with increasing yield anomaly SD (20%).\u003c/p\u003e \u003cp\u003eOn the other hand, the yield-stabilizing effect of R\u0026amp;D was not found to be significant in a dominant proportion (0.85) of the area where an increase in yield anomaly SD occurred. This value (0.85) was derived from the area with a non-significant R\u0026amp;D effect of 17% and the area with increasing yield anomaly SD of 20% (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003eE), which is 1.7 times greater than the value (0.50) calculated from the area with a non-significant R\u0026amp;D effect (8%) and the area with decreasing yield anomaly SD (16%). Importantly, the area with a significant R\u0026amp;D effect became smaller when estimated over a more recent, shorter period than over the entire study period (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003eF).\u003c/p\u003e \u003cp\u003e[Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003e here]\u003c/p\u003e \u003cp\u003e[Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003e here]\u003c/p\u003e\n\u003ch3\u003eSensitivity to different R\u0026D assumptions\u003c/h3\u003e\n\u003cp\u003eThe estimated area shares for the climate and R\u0026amp;D effect categories were found to be robust to the different assumptions used in the R\u0026amp;D variable calculation. We took into account whether the cumulative R\u0026amp;D expenditure was for agriculture or for the four crops, whether the time lag between research and adoption of technologies by producers was relatively shorter (6 years) or longer (12 years), and whether the technology obsolescence was relatively slower (10% per year) or faster (20% per year) (\u003cb\u003eMaterials and methods\u003c/b\u003e). The main results described in the former subsections (specifically, the fact that R\u0026amp;D is the primary explanatory variable explaining yield anomaly, and unfavorable agro-climatic conditions and/or non-significant R\u0026amp;D effect are proportionately more prevalent in the cropland area experiencing a decrease in yield stability (higher SD)) were robust to the different R\u0026amp;D assumptions (figs. S9\u0026ndash;S10).\u003c/p\u003e \u003cdiv id=\"Sec11\" class=\"Section2\"\u003e \u003ch2\u003eRegion and crop-specific characteristics\u003c/h2\u003e \u003cp\u003eWe identified the regions and crop species for which R\u0026amp;D were able to compensate for the negative impacts of climate change on yield stability after 2000. For the four crops taken together (479.8 Mha), R\u0026amp;D had no significant effect in 75% of this area. However, R\u0026amp;D was found to have a significant stabilizing effect on large shares of croplands in several specific crop species and regions, especially for soybean [42% (significant R\u0026amp;D effect) vs. 58% (non-significant)], in Sub-Saharan Africa (41% vs. 59%), East Asia (35% vs. 65%) and North America (29% vs. 71%) (Fig.\u0026nbsp;\u003cspan refid=\"Fig8\" class=\"InternalRef\"\u003e8\u003c/span\u003eA).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eClimate change had no significant effect on yield stability in almost half of this area (43%), but led to unfavorable agro-climatic conditions (higher SD) and favorable agro-climatic conditions (lower SD) in 36% and 21% of the area, respectively. The negative effect of climate change on yield stability was particularly pronounced for maize [45% (unfavorable) vs. 22% (favorable)] and soybean (57% vs. 23%), and in North America (81% vs. 3%), followed by East Asia (33% vs. 18%). On the contrary, the shift to favorable agro-climatic conditions was relatively prominent in South Asia (10% vs. 36%) (Fig.\u0026nbsp;\u003cspan refid=\"Fig8\" class=\"InternalRef\"\u003e8\u003c/span\u003eB).\u003c/p\u003e \u003cp\u003e[Fig.\u0026nbsp;\u003cspan refid=\"Fig8\" class=\"InternalRef\"\u003e8\u003c/span\u003e here]\u003c/p\u003e \u003c/div\u003e"},{"header":"Discussion","content":"\u003cp\u003eYield stability has increased in 16% of the global cropland. Although not mutually exclusive, nearly half of this area is associated with a significant stabilizing effect of R\u0026amp;D, and a quarter of this area is associated with a shift to favorable agro-climatic conditions. In contrast, yield stability has decreased in 20% of the global cropland. Half of this area is associated with a shift to unfavorable agro-climatic conditions, while about nine-tenths of this area is associated with a lack of significant R\u0026amp;D effect, although these areas may overlap. These results provide insights into the association between investment in agricultural R\u0026amp;D, climate change, and resulting changes in yield stability. In addition, based on these results, we could formulate the hypothesis that, in the global cropland area where yield stability has decreased, the progress made in R\u0026amp;D since 2000 has lagged behind the negative impacts of climate change. This question needs to be addressed in future research.\u003c/p\u003e \u003cp\u003eAlthough R\u0026amp;D effect is relatively more prevalent in the area with increasing yield stability, significant R\u0026amp;D effect is also detected in some part of the global cropland with unchanged and decreasing yield stability (8% and 3%, respectively). In this study, the R\u0026amp;D effect is considered significant when R\u0026amp;D expenditures have contributed to a decrease of yield anomaly SD over the study period. Even if the R\u0026amp;D effect is found to be significant, decreasing or unchanged yield stability may occur if the R\u0026amp;D effect is not strong enough to reverse the trend in yield anomaly SD. However, without any R\u0026amp;D effect, the yield stability status in these areas would probably deteriorate further.\u003c/p\u003e \u003cp\u003eIt should be noted that the R\u0026amp;D effect on yield stability identified here is empirical rather than causal, as some confounding factors potentially correlated with R\u0026amp;D expenditure were not considered in this study. Process-resolving analyses would be useful to determine the causal effects of individual technologies and practices, and confirm our empirical evidence of the yield-stabilizing effect of R\u0026amp;D, which could ultimately support the implementation of effective interventions. Although applied to study average yield and not yield stability, a relevant example is presented in (\u003cspan citationid=\"CR25\" class=\"CitationRef\"\u003e25\u003c/span\u003e), which decomposes the historical yield increases into individual drivers, including increases in the use of land, water, and other inputs, and increases in the efficiency of input use, with an analysis of policy implications.\u003c/p\u003e \u003cdiv id=\"Sec13\" class=\"Section2\"\u003e \u003ch2\u003eImportant factors of yield anomaly\u003c/h2\u003e \u003cp\u003eThe yield-stabilizing effect of R\u0026amp;D is not totally surprising. For instance, new varieties bred through R\u0026amp;D generally have greater tolerance to sub-optimal growing conditions and higher yield potential than earlier varieties (\u003cspan additionalcitationids=\"CR27\" citationid=\"CR26\" class=\"CitationRef\"\u003e26\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR28\" class=\"CitationRef\"\u003e28\u003c/span\u003e). In addition to breeding, R\u0026amp;D informs producers about practices to reduce climate risks, such as sowing date, soil, water and fertilizer management, and pests, diseases and weed control.\u003c/p\u003e \u003cp\u003eThe temperatures in the growing season, as well as occurrences of dry and heavy rainfall days are identified to be important climatic factors explaining yield anomaly. This result is consistent with previous findings (\u003cspan additionalcitationids=\"CR30\" citationid=\"CR29\" class=\"CitationRef\"\u003e29\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR31\" class=\"CitationRef\"\u003e31\u003c/span\u003e). Under warming, a decrease in the frequency of near-optimal temperatures is accompanied by an increase in the frequency of supra-optimal temperatures, which decreases yield stability (\u003cspan citationid=\"CR32\" class=\"CitationRef\"\u003e32\u003c/span\u003e). The effects on yield stability from heavily wet conditions identified here are more uncertain than that the effects of dry conditions, although there is a growing body of literature reporting the negative impacts of waterlogging and resulting excessive soil water (\u003cspan citationid=\"CR33\" class=\"CitationRef\"\u003e33\u003c/span\u003e).\u003c/p\u003e \u003cp\u003eThe importance of N application in explaining yield anomaly can be explained by previous findings that abiotic stresses due to growing season climate become increasingly important stressors as N inputs increase (\u003cspan citationid=\"CR34\" class=\"CitationRef\"\u003e34\u003c/span\u003e, \u003cspan citationid=\"CR35\" class=\"CitationRef\"\u003e35\u003c/span\u003e). The importance of SOC found in this study confirms previous findings that drought damage is mitigated by SOC build-up (\u003cspan citationid=\"CR36\" class=\"CitationRef\"\u003e36\u003c/span\u003e). Despite the initial expectation that irrigation would increase yield stability, the contribution of irrigation appears to be limited in this study. The fact that the existence of irrigation equipment does not necessarily ensure water availability during droughts (\u003cspan citationid=\"CR37\" class=\"CitationRef\"\u003e37\u003c/span\u003e), and that yield sensitivity at regional levels to heat no longer decreases as irrigated area saturates (\u003cspan citationid=\"CR38\" class=\"CitationRef\"\u003e38\u003c/span\u003e), partly explain the counterintuitive outcome.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec14\" class=\"Section2\"\u003e \u003ch2\u003eLimitations and uncertainties\u003c/h2\u003e \u003cp\u003eThere are a few limitations to this study. First, we only consider public expenditures and do not include expenditures by private sector and international organizations. Therefore, the R\u0026amp;D effect estimated in this study may be underestimated to some extent, given that the private sector's share of total agricultural R\u0026amp;D expenditure can be substantial, particularly in high-income countries (30%) compared with low- and middle-income countries (5%) (\u003cspan citationid=\"CR39\" class=\"CitationRef\"\u003e39\u003c/span\u003e). Moreover, the role of international organizations for agricultural R\u0026amp;D is crucial in low-income countries, and investments in agricultural R\u0026amp;D for one country may generate spillover benefits for neighboring countries and regions with similar agro-climatic conditions. Data on private R\u0026amp;D expenditure are difficult to obtain. Even for public expenditure, the information on the technology on which R\u0026amp;D expenditure was spent is rarely available (\u003cspan citationid=\"CR40\" class=\"CitationRef\"\u003e40\u003c/span\u003e). Future research should fill these data gaps to distinguish the relative contributions of individual technologies and management practices.\u003c/p\u003e \u003cp\u003eSecond, a direct comparison with previous estimates of climate-induced changes in yield stability is not possible due to a lack of relevant literature. There is literature showing changes in average yield (\u003cspan citationid=\"CR20\" class=\"CitationRef\"\u003e20\u003c/span\u003e, \u003cspan citationid=\"CR23\" class=\"CitationRef\"\u003e23\u003c/span\u003e). However, none of these studies assessed changes in yield stability and the R\u0026amp;D effect.\u003c/p\u003e \u003cp\u003eThe remaining limitations and uncertainties are mostly related to the underlying data. Third, this study is based on a single grid yield dataset. Different yield datasets often produce different results. Subnational census statistics are expected to be more reliable than satellite- or model-based yield datasets, including one adopted here. However, the spatial resolution of census statistics is coarser than the resolution of d4PDF climate simulation, which compromises the ability to capture the yield impacts from extreme agro-climatic conditions. In addition, some recent grid yield datasets (\u003cspan citationid=\"CR41\" class=\"CitationRef\"\u003e41\u003c/span\u003e) use temperature and precipitation as inputs in their yield estimation procedure, positing the circularity issues because we relate yield anomaly to temperature and precipitation. Fourth, attribution outcomes may be influenced by the choice of meteorological forcing datasets and climate models (\u003cspan citationid=\"CR42\" class=\"CitationRef\"\u003e42\u003c/span\u003e, \u003cspan citationid=\"CR43\" class=\"CitationRef\"\u003e43\u003c/span\u003e). Fifth, our analyses are based on the harvested area in 2010, although the area for the four crops in 2022 increased 12% since 2010. This increase in area would include not only a response to increasing demand but also a shift of areas suitable for cultivation in response to the observed warming (\u003cspan citationid=\"CR44\" class=\"CitationRef\"\u003e44\u003c/span\u003e, \u003cspan citationid=\"CR45\" class=\"CitationRef\"\u003e45\u003c/span\u003e). However, we used the grid harvested area map for 2010, the midpoint of the study period, as no annual crop-specific grid harvested area map is available at the global scale. Last, future research should consider more management variables than those considered here, such as pesticide application rates, which are thought to increase yield stability (\u003cspan citationid=\"CR46\" class=\"CitationRef\"\u003e46\u003c/span\u003e).\u003c/p\u003e \u003c/div\u003e"},{"header":"Materials and methods","content":"\u003cdiv id=\"Sec16\" class=\"Section2\"\u003e \u003ch2\u003eExperimental design\u003c/h2\u003e \u003cp\u003eThis section describes the data on crop yields, climate (actual, factual and counterfactual) and agricultural R\u0026amp;D expenditure, as well as the method used to detect changes in yield stability. Also included are the machine learning models and the yield impact attribution analyses.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec17\" class=\"Section2\"\u003e \u003ch2\u003eCrop yields\u003c/h2\u003e \u003cp\u003eWe obtained country annual yield data from the FAO statistical database for maize, rice, wheat and soybean (\u003cspan citationid=\"CR47\" class=\"CitationRef\"\u003e47\u003c/span\u003e). In addition, we used the 0.5\u0026deg; grid global dataset of historical yields (GDHY) (\u003cspan citationid=\"CR48\" class=\"CitationRef\"\u003e48\u003c/span\u003e). In the GHDY dataset, yield data are available for two seasons for maize and rice (major and secondary) and for wheat (winter and spring). For soybean, data are available for a single main season. The GDHY dataset is a hybrid of satellite vegetation index and national census statistics reported by FAO. The information on accumulated net primary production over the crop-specific growing period, calculated using remotely-sensed leaf area index and fraction of photosynthetically-active radiation as well as reanalysis downward shortwave radiation and crop-specific radiation-use efficiency, is combinedly used to spatially disaggregate FAO country average yields to the grid level. Grid yield values are available in 76\u0026ndash;92% of the global harvested area with variations by crop and lacking when crop calendar information (\u003cspan citationid=\"CR49\" class=\"CitationRef\"\u003e49\u003c/span\u003e) is missing.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec18\" class=\"Section2\"\u003e \u003ch2\u003eClimate\u003c/h2\u003e \u003cp\u003eWe used 0.5\u0026deg; grid daily data on maximum and minimum air temperatures and precipitation, covering three distinct climate conditions in 2000\u0026ndash;2019: \u003cem\u003eactual\u003c/em\u003e, \u003cem\u003efactual\u003c/em\u003e and \u003cem\u003ecounterfactual\u003c/em\u003e. The \u003cem\u003eactual\u003c/em\u003e climate refers to observed conditions. The \u003cem\u003efactual\u003c/em\u003e climate represents a modelled approximation of the actual climate, influenced by both human activities and natural forces (e.g., volcanic eruptions and solar activity). In contrast, the \u003cem\u003ecounterfactual\u003c/em\u003e climate is a modeled preindustrial climate without any appreciable human impacts on the global climate.\u003c/p\u003e \u003cp\u003eThe actual climate data were obtained from the global retrospective meteorological forcing data. We used the JRA55-CDFDM-S14FD forcing data (called JCS for simplicity), which is a bias-corrected Japanese 55-year Reanalysis (JRA55) using the cumulative distribution function-based downscaling and bias-correction method (CDFDM) and the reference climatology from S14FD forcing dataset (\u003cspan citationid=\"CR50\" class=\"CitationRef\"\u003e50\u003c/span\u003e).\u003c/p\u003e \u003cp\u003eFor the factual and counterfactual climates, we used a large ensemble of long-term climate simulations with and without historical trends in external forcing to the global climate, called d4PDF (the Database for Policy Decision Making for Future Climate Change) (\u003cspan citationid=\"CR17\" class=\"CitationRef\"\u003e17\u003c/span\u003e, \u003cspan citationid=\"CR51\" class=\"CitationRef\"\u003e51\u003c/span\u003e). The climate simulations were performed using the Meteorological Research Institute Atmospheric General Circulation Model, version 3.2 (MRI-AGCM3.2) with a grid interval of 60 km (\u003cspan citationid=\"CR52\" class=\"CitationRef\"\u003e52\u003c/span\u003e). Each of the factual and counterfactual climate simulations has 100 ensemble members associated with slightly different initial conditions and small perturbations in sea surface temperatures that represent observational uncertainties. Comparisons between the factual and counterfactual climate simulations with large ensemble members enable researchers to conduct a robust attribution analysis of the effects of human activities on average climate conditions, extreme climate events (\u003cspan additionalcitationids=\"CR54 CR55\" citationid=\"CR53\" class=\"CitationRef\"\u003e53\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR56\" class=\"CitationRef\"\u003e56\u003c/span\u003e) and crop yields (\u003cspan citationid=\"CR57\" class=\"CitationRef\"\u003e57\u003c/span\u003e, \u003cspan citationid=\"CR58\" class=\"CitationRef\"\u003e58\u003c/span\u003e). The daily outputs from the climate model were interpolated to the 0.5\u0026deg; resolution and bias-corrected using the CDFDM method along with the reference climatology from S14FD forcing dataset.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec19\" class=\"Section2\"\u003e \u003ch2\u003eAgricultural R\u0026amp;D expenditure\u003c/h2\u003e \u003cp\u003ePublic agricultural R\u0026amp;D includes government, higher education and nonprofit entities, but excludes the private for-profit sector (\u003cspan citationid=\"CR39\" class=\"CitationRef\"\u003e39\u003c/span\u003e). It includes salaries, operating and program costs and capital investments for all entities, excluding the private for-profit sector, involved in agricultural R\u0026amp;D expenditure in a country.\u003c/p\u003e \u003cp\u003eFor the Organisation for Economic Co-operation and Development (OECD) member economies, we sourced annual data on gross domestic expenditure on R\u0026amp;D expressed as a percentage of gross domestic production (GDP) from the Main Science and Technology Indicators (MSTI) database (\u003cspan citationid=\"CR59\" class=\"CitationRef\"\u003e59\u003c/span\u003e). Agricultural R\u0026amp;D expenditure for OECD member economies was calculated by multiplying the GDP share of domestic R\u0026amp;D expenditure with the agricultural share of GDP and the total GDP, both obtained from World Bank (\u003cspan citationid=\"CR60\" class=\"CitationRef\"\u003e60\u003c/span\u003e). For 88 low and middle-income countries, data were obtained from the International Food Policy Research Institute (IFPRI)\u0026rsquo;s Agricultural Science and Technology Indicators (ASTI) database (\u003cspan citationid=\"CR39\" class=\"CitationRef\"\u003e39\u003c/span\u003e). The country annual R\u0026amp;D expenditure data collected were adjusted for inflation to be expressed in 2015 USD.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec20\" class=\"Section2\"\u003e \u003ch2\u003eDetection of changes in yield stability\u003c/h2\u003e \u003cp\u003eChanges in yield stability were detected based on the method of Iizumi and Ramankutty (\u003cspan citationid=\"CR8\" class=\"CitationRef\"\u003e8\u003c/span\u003e). For each grid cell, an annual yield time series (1981\u0026ndash;2019) was detrended with a cubic smoothing spline representing time-evolving average yield (fig. S12 A). Although the present study analyzed changes in yield stability between 2000 and 2019, a longer yield time series was used for a more robust identification of yield trends. Yield anomalies (in tonnes per hectare) were computed as the differences between the original yield time series and the fitted trends (fig. S12 B). By applying a 9-year centered moving window, an annual time series of the standard deviation (SD) of yield anomalies was derived. A regression line was fitted to the recent portion (2004\u0026ndash;2015) of the yield anomaly SD series, and the estimated slope was tested for a monotonic trend using the two-sided Mann-Kendall trend test (\u003cspan citationid=\"CR61\" class=\"CitationRef\"\u003e61\u003c/span\u003e) (fig. S12 C). Note that the yield stability trends analyzed here were based on the yield anomalies over the 2000\u0026ndash;2019 period (fig. S12 C). A significant positive (negative) slope indicated a decreasing (increasing) yield stability. Monotonic trends in average yield and the coefficient of variation (CV) time series of yield anomaly in the same period were tested in the same manner as the SD time series.\u003c/p\u003e \u003cp\u003eTo aid visualizing multi-season and multi-crop results as a single map, the detected changes in yield stability for individual seasons and crops were combined according to the following rules. The major season result was selected as the representative when both major (winter) and secondary (spring) seasons were operated. When only one season was operated, the result for that season was considered even when it was the secondary (spring) season. Then the number of crops with specific yield change (either of a significant increase in yield anomaly CV, significant increase in yield anomaly SD and significant increase in average yield), out of four (maize, wheat, rice and soybean), was counted and divided by the number of crops grown for a given location (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003eD).\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec21\" class=\"Section2\"\u003e \u003ch2\u003eRandom forest yield anomaly models\u003c/h2\u003e \u003cp\u003eMachine learning models were developed to relate yield anomaly to climate, soil, management and R\u0026amp;D expenditure. We used random forest regression models (\u003cspan citationid=\"CR62\" class=\"CitationRef\"\u003e62\u003c/span\u003e) and did not explore alternative machine learning approaches as random forests are straightforward to implement, and comparable to other approaches in many applications (\u003cspan citationid=\"CR63\" class=\"CitationRef\"\u003e63\u003c/span\u003e, \u003cspan citationid=\"CR64\" class=\"CitationRef\"\u003e64\u003c/span\u003e). Random forests are generally insensitive to correlations between explanatory variables in terms of predictive accuracy and overfitting (\u003cspan citationid=\"CR63\" class=\"CitationRef\"\u003e63\u003c/span\u003e), which further supports their use in this study. A random forest algorithm deals with collinearity by randomly sampling a small subset of input variables to build each splitting node of each tree. As a result, the learning algorithm can only choose among a small number of inputs at each splitting node, thus reducing the risk of collinearity.\u003c/p\u003e \u003cp\u003eWe calculated the z-score values of yield anomaly by dividing annual yield anomaly by the SD of yield anomaly for the same period as that used to estimate yield trends (fig. S12 D). The z-scored yield anomalies were related to climate, soil, management and R\u0026amp;D variables (Table\u0026nbsp;1). Details of the explanatory variables are provided in the subsequent sections.\u003c/p\u003e \u003cp\u003eA specific model was constructed for each crop. Data from major and secondary maize (rice) seasons were combined and used to develop a single maize (rice) model. However, the winter and spring wheat models were developed separately due to potential differences in the yield responses to environmental conditions because winter and spring wheat has different growing seasons and physiological characteristics. Model fitting was conducted using the statistical software R (\u003cspan citationid=\"CR65\" class=\"CitationRef\"\u003e65\u003c/span\u003e), using the randomForest package (\u003cspan citationid=\"CR66\" class=\"CitationRef\"\u003e66\u003c/span\u003e) with the following settings: number of trees (ntree\u0026thinsp;=\u0026thinsp;500), number of predictors sampled at each split (mtry\u0026thinsp;=\u0026thinsp;3), and minimum size of terminal nodes (nodesize\u0026thinsp;=\u0026thinsp;5). These hyperparameter settings were selected based on the package\u0026rsquo;s default values.\u003c/p\u003e \u003cp\u003eFor model validation purposes, we used the leave-one-year-out cross validation technique and assessed the Pearson\u0026rsquo;s correlation coefficients and the root mean squared errors calculated between the actual and simulated yield anomalies. Furthermore, the model evaluation was done to see whether the RF models with d4PDF factual climate simulation could reproduce the actual global pattern of the changes in yield anomaly SD derived from the grid yield dataset. We also checked whether values of d4PDF factual and counterfactual climate simulations did not exceed the ranges covered by the training data (i.e., the JCS retrospective meteorological forcing) to ensure the reliability of the RF models. Last, we assessed the importance of the explanatory variables. We used two methods to rank model input variables according to their importance; the permutation method and the node purity method that rely on the residual sum of squares and the Gini index as importance metrics, respectively.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec22\" class=\"Section2\"\u003e \u003ch2\u003eClimate variables\u003c/h2\u003e \u003cp\u003eTemperature and precipitation were considered to be the primary climatic drivers of yield anomalies. To account the effects of short-term extreme temperatures, relative counts of hourly temperature exposure for a growing season were calculated using 9 bins with constant 5\u0026deg;C intervals (T\u0026thinsp;\u0026lt;\u0026thinsp;0, 0\u0026thinsp;\u0026le;\u0026thinsp;T\u0026thinsp;\u0026lt;\u0026thinsp;5, 5\u0026thinsp;\u0026le;\u0026thinsp;T\u0026thinsp;\u0026lt;\u0026thinsp;10, \u0026hellip;, 25\u0026thinsp;\u0026le;\u0026thinsp;T\u0026thinsp;\u0026lt;\u0026thinsp;30, 30\u0026thinsp;\u0026le;\u0026thinsp;T\u0026thinsp;\u0026lt;\u0026thinsp;35, and 35\u0026thinsp;\u0026le;\u0026thinsp;T). Hourly temperature values were obtained by fitting a sine curve to daily maximum and minimum temperatures. Relative counts of daily precipitation exposure during the growing season were also calculated using 5 bins with unequal intervals (P\u0026thinsp;\u0026lt;\u0026thinsp;1, 1\u0026thinsp;\u0026le;\u0026thinsp;P\u0026thinsp;\u0026lt;\u0026thinsp;5, 5\u0026thinsp;\u0026le;\u0026thinsp;P\u0026thinsp;\u0026lt;\u0026thinsp;20, 20\u0026thinsp;\u0026le;\u0026thinsp;P\u0026thinsp;\u0026lt;\u0026thinsp;30, and 30\u0026thinsp;\u0026le;\u0026thinsp;P). The growth periods were fixed based on the reported planting and harvesting dates (\u003cspan citationid=\"CR49\" class=\"CitationRef\"\u003e49\u003c/span\u003e). For winter wheat only, the period before the completion of vernalization was excluded from the calculation by identifying the first date at which the fraction of growing-degree days to the crop total thermal requirements exceeds 0.1; this was estimated using the crop phenology model with the crop-specific base (0\u0026deg;C) and maximum temperature (26\u0026deg;C) for winter wheat. These climate variables were calculated for each of the actual, factual and counterfactual climate conditions.\u003c/p\u003e \u003cdiv id=\"Sec23\" class=\"Section3\"\u003e \u003ch2\u003eSoil variable\u003c/h2\u003e \u003cp\u003eWe focused on soil organic carbon (SOC) content as a main soil variable, given the fact that in arid and semiarid regions of the world, high levels of SOC improve soil water holding capacity and moderate drought damage (\u003cspan citationid=\"CR36\" class=\"CitationRef\"\u003e36\u003c/span\u003e). The topsoil (0\u0026ndash;30 cm) SOC data were obtained from the Regridded Harmonized World Soil Database v1.2 (\u003cspan citationid=\"CR67\" class=\"CitationRef\"\u003e67\u003c/span\u003e) and aggregated to the 0.5\u0026deg; resolution. In addition to SOC, we also used fractions of clay, silt and sand, bulk density, cation exchange capacity of the clay fraction, pH, and available water storage capacity, all at the topsoil, to consider possible effects of soil characteristics on yield anomaly.\u003c/p\u003e \u003c/div\u003e \u003c/div\u003e \u003cdiv id=\"Sec24\" class=\"Section2\"\u003e \u003ch2\u003eManagement variables\u003c/h2\u003e \u003cp\u003eWe considered two variables related to producer\u0026rsquo;s crop management, i.e., irrigation and N fertilizer. Irrigation mitigates drought damages and heat stress through evaporative cooling of the crop canopy (\u003cspan citationid=\"CR68\" class=\"CitationRef\"\u003e68\u003c/span\u003e). The irrigation intensity used in this study considered changes due to the expansion of the irrigation-equipped area based on the global historical irrigation dataset for the period of 1900\u0026ndash;2005 (HID) (\u003cspan citationid=\"CR69\" class=\"CitationRef\"\u003e69\u003c/span\u003e) and crop-specific variations based on the global monthly irrigated and rainfed crop areas around the year 2000 (MIRCA2000) (\u003cspan citationid=\"CR70\" class=\"CitationRef\"\u003e70\u003c/span\u003e). The data were aggregated to the 0.5\u0026deg; resolution. Due to the lack of data, irrigation intensity values for the period 2006\u0026ndash;2019 were extrapolated using a linear regression using the data from 1998 to 2005, although this assumption may not be totally accurate for some regions experiencing decreases in irrigation area, such as Russia (\u003cspan citationid=\"CR37\" class=\"CitationRef\"\u003e37\u003c/span\u003e). In addition to irrigation, we used the grid annual N application rate for each crop (\u003cspan citationid=\"CR71\" class=\"CitationRef\"\u003e71\u003c/span\u003e) as it is known that abiotic stresses become important as N inputs increase (\u003cspan citationid=\"CR34\" class=\"CitationRef\"\u003e34\u003c/span\u003e, \u003cspan citationid=\"CR35\" class=\"CitationRef\"\u003e35\u003c/span\u003e).\u003c/p\u003e \u003cdiv id=\"Sec25\" class=\"Section3\"\u003e \u003ch2\u003eR\u0026amp;D variable\u003c/h2\u003e \u003cp\u003eWe incorporated the effect of agriculture R\u0026amp;D and resulting technologies and management practices that help stabilize yields into the RF models. To do so, an economic indicator was used. The indicator, \u003cem\u003eR\u003c/em\u003e\u003csub\u003e\u003cem\u003ei\u003c/em\u003e,\u003cem\u003et\u003c/em\u003e\u003c/sub\u003e, is the accumulated agricultural R\u0026amp;D expenditure for country \u003cem\u003ei\u003c/em\u003e from 1995 to year \u003cem\u003et\u003c/em\u003e and expressed in constant 2015 USD. The accumulation took into account the duration of research and the obsolescence of technology:\u003cdiv id=\"Equ1\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equ1\" name=\"EquationSource\"\u003e\n$$\\:{R}_{i,\\:t}={E}_{i,t-6}+\\left(1-\\delta\\:\\right){R}_{i,t-1}$$\u003c/div\u003e\u003cdiv class=\"EquationNumber\"\u003e1\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003ewhere \u003cem\u003eE\u003c/em\u003e is the country annual agricultural R\u0026amp;D expenditure (constant 2015 USD); and \u003cem\u003eδ\u003c/em\u003e is the annual rate of technology obsolescence (=\u0026thinsp;0.1 or 10% per year) (\u003cspan citationid=\"CR71\" class=\"CitationRef\"\u003e71\u003c/span\u003e). The lag time between research and technology adoption by producers was set at six years, although lag time may vary by country and over time (\u003cspan citationid=\"CR72\" class=\"CitationRef\"\u003e72\u003c/span\u003e). As the R\u0026amp;D variable is of key importance in this study, we performed the sensitivity analysis on this variable, as described in the latter subsection.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec26\" class=\"Section3\"\u003e \u003ch2\u003eAttribution analyses\u003c/h2\u003e \u003cp\u003eSimulation experiments were conducted using the RF models and d4PDF factual and counterfactual climate simulations. As shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003e, the \u003cem\u003efc\u003c/em\u003e and \u003cem\u003ect\u003c/em\u003e runs compared the factual and counterfactual climates and determine the climate effect to yield stability under actual R\u0026amp;D conditions. The \u003cem\u003efc\u003c/em\u003e and \u003cem\u003efc.em\u003c/em\u003e runs compared yield stability with and without R\u0026amp;D after 2000 under climate change.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec27\" class=\"Section3\"\u003e \u003ch2\u003eClimate effect\u003c/h2\u003e \u003cp\u003eThe changes in yield stability associated with human-induced climate change was detected through a comparison of \u003cem\u003efc\u003c/em\u003e and \u003cem\u003ect\u003c/em\u003e runs (Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003e). Using 100 samples of the slope value of yield anomaly SD derived from the RF models forced by d4PDF climate simulation, a histogram was derived for each of the factual and counterfactual yields (fig. S13. When the average slope value of the factual yield anomaly SD is higher (lower) than that of the counterfactual yield anomaly SD, it indicates that climate change contributed to decreasing (increasing) yield stability in 2000\u0026ndash;2019. The two-sided Wilcoxon rank sum test (\u003cspan citationid=\"CR73\" class=\"CitationRef\"\u003e73\u003c/span\u003e) was used to test the null hypothesis that the two types of slopes were equal.\u003c/p\u003e \u003c/div\u003e \u003c/div\u003e \u003cdiv id=\"Sec28\" class=\"Section2\"\u003e \u003ch2\u003eR\u0026amp;D effect\u003c/h2\u003e \u003cp\u003eThe R\u0026amp;D effect after 2000 to yield stability under climate change was identified through a comparison of the \u003cem\u003efc\u003c/em\u003e and \u003cem\u003efc.em\u003c/em\u003e runs (Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003e). We assumed that the \u003cem\u003efc\u003c/em\u003e run, which included the time-dependent R\u0026amp;D and management variables, reflect advances in R\u0026amp;D and associated practices. Then, we performed a counterfactual run in which advances in R\u0026amp;D and management practices are lacking to varying degrees. To do this, we fixed the R\u0026amp;D and management variables to constant values. In the 1st model specification (M1), values of the R\u0026amp;D and management variables were set at the 2000 levels to calculate yield anomalies from 2000 to 2019, while in the 5th model specification (M5), they were set equal to their yearly values until 2004 and were held constant thereafter (fig. S14). This model setting resulted in 20 model specifications per crop (M1 to M20). For each year in which advances in R\u0026amp;D and practices were eliminated in the RF models, a histogram of the slope value for the corresponding model specification was derived and compared with a histogram of the slope value for actual R\u0026amp;D and practices, each consisting of 100 members from the RF models forced by d4PDF factual climate simulation (fig. S15). When the average slope value for actual R\u0026amp;D and practice was significantly smaller than that without advances in R\u0026amp;D and practice, it indicates that changes in R\u0026amp;D and practices after 2000 increased yield stability. The absence of a significant difference indicates that changes in R\u0026amp;D and practices since 2000 did not have a significant effect on increasing yield stability under climate change.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec29\" class=\"Section2\"\u003e \u003ch2\u003eSensitivity analysis\u003c/h2\u003e \u003cp\u003eWe examined the sensitivity of variable importance and attribution results to the different assumptions used in the R\u0026amp;D variable calculation. The country annual R\u0026amp;D expenditure data collected were for agriculture and not limited to the four crops considered here. Therefore, we first examined if the use of R\u0026amp;D expenditure for these crops, not for agriculture, affected results. To that end, we assumed that annual R\u0026amp;D expenditure designated for the four crops was proportional to the share of value derived from production of the four crops to total value of agricultural production. We collected country annual data on value of production for the four crops and entire agriculture from FAO (\u003cspan citationid=\"CR47\" class=\"CitationRef\"\u003e47\u003c/span\u003e) and derived them for each country and year. Second, we considered shorter (6 years) and longer (12 years) time lags between research and adoption of technologies by producers. Finally, we considered relatively slower (10% per year) and faster (20% per year) technology obsolescence by setting \u003cem\u003eδ\u003c/em\u003e value in Eq.\u0026nbsp;(\u003cspan refid=\"Equ1\" class=\"InternalRef\"\u003e1\u003c/span\u003e) to 0.1 and 0.2, respectively.\u003c/p\u003e \u003cp\u003eWe performed five experiments in the sensitivity analysis: Ex. 1\u0026mdash;R\u0026amp;D expenditure for agriculture (c0), shorter time lag between research and adoption (r6) and slower technology obsolescence (o1) (denoted \u0026ldquo;c0r6o1\u0026rdquo;); Ex. 2\u0026mdash;R\u0026amp;D expenditure for agriculture, longer time lag (r12) and slower obsolescence (\u0026ldquo;c0r12o1\u0026rdquo;); Ex. 3\u0026mdash;R\u0026amp;D expenditure for the four crops (c1), shorter time lag and slower obsolescence (\u0026ldquo;c1r6o1\u0026rdquo;); Ex. 4\u0026mdash;R\u0026amp;D expenditure for the four crops, longer time lag and slower obsolescence (\u0026ldquo;c1r12o1\u0026rdquo;); and Ex. 5\u0026mdash;R\u0026amp;D expenditure for the four crops, shorter time lag and faster obsolescence (\u0026ldquo;c1r6o2\u0026rdquo;). For each experiment, the RF models were built and then the variable importance and the area shares for the climate and R\u0026amp;D effect categories were estimated.\u003c/p\u003e \u003cp\u003eIn the main text, we presented the results from Ex. 1 (c0r6o1). We found that the assumptions used in the R\u0026amp;D variable calculation affected the level of but not the pattern of cumulative R\u0026amp;D expenditure (fig. S16). The difference in the level of cumulative R\u0026amp;D expenditure was absorbed in the RF models, as ordinary regression models would do, and did not affect the variable importance and attribution results.\u003c/p\u003e \u003c/div\u003e"},{"header":"Declarations","content":"\u003cp\u003e\u003cstrong\u003eFunding\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eTI was supported by the Environment Research and Technology Development Fund (Grant Number: JPMEERF23S21120) of the Environmental Restoration and Conservation Agency provided by the Ministry of the Environment of Japan; and Grants-in-Aid for Scientific Research were received from the Japan Society for the Promotion of Science (Grant Numbers: 22H00577 and 23H00351). YI and HS were supported by the Advanced Studies of Climate Change Projection (SENTAN, Grant Number: JPMXD0722680395) of the Ministry of Education, Culture, Sports, Science and Technology of Japan. DM and TI were supported by INRAE-NARO Linkage Call 2023.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAuthor contributions\u0026nbsp;\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eTI conceptualized, developed methodology, and conducted formal analysis; TS, YM, KO, and TT curated crop data and assisted with the investigation; HS and YI curated climate data and assisted with the investigation; DM assisted with methodology development; all authors contributed to writing the original draft.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eData availability\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eThe country yield and area harvested data from the FAO statistical database can be obtained from https://www.fao.org/faostat/en/#home. Gross domestic R\u0026amp;D spending data for some developed countries can be obtained from the OECD’s MSTI database: https://data-explorer.oecd.org/. Agricultural R\u0026amp;D expenditure data for some developing countries are accessible from the IFPRI’s ASTI database: https://doi.org/10.7910/DVN/9OXBIB. The grid yield dataset is accessible at the Data Integration and Analysis System (DIAS): https://doi.org/10.20783/DIAS.564. The actual climate data and the bias-corrected d4PDF factual and counterfactual climate simulations can be respectively obtained from https://doi.org/10.20783/DIAS.523 and https://doi.org/10.20783/DIAS.544.\u003c/p\u003e"},{"header":"References","content":"\u003col\u003e\n \u003cli\u003eFAO, IFAD, UNICEF, WFP, WHO. \u003cem\u003eThe State of Food Security and Nutrition in the World 2020\u003c/em\u003e (FAO, 2020). https://doi.org/10.4060/ca9692en\u003c/li\u003e\n \u003cli\u003eN\u0026oacute;ia J\u0026uacute;nior, R. de S., Ewert, F., Webber, H., Martre, P., Hertel, T.W., van Ittersum, M.K., Asseng, S. Needed global wheat stock and crop management in response to the war in Ukraine. \u003cem\u003eGlob. Food Secur.\u003c/em\u003e\u003cstrong\u003e35\u003c/strong\u003e, 100662 (2022). https://doi.org/10.1016/j.gfs.2022.100662.\u003c/li\u003e\n \u003cli\u003eAday, S., Aday, M.S. Impact of COVID-19 on the food supply chain. \u003cem\u003eFood Qual. Saf.\u0026nbsp;\u003c/em\u003e\u003cstrong\u003e4\u003c/strong\u003e, 2020, 167\u0026ndash;180. https://doi.org/10.1093/fqsafe/fyaa024\u003c/li\u003e\n \u003cli\u003eGaupp, F. Extreme Events in a Globalized Food System. \u003cem\u003eOne Earth\u003c/em\u003e\u003cstrong\u003e2\u003c/strong\u003e, 518\u0026ndash;521 (2020). https://doi.org/10.1016/j.oneear.2020.06.001\u003c/li\u003e\n \u003cli\u003eZhu, P., Burney, J., Chang, J. et al. Warming reduces global agricultural production by decreasing cropping frequency and yields. \u003cem\u003eNat. Clim. Chang\u003c/em\u003e\u003cstrong\u003e12\u003c/strong\u003e, 1016\u0026ndash;1023 (2022). https://doi.org/10.1038/s41558-022-01492-5\u003c/li\u003e\n \u003cli\u003eWei, D. et al. Key role of planted and harvested area fluctuations in US crop production shocks. \u003cem\u003eNat. Sustain.\u003c/em\u003e\u003cstrong\u003e6\u003c/strong\u003e, 1177\u0026ndash;1185 (2023). https://doi.org/10.1038/s41893-023-01152-2\u003c/li\u003e\n \u003cli\u003eBezner Kerr, R. et al. \u0026ldquo;Food, Fibre, and Other Ecosystem Products\u0026rdquo; in \u003cem\u003eClimate Change 2022: Impacts, Adaptation and Vulnerability\u003c/em\u003e (Cambridge Univ. Press, 2022). doi:10.1017/9781009325844.007\u003c/li\u003e\n \u003cli\u003eIizumi, T., Ramankutty, N. Changes in yield variability of major crops for 1981\u0026ndash;2010 explained by climate change. \u003cem\u003eEnviron. Res. Lett.\u003c/em\u003e\u003cstrong\u003e11\u003c/strong\u003e, 034003 (2016). https://dx.doi.org/10.1088/1748-9326/11/3/034003\u003c/li\u003e\n \u003cli\u003eBen-Ari, T., Bo\u0026eacute;, J., Ciais, P. et al. Causes and implications of the unforeseen 2016 extreme yield loss in the breadbasket of France. \u003cem\u003eNat Commun.\u003c/em\u003e\u003cstrong\u003e9\u003c/strong\u003e, 1627 (2018). https://doi.org/10.1038/s41467-018-04087-x\u003c/li\u003e\n \u003cli\u003eThirumalai, K. et al. Extreme temperatures in Southeast Asia caused by El Ni\u0026ntilde;o and worsened by global warming. \u003cem\u003eNat. Commun.\u003c/em\u003e\u003cstrong\u003e8\u003c/strong\u003e, 15531 (2017). https://doi.org/10.1038/ncomms15531\u003c/li\u003e\n \u003cli\u003eBeillouin, D., Schauberger, B., Bastos, A., Ciais, P., Makowski, D. Impact of extreme weather conditions on European crop production in 2018. \u003cem\u003ePhil. Trans. R. Soc. B\u003c/em\u003e\u003cstrong\u003e375\u003c/strong\u003e, 20190510 (2020). https://doi.org/10.1098/rstb.2019.0510\u003c/li\u003e\n \u003cli\u003eWhite, R. H. et al. The unprecedented Pacific Northwest heatwave of June 2021. \u003cem\u003eNat. Commun.\u003c/em\u003e\u003cstrong\u003e14\u003c/strong\u003e, 727 (2023). https://doi.org/10.1038/s41467-023-36289-3\u003c/li\u003e\n \u003cli\u003eCooper, M., Messina, C. D. Breeding crops for drought-affected environments and improved climate resilience. \u003cem\u003eThe Plant Cell\u003c/em\u003e\u003cstrong\u003e35\u003c/strong\u003e, 162\u0026ndash;186 (2023). https://doi.org/10.1093/plcell/koac321\u003c/li\u003e\n \u003cli\u003eBerrang-Ford, L. et al. A systematic global stocktake of evidence on human adaptation to climate change. \u003cem\u003eNat. Clim. Chang.\u003c/em\u003e\u003cstrong\u003e11\u003c/strong\u003e, 989\u0026ndash;1000 (2021). https://doi.org/10.1038/s41558-021-01170-y\u003c/li\u003e\n \u003cli\u003eLobell, D. B. et al. Greater Sensitivity to Drought Accompanies Maize Yield Increase in the U.S. Midwest. \u003cem\u003eScience\u003c/em\u003e\u003cstrong\u003e344\u003c/strong\u003e, 516\u0026ndash;519 (2014). doi:10.1126/science.1251423\u003c/li\u003e\n \u003cli\u003eMiller, S. J., Clarke, E., Mathews, S. L. Little evidence of avoided yield loss in US corn when short-term forecasts correctly predict extreme heat. \u003cem\u003eEnviron. Res. Lett.\u003c/em\u003e\u003cstrong\u003e18\u003c/strong\u003e, 124041 (2023). https://dx.doi.org/10.1088/1748-9326/ad0bd5\u003c/li\u003e\n \u003cli\u003eMizuta, R. et al. Over 5,000 years of ensemble future climate simulations by 60-km global and 20-km regional atmospheric models. \u003cem\u003eBull. Amer. Meteor. Soc.\u003c/em\u003e\u003cstrong\u003e98\u003c/strong\u003e, 1383\u0026ndash;1398 (2017). https://doi.org/10.1175/BAMS-D-16-0099.1\u003c/li\u003e\n \u003cli\u003eOtto, F. E. L. Attribution of weather and climate events.\u0026nbsp;\u003cem\u003eAnnu. Rev. Environ. Resour.\u0026nbsp;\u003c/em\u003e\u003cstrong\u003e42\u003c/strong\u003e, 627\u0026ndash;646 (2017). https://doi.org/10.1146/annurev-environ-102016-060847\u003c/li\u003e\n \u003cli\u003eGillett, N. P. et al. Constraining human contributions to observed warming since the pre-industrial period. \u003cem\u003eNat. Clim. Chang\u003c/em\u003e\u003cstrong\u003e11\u003c/strong\u003e, 207\u0026ndash;212 (2021). https://doi.org/10.1038/s41558-020-00965-9\u003c/li\u003e\n \u003cli\u003eRomanovska, P., Undorf, S., Schauberger, B., Duisenbekova, A., Gornott, C. 2024. Human-induced climate change has decreased wheat production in northern Kazakhstan. \u003cem\u003eEnvironmental Research: Climate\u003c/em\u003e\u003cstrong\u003e3\u003c/strong\u003e, 031005. https://dx.doi.org/10.1088/2752-5295/ad53f7\u003c/li\u003e\n \u003cli\u003eFunk, C. et al. Anthropogenic enhancement of moderate-to-strong El Ni\u0026ntilde;o events likely contributed to drought and poor harvests in southern Africa during 2016. \u003cem\u003eBull. Amer. Meteor. Soc.\u003c/em\u003e\u003cstrong\u003e99\u003c/strong\u003e, S91\u0026ndash;S96 (2018). https://doi.org/10.1175/BAMS-D-17-0112.1\u003c/li\u003e\n \u003cli\u003eVerschuur, J. et al. Climate change as a driver of food insecurity in the 2007 Lesotho-South Africa drought. \u003cem\u003eSci. Rep.\u003c/em\u003e\u003cstrong\u003e11\u003c/strong\u003e, 3852 (2021). https://doi.org/10.1038/s41598-021-83375-x\u003c/li\u003e\n \u003cli\u003eFAO. \u003cem\u003eThe Impact of Disasters on Agriculture and Food Security 2023 \u0026ndash; Avoiding and reducing losses through investment in resilience\u003c/em\u003e (FAO, 2023). https://doi.org/10.4060/cc7900en\u003c/li\u003e\n \u003cli\u003eRay, D.K., Mueller, N.D., West, P.C., Foley, J.A. (2013). Yield Trends Are Insufficient to Double Global Crop Production by 2050. PLoS ONE \u003cstrong\u003e8\u003c/strong\u003e, e66428. https://doi.org/10.1371/journal.pone.0066428\u003c/li\u003e\n \u003cli\u003eFuglie, K.O., Morgan, S., Jelliffe, J. (2024). World agricultural production, resource use, and productivity, 1961\u0026ndash;2020 (Report No. EIB-268). U.S. Department of Agriculture, Economic Research Service.\u003c/li\u003e\n \u003cli\u003eDeBruin, J.L., Schussler, J.R., Mo, H. and Cooper, M. (2017), Grain Yield and Nitrogen Accumulation in Maize Hybrids Released during 1934 to 2013 in the US Midwest. Crop Science, 57: 1431-1446. https://doi.org/10.2135/cropsci2016.08.0704\u003c/li\u003e\n \u003cli\u003eH.E. Ahrends et al. 2018 Genetic yield gains of winter wheat in Germany over more than 100 years (1895\u0026ndash;2007) under contrasting fertilizer applications. Environmental Research Letters 13, 104003. https://dx.doi.org/10.1088/1748-9326/aade12\u003c/li\u003e\n \u003cli\u003eShew, A.M., Tack, J.B., Nalley, L.L. et al. Yield reduction under climate warming varies among wheat cultivars in South Africa. Nat Commun 11, 4408 (2020). https://doi.org/10.1038/s41467-020-18317-8\u003c/li\u003e\n \u003cli\u003eSchmidt, M., Felsche, E. The effect of climate change on crop yield anomaly in Europe. \u003cem\u003eClim. Resilience Sustainability\u003c/em\u003e\u003cstrong\u003e3\u003c/strong\u003e, e61 (2024). https://doi.org/10.1002/cli2.61\u003c/li\u003e\n \u003cli\u003eKim, W., Iizumi, T., Nishimori, M. Global patterns of crop production losses associated with droughts from 1983 to 2009. \u003cem\u003eJ. Appl. Meteorol.\u0026nbsp;\u003c/em\u003e\u003cem\u003eClimatol.\u003c/em\u003e\u003cstrong\u003e58\u003c/strong\u003e, 1233-1244 (2019). https://doi.org/10.1175/JAMC-D-18-0174.1\u003c/li\u003e\n \u003cli\u003eVogel, E. et al. The effects of climate extremes on global agricultural yields. \u003cem\u003eEnviron. Res. Lett.\u003c/em\u003e\u003cstrong\u003e14\u003c/strong\u003e, 054010 (2019). https://dx.doi.org/10.1088/1748-9326/ab154b\u003c/li\u003e\n \u003cli\u003eTigchelaar, M., Battisti, D.S., Naylor, R.L., Ray, D.K. Future warming increases probability of globally synchronized maize production shocks. \u003cem\u003eProc. Natl. Acad. Sci. USA\u003c/em\u003e\u003cstrong\u003e115\u003c/strong\u003e, 6644\u0026ndash;6649. https://doi.org/10.1073/pnas.1718031115\u003c/li\u003e\n \u003cli\u003eKim, Y.-U., Webber, H., Adiku, S.G.K., N\u0026oacute;ia J\u0026uacute;nior, R.d.S., Deswarte, J.-C., Asseng, S., Ewert, F., 2024. Mechanisms and modelling approaches for excessive rainfall stress on cereals: Waterlogging, submergence, lodging, pests and diseases. Agric. For. Meteorol. 344, 109819. https://doi.org/10.1016/j.agrformet.2023.109819\u003c/li\u003e\n \u003cli\u003eM\u0026uuml;ller, C., et al. Global patterns of crop yield stability under additional nutrient and water inputs. \u003cem\u003ePLoS ONE\u003c/em\u003e\u003cstrong\u003e13\u003c/strong\u003e, e0198748 (2018). https://doi.org/10.1371/journal.pone.0198748\u003c/li\u003e\n \u003cli\u003eAhrends, H. E. et al. Nutrient supply affects the yield stability of major European crops\u0026mdash;a 50 year study. \u003cem\u003eEnviron. Res. Lett.\u003c/em\u003e\u003cstrong\u003e16\u003c/strong\u003e, 014003 (2021). https://dx.doi.org/10.1088/1748-9326/abc849\u003c/li\u003e\n \u003cli\u003eQiao, L. et al. Soil quality both increases crop production and improves resilience to climate change. \u003cem\u003eNat. Clim. Chang.\u003c/em\u003e\u003cstrong\u003e12\u003c/strong\u003e, 574\u0026ndash;580 (2022). https://doi.org/10.1038/s41558-022-01376-8\u003c/li\u003e\n \u003cli\u003eMehta, P. et al. Half of twenty-first century global irrigation expansion has been in water-stressed regions. \u003cem\u003eNat. Water\u003c/em\u003e\u003cstrong\u003e2\u003c/strong\u003e, 254\u0026ndash;261 (2024). https://doi.org/10.1038/s44221-024-00206-9\u003c/li\u003e\n \u003cli\u003eZaveri, E., Lobell, D. B. The role of irrigation in changing wheat yields and heat sensitivity in India.\u0026nbsp;\u003cem\u003eNat. Commun.\u003c/em\u003e\u003cstrong\u003e10\u003c/strong\u003e, 4144 (2019). https://doi.org/10.1038/s41467-019-12183-9\u003c/li\u003e\n \u003cli\u003eIFPRI. \u003cem\u003eAgricultural Science and Technology Indicators: 2019 Global Food Policy Report Annex Table 1\u003c/em\u003e (IFPRI, 2019). https://doi.org/10.7910/DVN/9OXBIB\u003c/li\u003e\n \u003cli\u003eNalley L, Dixon B, Chaminuka P, Naledzani Z, Coale MJ (2018) The role of public wheat breeding in reducing food insecurity in South Africa. \u003cem\u003ePLoS ONE\u003c/em\u003e\u003cstrong\u003e13\u003c/strong\u003e, e0209598. https://doi.org/10.1371/journal.pone.0209598\u003c/li\u003e\n \u003cli\u003eWu, H. et al. AsiaRiceYield4km: seasonal rice yield in Asia from 1995 to 2015. \u003cem\u003eEarth Syst. Sci. Data\u003c/em\u003e\u003cstrong\u003e15\u003c/strong\u003e, 791\u0026ndash;808 (2023). https://doi.org/10.5194/essd-15-791-2023\u003c/li\u003e\n \u003cli\u003eRuane, A.C. et al. 2021. Strong regional influence of climatic forcing datasets on global crop model ensembles. Agricultural and Forest Meteorology, 300, 108313. https://doi.org/10.1016/j.agrformet.2020.108313\u003c/li\u003e\n \u003cli\u003eMengel, M., Treu, S., Lange, S., and Frieler, K.: ATTRICI v1.1 \u0026ndash; counterfactual climate for impact attribution, Geosci. Model Dev., 14, 5269\u0026ndash;5284, https://doi.org/10.5194/gmd-14-5269-2021, 2021.\u003c/li\u003e\n \u003cli\u003eSloat, L.L., Davis, S.J., Gerber, J.S. et al. Climate adaptation by crop migration. Nat Commun 11, 1243 (2020). https://doi.org/10.1038/s41467-020-15076-4\u003c/li\u003e\n \u003cli\u003eG. Zhang et al. 2017 Spatiotemporal patterns of paddy rice croplands in China and India from 2000 to 2015. Science of The Total Environment, 579, 82-92. https://doi.org/10.1016/j.scitotenv.2016.10.223.\u003c/li\u003e\n \u003cli\u003eMaggi, F. et al. PEST-CHEMGRIDS, global gridded maps of the top 20 crop-specific pesticide application rates from 2015 to 2025. \u003cem\u003eSci. Data\u003c/em\u003e\u003cstrong\u003e6\u003c/strong\u003e, 170 (2019). https://doi.org/10.1038/s41597-019-0169-4\u003c/li\u003e\n \u003cli\u003eFAO. \u003cem\u003eFAOSTAT\u003c/em\u003e (FAO, 2023). https://www.fao.org/faostat/en/#home\u003c/li\u003e\n \u003cli\u003eIizumi, T., Sakai, T. The global dataset of historical yields for major crops 1981\u0026ndash;2016.\u0026nbsp;\u003cem\u003eSci. Data\u003c/em\u003e\u003cstrong\u003e7\u003c/strong\u003e, 97 (2020). https://doi.org/10.1038/s41597-020-0433-7\u003c/li\u003e\n \u003cli\u003eSacks, W. J., Deryng, D., Foley, J. A. Ramankutty, N. Crop planting dates: an analysis of global patterns. \u003cem\u003eGlob. Ecol. Biogeogr.\u003c/em\u003e\u003cstrong\u003e19\u003c/strong\u003e, 607-620 (2010). https://doi.org/10.1111/j.1466-8238.2010.00551.x\u003c/li\u003e\n \u003cli\u003eIizumi, T. et al. Increasing heavy rainfall events and associated excessive soil water threaten a protein-source legume in dry environments of West Africa. \u003cem\u003eAgric. For. Meteorol.\u003c/em\u003e\u003cstrong\u003e344\u003c/strong\u003e, 109783 (2024). https://doi.org/10.1016/j.agrformet.2023.109783\u003c/li\u003e\n \u003cli\u003eIshii, M., Mori, N. d4PDF: large-ensemble and high-resolution climate simulations for global warming risk assessment. \u003cem\u003eProg. Earth Planet. Sci.\u003c/em\u003e\u003cstrong\u003e7\u003c/strong\u003e, 58 (2020). https://doi.org/10.1186/s40645-020-00367-7\u003c/li\u003e\n \u003cli\u003eMizuta, R. et al. Climate simulations using MRI-AGCM3.2 with 20-km grid. \u003cem\u003eJ. Meteorol. Soc. Jpn\u003c/em\u003e\u003cstrong\u003e90A\u003c/strong\u003e, 233-258 (2012). https://doi.org/10.2151/jmsj.2012-A12\u003c/li\u003e\n \u003cli\u003eShiogama, H. et al. Attributing historical changes in probabilities of record-breaking daily temperature and precipitation extreme events. \u003cem\u003eSOLA\u003c/em\u003e\u003cstrong\u003e12\u003c/strong\u003e, 225-231 (2016). https://doi.org/10.2151/sola.2016-045\u003c/li\u003e\n \u003cli\u003eImada, Y., Watanabe, M., Kawase, H., Shiogama, H., Arai, M. The July 2018 high temperature event in Japan could not have happened without human-induced global warming.\u0026nbsp;\u003cem\u003eSOLA\u003c/em\u003e\u003cstrong\u003e15A\u003c/strong\u003e, 8-12 (2019). https://doi.org/10.2151/sola.15A-002\u003c/li\u003e\n \u003cli\u003eKawase, H. et al. Contribution of historical global warming to local-scale heavy precipitation in western Japan estimated by large ensemble high-resolution simulations. \u003cem\u003eJ. Geophys. Res. Atmos.\u003c/em\u003e\u003cstrong\u003e124\u003c/strong\u003e, 6093\u0026ndash;6103 (2019). https://doi.org/10.1029/2018JD030155\u003c/li\u003e\n \u003cli\u003eImada, Y. et al. Advanced risk-based event attribution for heavy regional rainfall events. \u003cem\u003enpj Clim.\u0026nbsp;\u003c/em\u003e\u003cem\u003eAtmos. Sci.\u003c/em\u003e\u003cstrong\u003e3\u003c/strong\u003e, 37 (2020). https://doi.org/10.1038/s41612-020-00141-y\u003c/li\u003e\n \u003cli\u003eIizumi, T. et al. Crop production losses associated with anthropogenic climate change for 1981\u0026ndash;2010 compared with preindustrial levels. \u003cem\u003eInt. J. Climatol.\u003c/em\u003e\u003cstrong\u003e38\u003c/strong\u003e, 5405\u0026ndash;5417 (2018). https://doi.org/10.1002/joc.5818\u003c/li\u003e\n \u003cli\u003eSultan, B., Defrance, D., Iizumi, T. Evidence of crop production losses in West Africa due to historical global warming in two crop models. \u003cem\u003eSci. Rep.\u003c/em\u003e\u003cstrong\u003e9\u003c/strong\u003e, 12834 (2019). https://doi.org/10.1038/s41598-019-49167-0\u003c/li\u003e\n \u003cli\u003eOECD. \u003cem\u003eMSTI database\u003c/em\u003e (OECD, 2023). https://stats.oecd.org/Index.aspx?DataSetCode=MSTI_PUB\u003c/li\u003e\n \u003cli\u003eWorld Bank. \u003cem\u003eWorld Development Indicators\u003c/em\u003e (World Bank, 2023) https://databank.worldbank.org/reports.aspx?source=world-development-indicators#\u003c/li\u003e\n \u003cli\u003eHipel, K. W., McLeod, A. I. \u003cem\u003eTime Series Modelling of Water Resources and Environmental Systems\u003c/em\u003e (Elsevier, 1994).\u003c/li\u003e\n \u003cli\u003eBreiman, L. Random forests. \u003cem\u003eMach Learn.\u0026nbsp;\u003c/em\u003e\u003cstrong\u003e45\u003c/strong\u003e, 5\u0026ndash;32 (2001). https://doi.org/ 10.1023/A: 10109 33404 324\u003c/li\u003e\n \u003cli\u003eFolberth, C. et al. Spatio-temporal downscaling of gridded crop model yield estimates based on machine learning. \u003cem\u003eAgric. For. Meteorol.\u003c/em\u003e\u003cstrong\u003e264\u003c/strong\u003e, 1\u0026ndash;15 (2019). https://doi.org/10.1016/j.agrformet.2018.09.021\u003c/li\u003e\n \u003cli\u003eChen, M., Guilpart, N., Makowski, D. Comparison of methods to aggregate climate data to predict crop yield: an application to soybean. Environmental Research Letters 19, 054049 (2024). https://dx.doi.org/10.1088/1748-9326/ad42b5\u003c/li\u003e\n \u003cli\u003eR Core Team. \u003cem\u003eR: a language and environment for statistical computing\u003c/em\u003e (R Foundation for Statistical Computing, 2023)\u0026nbsp;\u003c/li\u003e\n \u003cli\u003eLiaw, A., Wiener, M. Classification and regression by randomForest. \u003cem\u003eR News\u003c/em\u003e\u003cstrong\u003e2\u003c/strong\u003e, 18\u0026ndash;22 (2002). https://journal.r-project.org/articles/RN-2002-022/\u003c/li\u003e\n \u003cli\u003eWieder, W. R., Boehnert, J., Bonan, G. B. Evaluating soil biogeochemistry parameterizations in Earth system models with observations. \u003cem\u003eGlob. Biogeochem. Cycles\u003c/em\u003e\u003cstrong\u003e28\u003c/strong\u003e, 211\u0026ndash;222 (2014). https://doi.org/10.1002/2013GB004665\u003c/li\u003e\n \u003cli\u003eMcDermid, S. et al. Irrigation in the Earth system. \u003cem\u003eNat. Rev. Earth Environ.\u003c/em\u003e\u003cstrong\u003e4\u003c/strong\u003e, 435\u0026ndash;453 (2023). https://doi.org/10.1038/s43017-023-00438-5\u003c/li\u003e\n \u003cli\u003eSiebert, S. et al. A global data set of the extent of irrigated land from 1900 to 2005. \u003cem\u003eHydrol. Earth Syst. Sci.\u003c/em\u003e\u003cstrong\u003e19\u003c/strong\u003e, 1521\u0026ndash;1545 (2015). https://doi.org/10.5194/hess-19-1521-2015\u003c/li\u003e\n \u003cli\u003ePortmann, F. T., Siebert, S., D\u0026ouml;ll, P. MIRCA2000\u0026mdash;Global monthly irrigated and rainfed crop areas around the year 2000: A new high‐resolution data set for agricultural and hydrological modeling. \u003cem\u003eGlob. Biogeochem. Cycles\u003c/em\u003e\u003cstrong\u003e24\u003c/strong\u003e, GB1011 (2010). https://doi.org/10.1029/2008GB003435\u003c/li\u003e\n \u003cli\u003eIizumi, T. et al. Responses of crop yield growth to global temperature and socioeconomic changes. \u003cem\u003eSci. Rep.\u003c/em\u003e\u003cstrong\u003e7\u003c/strong\u003e, 7800 (2017). https://doi.org/10.1038/s41598-017-08214-4\u003c/li\u003e\n \u003cli\u003eBaldos, U.L.C., Fuglie, K.O., Hertel, T.W. The research cost of adapting agriculture to climate change: A global analysis to 2050. \u003cem\u003eAgric. Econ.\u0026nbsp;\u003c/em\u003e\u003cstrong\u003e51\u003c/strong\u003e, 207\u0026ndash;220 (2020). https://doi.org/10.1111/agec.12550\u003c/li\u003e\n \u003cli\u003eHollander, M., Wolfe, D. A., Chicken, E. \u003cem\u003eNonparametric Statistical Methods, Third Edition\u003c/em\u003e (John Wiley \u0026amp; Sons, 2015). doi:10.1002/9781119196037\u003c/li\u003e\n \u003cli\u003eYu, Q. et al. A cultivated planet in 2010 \u0026ndash; Part 2: The global gridded agricultural-production maps. \u003cem\u003eEarth Syst. Sci. Data\u003c/em\u003e\u003cstrong\u003e12\u003c/strong\u003e, 3545\u0026ndash;3572 (2020). https://doi.org/10.5194/essd-12-3545-2020\u003c/li\u003e\n\u003c/ol\u003e"},{"header":"Table","content":"\u003cp\u003eTable 1 is available in the Supplementary Files section\u003c/p\u003e"}],"fulltextSource":"","fullText":"","funders":[],"hasAdminPriorityOnWorkflow":false,"hasManuscriptDocX":true,"hasOptedInToPreprint":true,"hasPassedJournalQc":"","hasAnyPriority":false,"hideJournal":true,"highlight":"","institution":"","isAcceptedByJournal":true,"isAuthorSuppliedPdf":false,"isDeskRejected":"","isHiddenFromSearch":false,"isInQc":false,"isInWorkflow":false,"isPdf":false,"isPdfUpToDate":true,"isWithdrawnOrRetracted":false,"journal":{"display":true,"email":"[email protected]","identity":"researchsquare","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":true,"externalIdentity":"","sideBox":"","snPcode":"","submissionUrl":"/submission","title":"Research Square","twitterHandle":"researchsquare","acdcEnabled":true,"dfaEnabled":false,"editorialSystem":"","reportingPortfolio":"","inReviewEnabled":false,"inReviewRevisionsEnabled":true},"keywords":"","lastPublishedDoi":"10.21203/rs.3.rs-3895355/v2","lastPublishedDoiUrl":"https://doi.org/10.21203/rs.3.rs-3895355/v2","license":{"name":"CC BY 4.0","url":"https://creativecommons.org/licenses/by/4.0/"},"manuscriptAbstract":"\u003cp\u003eAgricultural research and development (R\u0026amp;D) has increased crop yields, but little is known about its ability to increase yield stability in the context of increasingly frequent extreme weather events. Using a grid yield dataset, we show that from 2000 to 2019, the standard deviation (SD) of yield anomalies for maize, rice, wheat and soybean, increased in 20% of the global harvested area. Based on random forest models relating yield anomaly to climate, soil, management and public R\u0026amp;D expenditure, we show that cumulative agricultural R\u0026amp;D expenditure, proportion of growing season exposed to optimal hourly temperatures, and dry and very wet days are key factors explaining crop yield variability. An attribution analysis based on large ensemble climate simulations with and without human influence on the global climate shows that unfavorable agro-climatic conditions due to climate change has increased SD, while higher R\u0026amp;D expenditure has led to more contrasting trends in SD over 2000\u0026ndash;2019. Although R\u0026amp;D has continued steadily in most countries, this study indicates that the progress made in R\u0026amp;D since 2000 may have lagged behind the unfavorable effect of climate change on yield variability.\u003c/p\u003e","manuscriptTitle":"Assessing the capacity of agricultural research and development to increase the stability of global crop yields under climate change","msid":"","msnumber":"","nonDraftVersions":[{"code":2,"date":"2025-03-18 09:10:28","doi":"10.21203/rs.3.rs-3895355/v2","editorialEvents":[{"type":"communityComments","content":0}],"status":"published","journal":{"display":true,"email":"[email protected]","identity":"researchsquare","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":true,"externalIdentity":"","sideBox":"","snPcode":"","submissionUrl":"/submission","title":"Research Square","twitterHandle":"researchsquare","acdcEnabled":true,"dfaEnabled":false,"editorialSystem":"","reportingPortfolio":"","inReviewEnabled":false,"inReviewRevisionsEnabled":true}},{"code":1,"date":"2024-02-27 15:30:55","doi":"10.21203/rs.3.rs-3895355/v1","editorialEvents":[{"type":"communityComments","content":0}],"status":"published","journal":{"display":true,"email":"[email protected]","identity":"researchsquare","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":true,"externalIdentity":"","sideBox":"","snPcode":"","submissionUrl":"/submission","title":"Research Square","twitterHandle":"researchsquare","acdcEnabled":true,"dfaEnabled":false,"editorialSystem":"","reportingPortfolio":"","inReviewEnabled":false,"inReviewRevisionsEnabled":true}}],"origin":"","ownerIdentity":"4e90951e-2a11-4e68-8d54-1e5b94a72fe7","owner":[],"postedDate":"March 18th, 2025","published":true,"recentEditorialEvents":[],"rejectedJournal":[],"revision":"","amendment":"","status":"posted","subjectAreas":[{"id":46314374,"name":"Agroecology"},{"id":46314375,"name":"Climate Analysis and Modeling"}],"tags":[],"updatedAt":"2025-03-27T16:22:04+00:00","versionOfRecord":{"articleIdentity":"rs-3895355","link":"https://doi.org/10.1093/pnasnexus/pgaf099","journal":{"identity":"pnas-nexus","isVorOnly":true,"title":"PNAS Nexus"},"publishedOn":"2025-03-22 00:00:00","publishedOnDateReadable":"March 22nd, 2025"},"versionCreatedAt":"2025-03-18 09:10:28","video":"","vorDoi":"10.1093/pnasnexus/pgaf099","vorDoiUrl":"https://doi.org/10.1093/pnasnexus/pgaf099","workflowStages":[]},"version":"v2","identity":"rs-3895355","journalConfig":"researchsquare"},"__N_SSP":true},"page":"/article/[identity]/[[...version]]","query":{"redirect":"/article/rs-3895355","identity":"rs-3895355","version":["v2"]},"buildId":"8U1c8b4HqxoKbykW_rLl7","isFallback":false,"isExperimentalCompile":false,"dynamicIds":[84888],"gssp":true,"scriptLoader":[]}