# Causes of variation among rice models in yield response to CO2 examined with Free-Air CO2 Enrichment and growth chamber experiments

Toshihiro Hasegawa, Tao Li, Xinyou Yin, Yan Zhu, Kenneth Boote, Jeffrey Baker, Simone Bregaglio, Samuel Buis, Roberto Confalonieri, Job Fugice, et al.

<figure><img src="image_1.png" /></figure>

## To cite this version:

Toshihiro Hasegawa, Tao Li, Xinyou Yin, Yan Zhu, Kenneth Boote, et al. Causes of variation among rice models in yield response to CO2 examined with Free-Air CO2 Enrichment and growth chamber experiments. Scientific Reports, 2017, 7 (1), 10.1038/s41598-017-13582-y. hal-01629858

HAL Id: hal-01629858

https://hal.science/hal-01629858

Submitted on 6 Nov 2017

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
SCIENTIFIC REPORTS

<figure><img src="image_3.png" /></figure>

OPEN

# Causes of variation among rice models in yield response to CO₂ examined with Free-Air CO₂ Enrichment and growth chamber experiments

Received: 30 June 2017

Accepted: 25 September 2017

Published online: 01 November 2017

Toshihiro Hasegawa¹, Tao Li², Xinyou Yin³, Yan Zhu⁴,⁵,⁶,⁷, Kenneth Boote⁸, Jeffrey Baker⁹, Simone Bregaglio¹⁰, Samuel Buis¹¹, Roberto Confalonieri¹², Job Fugice¹³, Tamon Fumoto¹⁴, Donald Gaydon¹⁵, Soora Naresh Kumar¹⁶, Tanguy Lafarge¹⁷,¹⁸, Manuel Marcaida III², Yuji Masutomi¹⁹, Hiroshi Nakagawa¹⁴, Philippe Oriol¹⁷,¹⁸, Françoise Ruget¹¹, Upendra Singh¹³, Liang Tang⁴,⁵,⁶,⁷, Fulu Tao²⁰,²¹, Hitomi Wakatsuki¹⁴, Daniel Wallach²², Yulong Wang²³, Lloyd Ted Wilson²⁴, Lianxin Yang²³, Yubin Yang²⁴, Hiroe Yoshida¹⁴, Zhao Zhang²⁵ & Jianguo Zhu²⁶

The CO₂ fertilization effect is a major source of uncertainty in crop models for future yield forecasts, but coordinated efforts to determine the mechanisms of this uncertainty have been lacking. Here, we studied causes of uncertainty among 16 crop models in predicting rice yield in response to elevated [CO₂] (E-[CO₂]) by comparison to free-air CO₂ enrichment (FACE) and chamber experiments. The model ensemble reproduced the experimental results well. However, yield prediction in response to E-[CO₂] varied significantly among the rice models. The variation was not random: models that overestimated at one experiment simulated greater yield enhancements at the others. The variation was not associated with model structure or magnitude of photosynthetic response to E-[CO₂] but was significantly associated with the predictions of leaf area. This suggests that modelled secondary effects of E-[CO₂] on morphological development, primarily leaf area, are the sources of model uncertainty. Rice morphological development

¹Tohoku Agricultural Research Center, National Agriculture and Food Research Organization, Morioka, Iwate, 020-0198, Japan. ²International Rice Research Institute, DAPO Box, 7777, Metro Manila, Philippines. ³Centre for Crop Systems Analysis, Wageningen University & Research, Wageningen, PO Box 430, The Netherlands. ⁴National Engineering and Technology Center for Information Agriculture, Nanjing Agricultural University, Nanjing, Jiangsu, 210095, China. ⁵Key Laboratory for Crop System Analysis and Decision Making, Ministry of Agriculture, Nanjing Agricultural University, Nanjing, Jiangsu 210095, China. ⁶Jiangsu Key Laboratory for Information Agriculture, Nanjing Agricultural University, Nanjing, Jiangsu 210095, China. ⁷Jiangsu Collaborative Innovation Center for Modern Crop Production, Nanjing Agricultural University, Nanjing, Jiangsu 210095, China. ⁸University of Florida, Gainesville, Florida, 32611-0500, USA. ⁹United States Department of Agriculture, Agricultural Research Service, Big Spring, Texas, 79720, USA. ¹⁰CREA, Research Center for Agriculture and Environment, Bologna, 40128, Italy. ¹¹INRA, UMR1114 EMMAH, F-84914, Avignon, France. ¹²University of Milan, Cassandra lab, Milan, 20133, Italy. ¹³International Fertilizer Development Center, Muscle Shoals, Alabama, 35662, USA. ¹⁴Institute for Agro-Environmental Sciences, National Agriculture and Food Research Organization, Tsukuba, Ibaraki, 305-8604, Japan. ¹⁵CSIRO Agriculture and Food, St Lucia, QLD, 4067, Australia. ¹⁶Indian Agricultural Research Institute, New Delhi, 110012, India. ¹⁷Cirad, UMR AGAP, F-34398, Montpellier, France. ¹⁸AGAP, Univ Montpellier, CIRAD, INRA, INRIA, Montpellier SupAgro, Montpellier, France. ¹⁹Ibaraki University, College of Agriculture, Inashiki, Ibaraki, 300-0393, Japan. ²⁰Chinese Academy of Sciences, Institute of Geographical Sciences and Natural Resources Research, Beijing, 100101, P.R. China. ²¹Natural Resources Institute Finland (Luke), FI-00790, Helsinki, Finland. ²²INRA, UMR AGIR, Castanet Tolosan, France. ²³Yangzhou University, Hanjiang, Yangzhou, Jiangsu, 225009, China. ²⁴Texas A&M AgriLife Research Center, Beaumont, Texas, 77701, USA. ²⁵State Key Laboratory of Earth Surface Processes and Resource Ecology, Beijing Normal University, Beijing, 100875, P.R. China. ²⁶State Key Laboratory of Soil and Sustainable Agriculture, Institute of Soil Sciences, Chinese Academy of Sciences, Nanjing, 210008, China. Correspondence and requests for materials should be addressed to T.H. (email: thase@affrc.go.jp) or T.L. (email: t.li@irri.org)
is conservative to carbon acquisition. Uncertainty will be reduced by incorporating this conservative nature of the morphological response to E-[CO₂] into the models. Nitrogen levels, particularly under limited situations, make the prediction more uncertain. Improving models to account for [CO₂] × N interactions is necessary to better evaluate management practices under climate change.

Climate change is a daunting challenge to world agriculture in meeting the demand for food and energy from a growing population¹. The production of rice, a staple food crop in Asia, is strongly affected by climate change². Atmospheric CO₂ concentration ([CO₂]) is an important driver for climate change. The most recent projections adopted by the Intergovernmental Panel on Climate Change show that [CO₂] may reach as high as about 540 µmol mol⁻¹ by 2050 and 940 µmol mol⁻¹ by 2100³ and will have significant impacts on global climate systems.

On the other hand, rising [CO₂] will have some positive effects on the productivity of C₃ crops like rice, by increasing photosynthesis, biomass and grain yield. This CO₂ fertilization will alleviate some of the negative effects of increasing temperatures⁴, but prior studies suggest that the magnitude of the CO₂ fertilization effect is highly uncertain⁵,⁶.

Crop models are pivotal tools for assessing the impact of climate change variables⁷, but models differ in their descriptions of processes and forcing variables⁸. Even for the major processes such as photosynthesis and biomass production that are directly affected by elevated [CO₂] (E-[CO₂]), some models use the coarse-grained concept of radiation use efficiency (RUE) whereas other models use a light response curve (LRC) of instantaneous leaf photosynthesis scaled to hourly or daily canopy photosynthesis and crop respiration. Some crop models simulate instantaneous leaf photosynthesis by use of the biochemical model of Farquhar, von Caemmerer & Berry (FvCB)⁹. Such diversity in model algorithms, combined with inconsistency in model parameterization procedures, can create a large range of uncertainties in model projections¹⁰⁻¹³.

Another direct effect of E-[CO₂] is on stomatal conductance and transpiration, which could affect crop production under water limited environments¹⁴. However, because rice is generally grown in fields where 5–10 cm of floodwater is maintained, the effect of E-[CO₂] on crop yield via stomatal conductance for irrigated rice is in general considered to be negligible. We therefore hypothesized that the difference in model types for the primary response of photosynthesis and biomass production to E-[CO₂] creates uncertainties in yield prediction.

Coordinated crop model intercomparison has recently become feasible through the Agricultural Model Intercomparison and Improvement Project (AgMIP)¹⁵. Our initial attempt has shown that rice crop models differ considerably in yield and biomass responses to [CO₂]¹², but the sources of uncertainties were not fully analyzed at that time, particularly lacking comparison with experimental observations under E-[CO₂]. There is an urgent need to identify the sources of uncertainty and to improve the methods of prediction of the effects of E-[CO₂]⁶.

Two types of experimental facilities have been exploited to determine the effects of E-[CO₂] on various plant traits; small-scale growth enclosures or chambers¹⁶ and large-scale field environments using FACE (Free-Air CO₂ Enrichment)¹⁷. There are arguments, however, that the experimentally observed CO₂ enhancement on yield and other plant traits differ between FACE and enclosed chambers¹⁸,¹⁹, suggesting that uncertainty also exists in the observations and that testing multiple models against multiple sources of experimental results is desirable.

The CO₂ fertilization effects on rice yield are not constant and also depend on other factors such as temperature²⁰,²¹ and nitrogen (N) availability²²,²³. Most rice models take account of these factors, but differences in how crop models simulate [CO₂] × N interactions are yet to be quantified. In this study, we compiled three years of FACE data conducted at two sites with contrasting climate conditions: one in northern Japan with a cool temperature climate²² and the other in central China with a subtropical climate. Each site had three N levels to test the [CO₂] × N interactions. We also compiled results from a series of Soil-Plant-Atmosphere-Research (SPAR) chamber experiments⁶,²⁴. The SPAR chamber experiments conducted in Gainesville, Florida, USA, included six levels of [CO₂], so are suitable to test model responses to various levels of [CO₂]. Another SPAR experiment with three levels of N was conducted to test the [CO₂] × N interactions. With these datasets and collectively 16 rice models, uncertainty analyses using a multi-model ensemble and multi-experimental facilities are possible.

Our objectives are: (i) to determine the variation of yield prediction among existing commonly used rice models of different types, (ii) to identify the processes/mechanisms/factors that cause this uncertainty, (iii) to analyze to what extent the models differ in simulating [CO₂] × N interactions, and (iv) to propose possible directions for model improvement.

## Results

**Modelled primary responses to changes in [CO₂].** All 16 models listed in Table 1 simulate leaf area index, canopy light interception, biomass production and grain yield¹²,¹³ and account for the direct effect of E-[CO₂] on photosynthetic rates or canopy radiation use efficiency (RUE). All models have a non-linear asymptotic response of leaf photosynthetic rate or canopy radiation use efficiency (RUE) to changes in [CO₂], but the type of algorithm used varies (Table 1, Methods) and the responses differ in magnitude (Fig. 1). At 567 µmol mol⁻¹ (averaged [CO₂] for all FACE treatments), increases in assimilation or RUE relative to that at 367 µmol mol⁻¹ (average daytime ambient [CO₂] in the FACE experiments) ranged from 14 to 30% for leaf photosynthesis and from 8 to 24% for canopy level RUE response (Fig. 1). Averaged responses were similar across different model types (16–20%), which agreed well with the summary of observations in the FACE and chamber experiments 20.

**Model simulations under ambient [CO₂].** Models were first run under ambient [CO₂] conditions and different N levels at the Shizukuishi (Japan) and Wuxi (China) FACE sites (Methods, Table S1). Observed yields of rice under ambient [CO₂] were generally greater at Wuxi than at Shizukuishi because of the higher N fertilization rates, longer duration, and higher yield potential of the cultivar used at Wuxi (Methods). Simulated yields of all models under ambient [CO₂] averaged 7.5 t ha⁻¹ and 9.3 t ha⁻¹ at Shizukuishi and Wuxi, respectively, in good
<table><thead><tr><th rowspan="3">Model</th><th colspan="4">Exercise simulated</th><th colspan="2">CO<sub>2</sub> response for primary production</th><th colspan="3">Leaf area increase</th><th colspan="3">Yield formation</th><th rowspan="3">Reference</th></tr><tr><th rowspan="2">FACE</th><th rowspan="2">SPAR</th><th colspan="2">Leaf-level</th><th rowspan="2">RUE</th><th rowspan="2">Other direct effects of [CO<sub>2</sub>]</th><th colspan="2">Resource-driven</th><th rowspan="2">Temperature-driven</th><th rowspan="2">Direct CO<sub>2</sub> effect on grain set or harvest index</th><th rowspan="2">Grain number</th><th rowspan="2">Partitioning coefficient</th></tr><tr><th>LRC</th><th>FvCB</th><th>Carbon</th><th>Nitrogen</th></tr></thead><tbody><tr><td>1. APSIM-ORYZA</td><td>●</td><td>●</td><td>●</td><td></td><td></td><td></td><td>●</td><td></td><td>●</td><td>●</td><td>●</td><td>●</td><td>47,48</td></tr><tr><td>2. CERES-RICE</td><td>●</td><td>●</td><td></td><td></td><td>●</td><td>Gs, Tr</td><td></td><td></td><td>●</td><td></td><td>●</td><td></td><td>49</td></tr><tr><td>3. DNDC-Rice</td><td>●</td><td>●</td><td></td><td>●</td><td></td><td>Gs, Tr</td><td></td><td>●</td><td></td><td></td><td></td><td>●</td><td>20,50,51</td></tr><tr><td>4. GECROS</td><td>●</td><td>●</td><td></td><td>●</td><td></td><td>Gs, Gm, Tr</td><td>●</td><td>●</td><td></td><td>●</td><td></td><td>●</td><td>52</td></tr><tr><td>5. GEMRICE</td><td>●</td><td>●</td><td>●</td><td></td><td></td><td>phenology, spikelet sterility</td><td></td><td></td><td>●</td><td>●</td><td>●</td><td></td><td>53–56</td></tr><tr><td>6. H/H</td><td>●</td><td>●</td><td></td><td>●</td><td></td><td>Gs</td><td></td><td>●</td><td></td><td></td><td>●</td><td></td><td>20,51,57</td></tr><tr><td>7. InfoCrop</td><td>●</td><td>●</td><td></td><td></td><td>●</td><td>Tr,</td><td>●</td><td>●</td><td></td><td></td><td>●</td><td></td><td>58</td></tr><tr><td>8. MATCRO<sup>§</sup></td><td>●</td><td>●</td><td></td><td>●</td><td></td><td>Gs</td><td></td><td></td><td></td><td></td><td>●</td><td></td><td>59–61</td></tr><tr><td>9. MCWLA<sup>§</sup></td><td>●</td><td>●</td><td></td><td>●</td><td></td><td>Gs, Tr</td><td></td><td></td><td>●</td><td></td><td>●</td><td></td><td>62</td></tr><tr><td>10. ORYZA2000</td><td>●</td><td>●</td><td>●</td><td></td><td></td><td></td><td>●</td><td></td><td></td><td></td><td>●</td><td></td><td>47</td></tr><tr><td>11. RiceGrow</td><td>●</td><td>●</td><td>●</td><td></td><td></td><td></td><td>●</td><td></td><td></td><td></td><td>●</td><td></td><td>63</td></tr><tr><td>12. RicePSM</td><td>●</td><td>●</td><td>●</td><td></td><td></td><td></td><td>●</td><td>●</td><td></td><td>●</td><td></td><td>●</td><td>64</td></tr><tr><td>13. SAMARA<sup>§</sup></td><td></td><td>●</td><td></td><td></td><td>●</td><td>Tr</td><td></td><td></td><td>●</td><td></td><td>●</td><td></td><td>65</td></tr><tr><td>14. SIMRIW<sup>§</sup></td><td>●</td><td>●</td><td></td><td></td><td>●</td><td>phenology, spikelet sterility</td><td></td><td></td><td>●</td><td>●</td><td>●</td><td></td><td>66,67</td></tr><tr><td>15. STICS</td><td>●</td><td></td><td></td><td></td><td>●</td><td>Gs, Tr</td><td></td><td></td><td>●</td><td></td><td>●</td><td></td><td>68,69</td></tr><tr><td>16. WARM<sup>§</sup></td><td>●</td><td>●</td><td></td><td></td><td>●</td><td></td><td></td><td></td><td></td><td></td><td>●</td><td></td><td>70,71</td></tr></tbody></table>

**Table 1.** Models used for the simulation exercise and their main characteristics. Gs, stomatal conductance; Gm, mesophyll conductance; Tr, transpiration. Leaf area increase; Resource-driven, dependent: LAI increases with the resource such as C and N allocated to the leaves; Temperature-driven, LAI increases without any effects of resource availability (only as a function of developmental stages or temperatures). Yield formation: Grain number, yield is calculated by grain number × individual grain weight; Partitioning coefficient, yield is calculated by biomass × harvest index; Direct CO₂ effect on grain set or harvest index, models that account for a direct effect of E-[CO₂] on grain set or harvest index. §The models that do not include the quantification of the effects of different nitrogen on crop growth and yield.

<figure><img src="image_4.png" /></figure>

**Figure 1.** [CO₂] response curves used for (a) leaf CO₂ assimilation rate or (b) radiation use efficiency used in 16 rice models (Primary CO₂ response). Values were scaled to that at 367 µmol mol⁻¹ (average daytime ambient [CO₂] in the FACE experiments). Each response was estimated under the following conditions: Photosynthetically active radiation, 2000 µmol m²·s⁻¹; Relative humidity, 70%, Air temperature, 25 °C; Wind speed, 1 m s⁻¹, Leaf N content, 2 gm⁻²; Leaf N concentration, 40 mg g⁻¹; Specific leaf mass, 200 cm² g⁻¹.
<figure><img src="image_5.png" /></figure>

**Figure 2.** Simulated and observed yield and biomass. Grain yield and biomass under ambient [CO₂] in (a) FACE and (b) SPAR experiments; % increase in yield and biomass in response to elevated [CO₂] in (c) FACE and (d) SPAR experiments. Each bar represents average across different N treatments and years at two FACE sites left (Shizukuishi and Wuxi) and in the SPAR chamber experiments (Exp. 1, 2 and 3). Error bars represent the maximum and minimum values of each measurement or simulation. Data sets used for calibration are not included. For % increase in elevated [CO₂] in the SPAR chambers (d), mean of values under 500 and 660 µmol mol⁻¹ relative to ambient are presented. Comparison by treatment in the FACE experiments are detailed in Figures S2 and S3.

agreement with the observed yields of 7.9 t ha⁻¹ and 9.4 t ha⁻¹ averaged over N treatments and years at each site (Fig. 2a). A pairwise T-test comparing simulations with observations for all data, excluding those used for model calibration, showed modestly significant differences in some instances, but the differences are less than the experimental standard errors (Detailed analyses are given in Tables S2 and S3). Model ensemble means agreed well with the observed means and variation in yield and biomass under different N and years (Figure S2). However, simulated yields of individual models for each N treatment ranged from 4.0 t ha⁻¹ to 11.1 t ha⁻¹ at Shizukuishi and from 5.8 t ha⁻¹ to 13.3 t ha⁻¹ at Wuxi (Fig. 2a). At both sites, yield and biomass varied substantially among years and N treatments but the model-to-model variation was the most important cause of variability, accounting for up to two-thirds of total variability (Table S5, Supplementary text for uncertainty analysis). The variance of the simulated yield and biomass among the models was 13–62 times greater than the residual variance of the experiments (Table S4).

Models were also run for reference [CO₂] (330 µmol mol⁻¹) for the SPAR chamber experiments conducted in Florida, USA (Methods). Simulated yield and biomass averaged over the three experiments (standard N only) differed from observed values by only 1.3 t/ha and 0.2 t/ha, respectively (Fig. 2b), and the differences were less than the standard error for the measurement mean (Table S5). The variance in simulated yield and biomass was mostly attributable to the variability between models (71% for grain yield and 57% for biomass, Table S5), and was 7–8 times greater than the variance between experiments (Tables S5, S6).

**Model simulations for the effects of elevated [CO₂].** Model simulations for the yield enhancement due to E-[CO₂]. Simulations were conducted under E-[CO₂] conditions at the two FACE sites. The mean [CO₂] in the FACE plots was 567 µmol mol⁻¹ (Methods). The observed increase due to E-[CO₂] on grain yield averaged 12.8% and was similar between the two sites: 12.4% at Shizukuishi and 13.3% at Wuxi (Fig. 2c). Simulated increase of grain yield due to E-[CO₂] averaged for all models was slightly greater than the measurements, by about 3 and 1 percent point at Shizukuishi and Wuxi, respectively. These differences are less than the observed residual standard error (see Figure S3 and Table S2). However, simulated yield increase due to E-[CO₂] varied considerably among models, treatments and years, ranging from -1 to 50% (Figs 2c, S4). Simulated increases in biomass in response to E-[CO₂] also varied greatly depending on model and treatment. The ranges of simulated and observed increases were greater at Shizukuishi than at Wuxi (Figs 2c and S4).

Between-model variation accounted for about 60-70% of the variation in simulated yield increase (Table S3). Simulated yield increases for individual models, averaged over treatments and years at each site, varied from 2
<figure><img src="image_6.png" /><figcaption>**Figure 3.** Comparison between yield response to elevated [CO₂] of 14 individual rice models (a) between two FACE sites and (b) between FACE and SPAR chamber experiments. Each point is an average over different N treatments and years for FACE. Red, LRC (light response curve-type photosynthesis model); Yellow, FvCB (Farquhar, von Caemmerer & Berry photosynthesis model); Green, RUE (Radiation use efficiency). Means of the observed and simulated are also shown. FACE = mean value for Shizukuishi and Wuxi. SPAR = mean value for 500 and 660 µmol mol⁻¹ (ex 1, 2 & 3). SD is the measure of variation within each model due to different conditions (N, years or site).</figcaption></figure>

to 38% at Shizukuishi (Fig. 3a, X-axis) and from 3 to 26% at Wuxi (Fig. 3a, Y-axis). Simulated yield increase due to E-[CO₂] at Shizukuishi was significantly correlated with that at Wuxi (*P* < 0.01, Fig. 3a) so that models that overestimated at one site also simulated greater yields at the other site, suggesting that the magnitude of simulated CO₂ fertilization effects by each model was consistent across the sites.

No clear difference was observed between model types, but substantial variation exists within each model type (Figs 3a and S4), suggesting that the model structure for primary CO₂ response is not the major source of variation.

Simulated yield and biomass by each model greatly differed among treatments (Figures S3, S4). However, the between-model variation accounted for about 30 to 70% of the total variation in simulated yields, suggesting that within-model variation is of major importance (Table S4).

Simulated yield and biomass increases in response to E-[CO₂] in the SPAR experiments (mean of values at 500 and 660 µmol mol⁻¹, comparable to the change of [CO₂] in the FACE studies), were also similar to the measured increases (Fig. 2d). The observed and simulated yield increases were generally greater in the SPAR chamber than in the FACE studies (Fig. 2c,d), although the simulations slightly overestimated the FACE response and underestimated the chamber response.

A significantly positive correlation exists in the modelled yield response between FACE and chamber studies (Fig. 3b, *P* < 0.05). This result supports the consistency of model performances under different experimental setups.

*Model simulations for the yield response to various [CO₂].* The SPAR chamber experiments tested six levels of [CO₂], ranging from 160 to 900 µmol mol⁻¹ for two growing seasons. The observed yields showed a typical asymptotic response to increasing [CO₂]: grain yield increased more with increasing [CO₂] within the sub-ambient [CO₂] range than the supra-ambient range, and above 500 µmol mol⁻¹ yield plateaued to a similar level of around 10 to 40% increase relative to yield at [CO₂] at 367 µmol mol⁻¹ (Fig. 4a). Model ensemble means agreed well with the response observed at different [CO₂] but with a large variation between different models. There was no clear difference between model types, although FvCB-type models tended to be more sensitive (responsive) to [CO₂]. Biomass followed a response pattern similar to grain yield with a similar variation among models (Fig. 4b).

*Model simulations for the yield response to various [CO₂] under different N levels.* Both FACE sites had three levels of N applications. However, because the N doses were based on the conventional practices at each site, they differed between the two sites (Methods). Nevertheless, the effect of [CO₂] fertilization on yield was similar at both sites, with the observed yield enhancement being smaller when N supply was limited (Fig. 5), as reported by Kim *et al.*²² for Shizukuishi and Yang *et al.*²³ for Wuxi. However, simulated yield enhancement did not show a sufficient degree of dependence on N levels (Fig. 5). At all N levels, variation in yield response among models was large but variation tended to be greater at lower N levels. There was no clear difference between the model types in magnitude or N dependence of the yield enhancements.

One of the SPAR experiments (Experiment 3) compared factorial combinations of two [CO₂] and three N levels. Enhancements due to E-[CO₂] were not different among N levels (Figure S5) for either simulated or observed yields. Model-to-model variation in yield response to E-[CO₂] was generally smaller at the highest N levels and increased as N decreased. No clear difference was observed between model types.

**Factors influencing the variation in simulated [CO₂] response.** Simulated yield response to E-[CO₂] was positively correlated with the simulated biomass response for both FACE (Figs 6a, S6a, *P* < 0.001) and SPAR
<figure><img src="image_7.png" /><figcaption>**Figure 4.** Simulated and observed response of (a) grain yield and (b) aboveground biomass to six [CO₂] levels. LRC, light response curve-type photosynthesis model; FvCB, Farquhar, von Caemmerer & Berry photosynthesis model; RUE, radiation use efficiency. Numbers in the figures represent the growth [CO₂] conditions. % changes are relative to the value at 367 µmol mol⁻¹. Solid curve represents the model ensemble mean.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 5.** Simulated (a) and observed (b) yield enhancements due to elevated [CO₂] under different N levels obtained at Shizukuishi FACE site in Japan and Wuxi site in China. LRC, light response curve-type photosynthesis model; FvCB, Farquhar, von Caemmerer & Berry photosynthesis model; RUE, radiation use efficiency. Note that five out 14 models do not take N limitation into account (Table 1) and simulated only under the highest N conditions at each site.</figcaption></figure>
<figure><img src="image_9.png" /><figcaption>**Figure 6.** Factors affecting the simulated increase in grain yield due to E-[CO₂] in the high N treatment in the FACE experiments at two sites. (a) grain yield increase versus biomass increase, (b) grain yield increase versus harvest index increase, (c) biomass increase versus primary [CO₂] increase (leaf CO₂ assimilation rate, CAR, or radiation use efficiency, RUE), and (d) biomass increase versus maximum LAI increase. The data from all N treatments and the chamber experiments are shown in Figures S6–10. LRC, light response curve-type photosynthesis model; FvCB, Farquhar, von Caemmerer & Berry photosynthesis model; RUE, radiation use efficiency. SZ = Shizukuishi, WX = Wuxi.</figcaption></figure>

chamber studies (Figure S7, P < 0.001). About 69% of the variation in the simulated yield increase in the FACE studies and 67–97% in the SPAR chamber experiments was attributed to the simulated response in biomass. On the other hand, variation among models in simulated harvest index was not or much less correlated with yield response than that in biomass response (Figs 6b, S6b, S8).

To identify factors affecting biomass increase, we examined the relationship between primary [CO₂] response of the model (Fig. 1) and simulated biomass increase (Figs 6c, S6c and S9), because this is the primary direct mechanism by which modelled biomass changes in response to different [CO₂]. However, variation in simulated primary response to different [CO₂] was not correlated with biomass response, except in sub-ambient [CO₂] treatments. A weak but significantly positive correlation was observed between primary response and biomass response to 160 and 250 µmol mol⁻¹ (Figure S9d and e).

Another determinant of biomass production is canopy light interception, which is largely determined by leaf area index (LAI). The variation in the biomass response to [CO₂] was closely and positively correlated with that of the maximum LAI response in both FACE and SPAR experiments and at all [CO₂] levels (Figs 6d, S6d and S10). Models differ in their LAI submodules and can be classified into four types (Table 1). Some models use a resource-driven approach for LAI growth, the classic approach being carbon (or biomass)-dependent modules. In this approach, carbon allocated to the leaves is converted to leaf area based on fixed area-to-weight ratios for leaf (C-dependent type). Other models use N as a key driver for leaf growth, in which leaf area increases in proportion to N uptake, or N allocated to the leaves (N-dependent type). In some models, leaf area growth is a function of temperature, where, for instance, leaf area increases in proportion to growing degree days, regardless of C or
<figure><img src="image_10.png" /><figcaption>**Figure 7.** Simulated response of maximum LAI, final biomass and grain yield at Shizukuishi (a) and Wuxi (b) for models with different approaches to calculating LAI. Data are expressed as a ratio to response at 367 µmol mol⁻¹. At each site, all variables are significantly different between model types (*P* < 0.05 for biomass at Wuxi and *P* < 0.001 for all other variables). The main effect of N or interactions between N and model type were not significant for any variables. Letters with different letters beside the bar indicate that the between-group is different at *P* = 0.05, by the Tukey b method (Wholly Significant Difference).</figcaption></figure>

N availability (Resource-independent, temperature-driven type). Some models use the carbon-driven approach combined with other limitations.

There were significant differences in the LAI response to E-[CO₂] among the model types. Carbon-dependent models showed as much as a 30% increase in maximum LAI in response to E-[CO₂], whereas no increase was recorded with the temperature-driven (resource-independent) type (Fig. 7). The other two types had intermediate responses. The large variation in maximum LAI response resulted in significant differences in responses in biomass and yield.

## Discussion

Crop response to increasing [CO₂] has been identified as a major source of uncertainty in predicting crop production 6. In this study, we confirmed substantial variation in yield prediction in response to E-[CO₂] among rice crop models, which was much greater than the experimental variation. Here we used data from FACE and SPAR chamber experiments to identify the sources of uncertainty. Experiments with small pots could limit the response to E-[CO₂]²⁵, but the SPAR chamber experiments used large containers containing soil of 50 cm depth, so the effect of the limited rooting zone was assumed to be minimized.

Despite the large differences between models, the mean computed over all model simulations was close to the observed yield and yield response to E-[CO₂], as previously reported¹⁰⁻¹²,²⁶. The observed yield increase in response to E-[CO₂] was greater in the chamber studies than in the FACE study, in agreement with the meta-analysis by Long *et al*.⁵. However, the mean over all models reproduced similar differences in yield response between FACE and chamber studies, although the mean slightly overestimated the FACE results and underestimated the chamber results. The reason why simulations followed the yield responses of both experimental setups could be partly due to growth conditions (temperature and N fertilization) in the SPAR chambers, with ample N under warmer air being beneficial for the CO₂ fertilization effects, compared to the two FACE studies. Such a difference could have caused the greater yield response to E-[CO₂] in the SPAR chambers and to some extent in the simulations. There are arguments as to whether results from FACE and enclosure studies are intrinsically different⁵,²⁷⁻²⁹, but the present study suggests that the differences between the two technologies are due to differences in environmental and management factors and that model ensembles can reproduce these differences.

Differences in predictions among models were relatively consistent between the two sites (Fig. 3a) and between FACE and chamber studies (Fig. 3b), unlike a case study on wheat²⁶. Models that tend to overestimate the increase in one environment are likely to overestimate in the other environments, which suggests that model uncertainties are not random, and that opportunities exist for reduction in model uncertainty through identifying the mechanisms underlying the consistent biases.

Yield is a result of biomass production and allocation to harvestable organs (harvest index). The present study showed that much of the variation in yield resulted from variation in total biomass. The primary response to [CO₂] is the direct mechanism by which biomass is influenced in different [CO₂] conditions, and the models differ in complexity or mechanisticity for the primary response. However, the model type (LRC, FvCB or RUE) did not account for differences in simulated biomass response to E-[CO₂] (Figs 6c and S9), which was not surprising since the magnitude of the primary response to E-[CO₂] did not differ among model types (17–21%, Figure S1).
Model primary response to E-[CO₂] (+200 µmol mol⁻¹) ranged from 8 to 30% among models (Fig. 1), but was not related to simulated biomass response to E-[CO₂]. On the other hand, simulated LAI in response to E-[CO₂] accounted for much of the model-to-model variation in the biomass response to E-[CO₂] (Figs 6d,S6d and S10 despite the fact that none of the models included direct mechanisms by which E-[CO₂] influences leaf area expansion. Nevertheless, increase in LAI in response to E-[CO₂] varied substantially among the models. This is in contrast to a recent study testing six wheat models against Australian FACE data³⁰, which found similar LAI increases among models of 7–24% due to E-[CO₂]. The 16 rice models tested in this study differed in LAI formulations: for instance, some models simulate LAI increase assuming leaf growth is limited by carbon allocated to the leaves, whereas some use a function of temperature only (Table 1). Because carbon gain is enhanced by E-[CO₂], LAI in the C-driven model is increased by E-[CO₂]. As a result, light interception increases, which enhances biomass production. This positive feedback results in a biomass increase that is much greater than the primary [CO₂] response. On the other hand, the advantage of greater light interception does not occur with resource-independent type models, so the biomass response is largely determined by the primary response to E-[CO₂]. In an earlier study with wheat, Ewert³¹ pointed out the importance of accurate estimation of LAI for simulating the effect of [CO₂] on canopy carbon gain. The present study shows that much more of the between-model variability in biomass in rice can be attributed to variability in simulated leaf area than to variability in primary response to E-[CO₂].

Variation in simulated leaf area among models can also be a source of uncertainty in estimating gas and energy exchange between land surface and atmosphere, which is one of the most important feedbacks to the climate system. Because agricultural land represents a large segment of terrestrial ecosystems, accurate estimation of gas and heat exchange in croplands is important. Crop canopy cover or LAI is one of the most important state variables for the gas and heat exchange process, but its response to climate change is a major source of uncertainty. Improving models to reduce the uncertainty in leaf area estimation is thus also important for the modelling terrestrial energy and water cycles.

Experimental observations in the FACE experiments showed that the response of maximum LAI to E-[CO₂] was generally below 10% (Fig. 6d)³²,³³, despite increased tiller numbers. Likewise in the SPAR chamber experiments, leaf area around flowering in E-[CO₂] (averaged for 500 and 660 µmol mol⁻¹) was greater than the control by about 6%¹⁶. For rice, laminar length becomes shorter under E-[CO₂]³⁴, which offsets the increase in the total leaf number. It is worth noting, however, that there are certain species that show substantially greater response of LAI to E-[CO₂], such as poplar³⁵ and some invasive weed species³⁶. Leaf area expansion for these species is likely driven by carbon supply. On the other hand, rice is generally conservative in leaf morphology and expansion as a function of [CO₂], and is mostly limited by N supply and temperature. Modern cultivars of major crops such as rice may have been bred to be less plastic morphologically in response to excess amounts of resources, in order to be adapted to increased fertilization and dense planting. In previous studies comparing old and new cultivars of wheat and oats³⁷,³⁸, morphological plasticity of modern cultivars in response to increased carbon gain has been reported to be less than that of cultivars released prior to the Green Revolution.

The weak correlation to biomass changes under sub-ambient [CO₂] conditions (Figure S9d,e) suggests the modest importance of carbon as a resource to determine biomass responses. Here we show a strong need to improve the simulation of leaf area production in response to E-[CO₂] to reduce uncertainties in predicting yield under E-[CO₂], which is more important than the types of photosynthetic responses to elevated [CO₂].

It is generally reported that N availability limits biomass response to elevated [CO₂]¹⁹,⁴⁰. In the rice FACE experiments, the observed yield response to E-[CO₂] was indeed limited under low N conditions²². Final biomass increased in E-[CO₂] in the FACE experiments, even in the low N treatment³², but harvest index decreased to reduce the positive effects of E-[CO₂] on grain yield⁴¹. The models failed to reproduce the limited response to E-[CO₂] at low N. Models differ in the way they treat the effect of N on biomass production and partitioning to grains, in particular under E-[CO₂], which resulted in a large variation in simulated yield response to E-[CO₂] under low N conditions both in the FACE and SPAR chamber experiments (Figs 5 and S5). There is a strong need to improve the accuracy of modeling [CO₂] and N interactions because N fertilization strategies are key issues for climate change adaptation. We need further experimental evidence on how dry matter partitioning is affected by E-[CO₂] under N-limited conditions to reduce model uncertainty.

This study shows that biomass response, and not harvest index, is the major source of between-model variability in yield response to E-[CO₂]. Previous experimental observations showed that response of harvest index to elevated [CO₂] varied depending on the extreme temperature conditions around flowering and grain filling²⁰,⁴². In the two FACE and one SPAR chamber experiments used in this study, temperatures at the critical period did not exceed the upper or lower thresholds triggering severe stresses on grain setting, which in turn could reduce harvest index. Currently, only a few models account for any potential direct effect of E-[CO₂] on grain-setting or harvest index (Table 1). We need to test model uncertainty under extreme temperatures, combined with different temperature and [CO₂] regimes, against different data sets that cover a range of temperature conditions.

In conclusion, predictions of rice yield and biomass in response to E-[CO₂] varied significantly among the 16 rice models - the variation was much greater than the experimental variation observed in the FACE and SPAR chamber experiments. Variation was not random; models that tended to overestimate the responses in one environment tended to overestimate in the other environments. We hypothesized that the variation between models could be accounted for by the type of algorithm used to describe photosynthetic or RUE response to E-[CO₂], because this is the only mechanism by which models take account of the direct effects of E-[CO₂]. However, the results showed that the variation in simulated yield was not associated with the primary response to E-[CO₂], but was significantly associated with the variation of modelled CO₂-related responses in LAI. This suggests that “modelled” secondary or indirect effects of E-[CO₂] on morphological development are the sources of uncertainty. In reality, rice morphology and LAI are conservative in response to changes in E-[CO₂]. Accounting for this more conservative nature in the models may reduce the uncertainty in biomass and yield prediction. Nitrogen
levels (particularly under limited N) make the prediction more uncertain. Improving models to account for the
[CO₂] × N interaction is therefore recommended to better evaluate management practices under climate change.

**Materials and Methods**

**Experimental data.** Data provided for the FACE studies were obtained from two rice FACE sites in Asia (Table S1); the Shizukuishi site (N39°38′, E140°57′) located in northern Japan and the Wuxi site (N 31°37′, E 120°28′) in central eastern China, both representing the typical rice growing regions in their respective countries. Details of experimental procedures are given in Okada *et al.*<sup>43</sup> and Kim *et al.*<sup>22</sup> for the Shizukuishi site and in Liu *et al.*<sup>44</sup> and Yang *et al.*<sup>33</sup> for the Wuxi site, also summarized in Table S1. Briefly, the FACE facilities were established in farmers' fields, and CO₂ concentration in octagonal plots (ring, hereafter) measuring 12 m in diameter was increased by fumigating CO₂ through emission tubes installed on the peripheries of the rings. Pure CO₂ was released from the windward side, with the proportional–integral–derivative algorithm controlling the pressure of CO₂ release<sup>43</sup>. Target [CO₂] at both sites was 200 µmol mol⁻¹ above the ambient [CO₂].

Within the CO₂ treatments, three levels of N fertilizer were included; low, standard and high, with the stand-
ard being similar to farmers' practice at each site. Because conventional (or standard) N input differed greatly
between the two sites, levels of N inputs also differed substantially; total N applied at the Shizukuishi site ranged
from 40 to 150 kg ha⁻¹, whereas applied N at the Wuxi site ranged from 150 to 350 kg ha⁻¹. Eighty-three per-
cent (in Shizukuishi) and 60% (in Wuxi) of the total N was applied as a basal dressing before transplanting
and as a side-dressing at early tillering, with the remaining N side-dressed at panicle initiation. Cultivars were
Akitakomachi at the Shizukuishi site and Wuxianjing 14 at the Wuxi site, which are both japonica types, planted
widely in the respective regions. Seedlings were transplanted at 4–5 leaf age and grown under flooded conditions
until 10–15 d before harvest. The planting density was 19 and 24 hills m⁻² at Shizukuishi and Wuxi, respectively,
with three seedlings per hill at both sites.

The chamber experimental data were from 1 m × 2 m sunlit, Soil Plant Atmosphere Research (SPAR) cham-
bers with high transmission of photosynthetically-active radiation<sup>45</sup> located in Gainesville, Florida, USA (N
29°37.8', W 82°22.2'). Radiation transmission was about 100%, especially considering that irradiance would also
have a higher diffuse fraction<sup>45</sup>. Details of the soil and cultivation practices are given by Baker *et al*.16. Briefly, the
soil was a Candler fine sand (sandy, siliceous, hyperthermic uncoated Typic Quartzsamments) placed in the soil
base 60 cm deep, with 50 cm of soil with a 5-cm paddy flood depth. There were no percolation losses of water or N,
as the chamber base was sealed. The crop was direct-seeded (over-seeded and thinned), grown under a constant
30 °C for the first 7–10 days until flood paddy and temperature-CO₂ treatments were established. The row spacing
was 18 cm.

Response to [CO₂] levels was tested with the following three chamber experiments using IR30, which is an *indica* cultivar. In Experiments 1 and 2¹⁶, plants were grown at six levels of [CO₂] at 160, 250, 330, 500, 660, and 900 µmol mol⁻¹ at a constant day/night air temperature of 31 °C and water temperature of 27 °C. Fertilization applications were aimed to be sufficient (up to 50% more than required for maximum growth). Nitrogen in the form of urea was applied in 4 splits; immediately prior to establishing the flood depth and driven into the soil by flood application (thus volatilization losses of N should be very limited). A total of 290 kg ha⁻¹ of N was applied for all treatments. In Experiment 3²⁴,⁴⁶, plants were grown at two levels of [CO₂] at 330 and 660 µmol mol⁻¹ at day/night air temperatures of 28 °C/21 °C and at a constant water temperature of 25 °C. At each [CO₂] level, three levels of N were imposed; 0, 100 and 200 kg ha⁻¹ in 3 splits.

**Crop models and model simulations.** Sixteen rice models were used in this study, but not all models sim-
ulated each experimental setup; 14 rice models were tested against the FACE data and 15 against the SPAR data
(Table S2). Detailed descriptions of the models are available in the references specified in the table and a summary
or classification of the models can be found in Li *et al.*¹² and Confalonieri *et al.*¹³. All models take into account
the effects of [CO₂] either on leaf-level [CO₂] assimilation or canopy-level conversion of radiation to biomass,
and can be classified into two groups; the first uses single-leaf photosynthesis with multiple canopy layers and the
second uses canopy-based radiation conversion or canopy photosynthesis by one big-leaf (Table 1). Models in
the single-leaf photosynthesis group can further be classified into two categories; one that uses an empirical light
response curve modified by [CO₂] (LRC) and the other that uses the Farquhar *et al.* model⁹ (FvCB). Many models
account for other direct effects of [CO₂] such as stomatal conductance and a few models also include direct effects
of [CO₂] on phenology and spikelet sterility.

Models also have large differences in the treatment of leaf area development, which is an important growth
process that strongly affects radiation interception and biomass production. In this study, we classified the leaf
area production submodules into two categories; (1) leaf area increase limited by resources such as carbon and/
or nitrogen allocated to the leaves and (2) leaf area increase independent of resource availability. The first type is
subdivided into C-limited and N-limited models, although some models include dependence on both C and N.
The second type often uses developmental stage or a temperature function to determine leaf area increase.

Specific simulation procedures for the FACE and SPAR chamber experimental conditions are described in
supplementary materials and outlined below.

*Model simulations for the FACE experiments.* Each modeling group received calibration data that covers all three N levels with both ambient and elevated [CO₂] treatments (Shizukuishi in 1999 and Wuxi in 2002). These data sets were used to test and improve model performance under E-[CO₂]. Each group could determine how to use the available data to perform calibration or modify the functions, and the final model equations and parameters were recorded and reported to the AgMIP-Rice team leader. Then the same final model, except for cultivar and soil parameters that were site dependent, was run for all N × [CO₂] combinations of field experiments conducted in 1998 and 2000 at Shizukuishi, and in 2001 and 2003 at Wuxi.
*Model simulations for the SPAR experiments.* Each modeling group received all required crop management, soil and weather data, phenology dates (anthesis and maturity), seasonal biomass accumulation and yield for the “ambient” [CO₂] (330 µmol mol⁻¹) treatment for all experiments. Modellers set life cycle and productivity traits for the ambient treatment, starting with the calibrated model based on FACE experiments. The model calibrated using the ambient SPAR treatment was then used to simulate all the [CO₂] levels and temperature × [CO₂] treatments of the SPAR study. For each modelling group, the model used to simulate the FACE and SPAR experiments used the same functions and parameters for the [CO₂] response.

FACE and SPAR experiments used different ambient control [CO₂]. To compare the enhancements by E-[CO₂] between the two types of experiments, response curves for the SPAR results were first estimated using non-linear regression for each output. Using those response curves, the output at the average ambient [CO₂] of the two FACE experimental sites (367 µmol mol⁻¹) was derived. Relative increases of the SPAR results were calculated based on those estimated reference values.

**Data availability.** The datasets analysed and all modifications made to individual models during the current study are available from the corresponding authors on reasonable request.

## References

1. Hatfield, J. L. et al. Climate impacts on agriculture: Implications for crop production. *Agron. J.* **103**, 351–370 (2011).

2. Wassmann, R. et al. Regional vulnerability of climate change impacts on Asian rice production and scope for adaptation. *Adv. Agron.* **102**, 91–133 (2009).

3. Meinshausen, M. et al. The RCP greenhouse gas concentrations and their extensions from 1765 to 2300. *Clim. Change* **109**, 213–241 (2011).

4. Parry, M. L., Rosenzweig, C., Iglesias, A., Livermore, M. & Fischer, G. Effects of climate change on global food production under SRES emissions and socio-economic scenarios. *Glob. Environ. Chang.* **14**, 53–67 (2004).

5. Long, S. P., Ainsworth, E. A., Leakey, A. D. B., Nosberger, J. & Ort, D. R. Food for thought: lower than expected crop yield stimulation with rising CO₂ concentrations. *Science* **(80–).**, **312**, 1918–1921 (2006).

6. Rosenzweig, C. et al. Assessing agricultural risks of climate change in the 21st century in a global gridded crop model intercomparison. *Proc. Natl. Acad. Sci. USA* **111**, 3268–73 (2014).

7. Rosenzweig, C. & Parry, M. L. Potential impact of climate change on world food supply. *Nature* **367**, 133–138 (1994).

8. White, J. W., Hoogenboom, G., Kimball, B. A. & Wall, G. W. Methodologies for simulating impacts of climate change on crop production. *Field Crop Res* **124**, 357–368 (2011).

9. Farquhar, G. D., von Caemmerer, S. & Berry, J. A. A biochemical model of photosynthesis CO₂ assimilation in leaves of C3 species. *Planta* **149**, 78–90 (1980).

10. Asseng, S. et al. Uncertainty in simulating wheat yields under climate change. *Nat. Clim. Chang.* **3**, 827–832 (2013).

11. Bassu, S. et al. How do various maize crop models vary in their responses to climate change factors? *Glob. Chang. Biol.* **20**, 2301–2320 (2014).

12. Li, T. et al. Uncertainties in predicting rice yield by current crop models under a wide range of climatic conditions. *Glob. Chang. Biol.* **21**, 1328–1341 (2015).

13. Confalonieri, R. et al. A taxonomy-based approach to shed light on the babel of mathematical models for rice simulation. *Environ. Model. Softw.* **85**, 332–341 (2016).

14. Leakey, A. D. B. et al. Elevated CO₂ effects on plant carbon, nitrogen, and water relations: six important lessons from FACE. *J. Exp. Bot.* **60**, 2859–76 (2009).

15. Rosenzweig, C. et al. The Agricultural Model Intercomparison and Improvement Project (AgMIP): Protocols and pilot studies. *Agric. For. Meteorol.* **170**, 166–182 (2013).

16. Baker, J. T., Allen, L. H. & Boote, K. J. Growth and yield responses of rice to carbon dioxide concentration. *J. Agric. Sci.* **115**, 313–320 (1990).

17. Kimball, B. A. et al. Productivity and water use of wheat under free-air CO₂ enrichment. *Glob. Chang. Biol.* **1**, 429–442 (1995).

18. de Graaff, M.-A., van Groenigen, K.-J., Six, J., Hungate, B. & Kessel, Van, C. Interactions between plant growth and soil nutrient cycling under elevated CO₂: a meta-analysis. *Glob. Chang. Biol* **12**, 2077–2091 (2006).

19. Ainsworth, E. A. Rice production in a changing climate: a meta-analysis of responses to elevated carbon dioxide and elevated ozone concentration. *Glob. Chang. Biol* **14**, 1642–1650 (2008).

20. Hasegawa, T. et al. In *Improving Modeling Tools to Assess Climate Change Effects on Crop Response* (eds Hatfield, J. L. & Fleisher, D.) 45–68, https://doi.org/10.2134/advagricsystmodel7.2014.0015 (American Society of Agronomy, 2016).

21. Cai, C. et al. Responses of wheat and rice to factorial combinations of ambient and elevated CO₂ and temperature in FACE experiments. *Glob. Chang. Biol.* **22**, 856–874 (2016).

22. Kim, H.-Y. et al. Effects of free-air CO₂ enrichment and nitrogen supply on the yield of temperate paddy rice crops. *Field Crop Res* **83**, 261–270 (2003).

23. Yang, L. et al. The impact of free-air CO₂ enrichment (FACE) and N supply on yield formation of rice crops with large panicle. *Field Crop Res* **98**, 141–150 (2006).

24. Baker, J. T., Allen, L. H., Boote, K. J. & Pickering, N. B. In *Carbon Dioxide and Terrestrial Ecosystems* 265–282, https://doi.org/10.1016/B978-012505295-5/50016-1 (Elsevier, 1996).

25. Arp, W. J. Effects of source-sink relations on photosynthetic acclimation to elevated CO₂. *Plant, Cell Environ.* **14**, 869–875 (1991).

26. Martre, P. et al. Multimodel ensembles of wheat growth: many models are better than one. *Glob. Chang. Biol.* **21**, 911–925 (2015).

27. Tubiello, F. N. et al. Crop response to elevated CO₂ and world food supply. A comment on 'Food for Thought' by Long *et al.*, Science **312**:1918–1921, 2006. *Eur J Agron* **26**, 215–223 (2007).

28. Bunce, J. A. Responses of cotton and wheat photosynthesis and growth to cyclic variation in carbon dioxide concentration. *Photosynthetica* **50**, 395–400 (2012).

29. Yin, X. Improving ecophysiological simulation models to predict the impact of elevated atmospheric CO₂ concentration on crop productivity. *Ann. Bot.* **112**, 465–75 (2013).

30. O'Leary, G. J. et al. Response of wheat growth, grain yield and water use to elevated CO₂ under a Free-Air CO₂ Enrichment (FACE) experiment and modelling in a semi-arid environment. *Glob. Chang. Biol* **21**, 2670–2686 (2015).

31. Ewert, F. Modelling Plant Responses to Elevated CO₂: How Important is Leaf Area Index? *Ann. Bot.* **93**, 619–627 (2004).

32. Kim, H.-Y., Lieffering, M., Kobayashi, K., Okada, M. & Miura, S. Seasonal changes in the effects of elevated CO₂ on rice at three levels of nitrogen supply: a free air CO₂ enrichment (FACE) experiment. *Glob. Chang. Biol.* **9**, 826–837 (2003).

33. Yang, L. et al. Seasonal changes in the effects of free-air CO₂ enrichment (FACE) on dry matter production and distribution of rice (Oryza sativa L.). *Field Crop Res* **98**, 12–19 (2006).
34. Hasegawa, T., Okada, M., Inoue, M. & Shimono, H. Does elevated CO2 concentration affect lamina length of rice cultivars? In *New directions for a diverse planet, Proceedings of the 4th International Crop Science Congress* at http://www.regional.org.au/au/asa/2004/poster/2/7/1/1172_hasegawat.htm (ed. Fischer, R.) (2005).

35. Ferris, R., Sabatti, M., Miglietta, F., Mills, R. F. & Taylor, G. Leaf area is stimulated in Populus by free air CO₂ enrichment (POPFACE), through increased cell expansion and production. *Plant, Cell Environ.* **24**, 305–315 (2001).

36. Ziska, L. H. Evaluation of the growth response of six invasive species to past, present and future atmospheric carbon dioxide. *J. Exp. Bot.* **54**, 395–404 (2003).

37. Ziska, L. H. Three-year field evaluation of early and late 20th century spring wheat cultivars to projected increases in atmospheric carbon dioxide. *Field Crop Res* **108**, 54–59 (2008).

38. Ziska, L. H. & Blumenthal, D. M. Empirical selection of cultivated oat in response to rising atmospheric carbon dioxide. *Crop Sci.* **47**, 1547–1552 (2007).

39. Reich, P. B. et al. Nitrogen limitation constrains sustainability of ecosystem response to CO₂. *Nature* **440**, 922–925 (2006).

40. Terrer, C., Vicca, S., Hungate, B. A., Phillips, R. P. & Prentice, I. C. Mycorrhizal association as a primary control of the CO₂ fertilization effect. *Science* **353**, 72–74 (2016).

41. Zhang, G. et al. The effects of free-air CO₂ enrichment (FACE) on carbon and nitrogen accumulation in grains of rice (*Oryza sativa* L.). *J. Exp. Bot.* **64**, 3179–3188 (2013).

42. Prasad, P. V. V., Boote, K. J., Allen, L. H., Sheehy, J. E. & Thomas, J. M. G. Species, ecotype and cultivar differences in spikelet fertility and harvest index of rice in response to high temperature stress. *Field Crop Res* **95**, 398–411 (2006).

43. Okada, M. et al. Free-air CO₂ enrichment (FACE) using pure CO₂ injection: system description. *New Phytol.* **150**, 251–260 (2001).

44. Liu, G. et al. Rice-wheat rotational FACE platform I. System structure and control. *Chinese J. Appl. Ecol.* **13**, 1253–1258 (2002).

45. Pierce Jones, J. W., Jones, L. H., Allen, Jr & Mishoe., J. W. Dynamic Computer Control of Closed Environmental Plant Growth Chambers. Design and Verification. *Trans. ASAE* **27**, 879–888 (1984).

46. Baker, J. T., Allen, L. H. J. & Boote, K. J. *In Climate Change and Agriculture: Analysis of Potential International Impacts* (eds Rosenzweig, C., Jones, J. W. & Allen, L. H.) 31–47, https://doi.org/10.1017/CBO9781107415324.004 (1995).

47. Bouman, B. A. M. & van Laar, H. H. Description and evaluation of the rice growth model ORYZA2000 under nitrogen-limited conditions. *Agric. Syst.* **87**, 249–273 (2006).

48. Gaydon, D. S. et al. Rice in cropping systems—Modelling transitions between flooded and non-flooded soil environments. *Eur. J. Agron.* **39**, 9–24 (2012).

49. Singh, U., Ritchie, J. T. & Godwin, D. C. *A Users Guide to CERES-Rice V2.10, Simulation manua*. IIFDC-SM-4, (IFDC, 1993).

50. Fumoto, T., Kobayashi, K., Li, C., Yagi, K. & Hasegawa, T. Revising a process-based biogeochemistry model (DNDC) to simulate methane emission from rice paddy fields under various residue management and fertilizer regimes. *Glob. Chang. Biol.* **14**, 382–402 (2007).

51. Yoshida, R. et al. Adaptation of rice to climate change through a cultivar-based simulation: a possible cultivar shift in eastern Japan. *Clim. Res.* **64**, 275–290 (2015).

52. Yin, X. & van Laar, H. H. *Crop Systems Dynamics: An Ecophysiological Simulation Model for Genotype-by-Environment Interactions*. (Wageningen Academic Publishers, 2005).

53. Yoshida, H., Horie, T., Katsura, K. & Shiraiwa, T. A model explaining genotypic and environmental variation in leaf area development of rice based on biomass growth and leaf N accumulation. *Field Crop Res* **102**, 228–238 (2007).

54. Yoshida, H., Horie, T. & Shiraiwa, T. A model for explaining genotypic and environmental variation in vegetative biomass growth in rice based on observed LAI and leaf nitrogen content. *Field Crop Res* **108**, 222–230 (2008).

55. Yoshida, H. & Horie, T. A process model for explaining genotypic and environmental variation in growth and yield of rice based on measured plant N accumulation. *Field Crop Res* **113**, 227–237 (2009).

56. Yoshida, H., Horie, T., Nakazono, K., Ohno, H. & Nakagawa, H. Simulation of the effects of genotype and N availability on rice growth and yield response to an elevated atmospheric CO₂ concentration. *Field Crop Res* **124**, 433–440 (2011).

57. Hasegawa, T. & Horie, T. In *Applications of systems approaches at the field level Vol. 2* (eds Kropff, M. J. et al.) 243–257 (Kluwer Academic Publishers, 1997).

58. Aggarwal, P. K., Kalra, N., Chander, S. & Pathak, H. InfoCrop: A dynamic simulation model for the assessment of crop yields, losses due to pests, and environmental impact of agro-ecosystems in tropical environments. I. Model description. *Agric. Syst.* **89**, 1–25 (2006).

59. Masutomi, Y., Ono, K., Mano, M., Maruyama, A. & Miyata, A. A land surface model combined with a crop growth model for paddy rice (MATCRO-Rice v. 1) – Part 1: Model description. *Geosci. Model Dev.* **9**, 4133–4154 (2016).

60. Masutomi, Y. et al. A land surface model combined with a crop growth model for paddy rice (MATCRO-Rice v. 1) – Part 2: Model validation. *Geosci. Model Dev.* **9**, 4155–4167 (2016).

61. Takata, K., Emori, S. & Watanabe, T. Development of the minimal advanced treatments of surface interaction and runoff. *Glob. Planet. Change* **38**, 209–222 (2003).

62. Tao, F., Zhang, Z., Liu, J. & Yokozawa, M. Modelling the impacts of weather and climate variability on crop productivity over a large area: A new super-ensemble-based probabilistic projection. *Agric. For. Meteorol.* **149**, 1266–1278 (2009).

63. Tang, L. et al. RiceGrow: A rice growth and productivity model. *NJAS - Wageningen J. Life Sci.* **57**, 83–92 (2009).

64. Wu, G. W. & Wilson, L. T. Parameterization, verification, and validation of a physiologically complex age-structured rice simulation model. *Agric. Syst.* **56**, 483–511 (1998).

65. Dingkuhn, M., Pasco, R. & Soulie, J. SAMARA V2.2 description. at http://umr-agap.cirad.fr/en/content/download/4886/38078/version/3/file/SAMARA+V22+description.pdf (2013).

66. Horie, T. A model for evaluating climatic productivity and water balance of irrigated rice and its application to Southeast Asia. *Southeast Asian Stud.* **25**, 62–74 (1987).

67. Horie, T., Nakagawa, H., Centeno, H. G. S. & Kropff, M. J. *In Modeling the impact of climate change on rice production in Asia* (eds Matthews, R. B., Kropff, M. J., Bachelet, D. & van Laar, H. H.) 51–56 (CAB Interanational in association with International Rice Research Institute, 1995).

68. Brisson, N. et al. STICS: a generic model for the simulation of crops and their water and nitrogen balances. I. Theory and parameterization applied to wheat and corn. *Agronomie* **18**, 311–346 (1998).

69. Brisson, N. et al. An overview of the crop model STICS. *Eur. J. Agron.* **18**, 309–332 (2003).

70. Confalonieri, R., Rosenmund, A. S. & Baruth, B. An improved model to simulate rice yield. *Agron. Sustain. Dev.* **29**, 463–474 (2009).

71. Stockle, C. O., Dyke, P. T., Williams, J. R., Jones, C. A. & Rosenberg, N. J. A method for estimating the direct and climatic effects of rising atmospheric carbon dioxide on growth and yield of crops: Part II-Sensitivity analysis at three sites in the Midwestern USA. *Agric. Syst.* **38**, 239–256 (1992).

## Acknowledgements

The authors are grateful to the AgMIP project, which has made the collaboration behind this study possible. YZ's work was supported by the National High-Tech Research and Development Program of China (2013AA100404).

FT is supported by the National Science Foundation of China (31561143003) and the Academy of Finland,
PLUMES project (277403). TH's participation was supported by MAFF and the Global Environment Research (S-10-2) of MOE, Japan. TL's work was supported by the project of Global Foresight and Strategic Plan. RC was partly supported by EC- FP7 under Grant Agreement No. 613817 (MODEXTREME). Nanjing Agricultural University, International Rice Research Institute, Ministry of Agriculture, Forestry and Fisheries and Ministry of the Environment of Japan, and Cassandra lab. of the University of Milan supported the annual research meetings of the AgMIP Rice Team. USDA is an equal opportunity provider and employer.

**Author Contributions**

T.H., T.Li. and X.Y. equally contributed to this work. T.Li., T.H., X.Y., K.B., Y.Z. designed the research. T.H., X.Y., T.Li. wrote the manuscript. T.H., T.Li., X.Y., Y.Z., S.B., S.B., R.C., J.F., T.F., D.G., S.N.K., T.La., M.M., Y.M., H.N., P.O., F.R., U.S., L.T., F.T., L.T.W., Y.Y., H.Y., Z.Z. conducted simulations. T.H., T.Li., H.W., D.W. analyzed the simulation output. T.H., J.B., W.Y., J.Z., L.Y. provided the experimental details and data, and J.B., W.Y., J.Z., L.Y. approved the analytical results. All authors have read and approved the manuscript.

**Additional Information**

**Supplementary information** accompanies this paper at https://doi.org/10.1038/s41598-017-13582-y.

**Competing Interests:** The authors declare that they have no competing interests.

**Publisher's note:** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_11.png" /></figure>

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.

© The Author(s) 2017