Revision of estimates of climate change impacts on rice yield and quality in Japan by considering the combined effects of temperature and CO 2 concentration

Rice is the most important cereal crop in Japan, and therefore the impact of projected climate change on its production and quality has been assessed using rice growth models accounting for the effects of rising temperature and atmospheric CO 2 concentration ( [CO 2 ] ) on important growth processes. Recent experimental studies, however, have shown some negative effects of interactions between [CO 2 ] and temperature on yield and quality of rice which were not accounted for by previous impact assessments. This study examined the importance of [CO 2 ] × temperature interactions in the nationwide impacts of climate change on grain yield and quality of rice in Japan by 2100. We introduced new functions accounting for the effects of interactions on yield. Then we adopted the acceleration by elevated [CO 2 ] in the estimation of the occurrence of chalky grains, an indicator of appearance quality of rice. We applied the modified model to Japan at a spatial resolution of 1 km using 10 climate scenarios ( 5 Global Circulation Models × 2 representative concentration pathways [RCPs] ) from 1981 to 2100. The effects of the newly introduced negative effects of [CO 2 ] × temperature were evaluated by comparing simulations with and without the interaction in each scenario. Nationwide production was estimated to decrease by up to 28 % and the percentage of white chalky grains to increase up to 16 % relative to the previous assessment results, especially in RCP8.5, in which larger increases were projected in both temperature and [CO 2 ]. The result suggests that the positive effect of elevated [CO 2 ], which had been expected to offset the negative effect of increased temperature on rice productivity, may be limited in the future, and rice quality degradation may be more severe than predicted previously.


Introduction
The negative impacts of projected climate change on agriculture are considered unavoidable e.g., IPCC, 2014 , and the introduction of adaptation measures based on impact assessments is critical for national food security e.g., FAO, 2008 . In Japan, the importance of such impact assessments is widely understood, as evidenced by the Climate Change Adaptation Act of 2018 Government of Japan, 2018 , which requires climate change adaptation plans to be updated every 5 years or so following new impact assessments. Therefore, impact assessment models need to be improved continuously by incorporating new scientific knowledge in order to reduce uncertainty in prediction and improve the reliability of results.
Rice is the most important cereal crop in Japan, and thus many studies have evaluated the impact of projected climate change on rice production and quality using rice growth models Horie et al., 1995;Seino, 1997;Yonemura et al., 1998;Hayashi et al., 2001;Iizumi et al., 2006Iizumi et al., , 2011Yokozawa et al., 2009 . Following those studies, our group conducted a comprehensive impact assessment of Japanese rice production and quality degradation risk based on multiple projected climate scenarios by using a process-based rice growth model Ishigooka et al., 2017 . The predicted rice yield would increase almost continuously within the 21st century, except under climate scenarios with an extremely large temperature increase, because, in many cases, the yield increase effect of elevated atmospheric CO 2 concentration [ CO 2 ] was greater than the yield decrease effect of high temperature shortening the growing period and increasing spikelet sterility. To predict rice quality, we introduced the HD_m26 heat stress index, calculated as the heat degree-days above 26 C during 20 days after the heading date, and related to the risk of rice quality degradation due to high temperature Ishigooka et al., 2011 . Since HD_m26 depends only on temperature, it increased greatly with the projected temperature increase Ishigooka et al., 2017 . Environmental control experiments such as Free-Air CO 2 Enrichment "FACE" experiments e.g., Hasegawa et al., 2013 ] is increased around target crops in an open field e.g., Kobayashi 2001 , revealed that the yield increase effect of elevated [ CO 2 ] is suppressed under high temperatures; that is, the positive effect of elevated [ CO 2 ] on yield may not eventuate Hasegawa et al., 2016 . Therefore, our predicted change in rice total production Ishigooka et al., 2017 may be overestimated. Studies of the combined effects of high temperature and elevated [ CO 2 ] on rice quality in FACE experiments in Ibaraki Prefecture in Japan Usui et al., 2014and in Jiangsu Province in China Yang et al., 2007 showed that elevated [ CO 2 ] significantly increased the emergence of chalky grains, which reduce the appearance quality of rice, even at the same air temperature. Therefore, the risk of rice quality degradation due to temperature in our previous study Ishigooka et al., 2017 may be exacerbated by the effects of elevated [ CO 2 ] .
The objective of this study was to present an advanced impact assessment of Japanese rice production and quality under climate change using an improved model based on these new findings of the combined effect of elevated [ CO 2 ] and high temperature on rice production and quality.
To evaluate rice quality practically, we adopted the percentage of chalky grains "PCG" as an indicator of rice quality associated with high temperature, instead of the HD_m26 heat stress index used in Ishigooka et al. 2017 , enabling us to analyze the effect of elevated [ CO 2 ] on grain quality by using the findings from the above FACE experiments. PCG has been used as an indicator of quality in analyses of the impact of past climate change on quality e.g., Morita et al., 2016;Takimoto et al., 2019 . Most models developed for estimating PCG Wakiyama et al., 2010;Masutomi et al., 2015Masutomi et al., , 2019Yoshida et al., 2016;Takimoto et al., 2019 were developed for specific cultivars and would thus not be suitable for a nationwide estimation. Nishimori et al. 2020 analyzed the relationships of various indices of rice quality in 15 major cultivars currently grown in Japan and covered in the Hasegawa-Horie H/H model Hasegawa and Horie, 1997 with meteorological indices at a large number of survey points around Japan to avoid regional bias, and proposed a statistical model for estimating PCG from HD_m26. Here, we used that model to prepare a scheme to estimate PCG by considering the combined effect of high temperature and elevated [ CO 2 ] .
Our calculations were based on a third-order mesh ~1-km ~1-km instead of the second-order mesh ~10-km ~10-km used in our previous study, for impact assessment at the regional scale in a complex geographical region Ishigooka et al., 2020 .

Model 2.1.1 Rice growth model for yield and phenology
We estimated rice yield and phenology with the H/H model Hasegawa and Horie, 1997;Fukui et al., 2015;Yoshida et al., 2015 which had been validated well for nationwide large-scale simulation as described by Ishigooka et al. 2017 . The H/H model has three major components: phenological development, biomass production, and yield formation. In the biomass production component, the effect of increasing yield caused by changes in photosynthetic rate due to the increase in [ CO 2 ] is presented by the combination of the Ball-Berry stomatal conductance model Collatz et al., 1991 with the Farquhar-von Caemmerer-Berry photosynthesis model Farquhar et al., 1980, as described in detail in Yoshida et al. 2015 . Since the effect of temperature on the yield increase due to elevated [ CO 2 ] is not considered, yield tends to be large at high [ CO 2 ] regardless of temperature Ishigooka et al., 2017 . Rice FACE experiments conducted in Shizukuishi, Iwate and Tsukubamirai, Ibaraki in Japan using a standard japonica-type cultivar Akitakomachi over 11 growing seasons have shown that the positive effects of elevated [ CO 2 ] on grain yield decreased with air temperature Hasegawa et al., 2016 . Here, we tested how the negative effects of increasing temperatures on grain yield enhancement by elevated [ CO 2 ] influence the national-scale yield projection by incorporating the empirical relation in Fig. 1a into the original H/H model. As shown in Fig. 1a, the relationship between the yield increase with elevated [ CO 2 ] and mean temperature during 30 days after the heading stage can be approximated by the following linear function within the range from ~20 C to ~30 C : where ∆Y is the yield change with elevated [ CO 2 ] here, +200 ppm above ambient [ CO 2 ] ≈ 380 ppm , a is the coefficient related to the temperature dependence of the CO 2 fertilization effect = -2.1, Table 1 , and T m30 is the mean temperature during 30 days after heading. Thus, the yield increase is maximum ~+20 at low temperature ~20 C , and decreases by 2.1 per 1-C increase to nil at high temperature ~30 C . Using this relationship, we introduced a reduction factor, RF, to modify the yield estimated in our previous study Ishigooka et al., 2017 : where a' is the temperature dependence of the CO 2 fertilization effect scaled to that at T m30 = 20 C Table 1 . Then, the modified yield, Y mod , can be given as: where Y o is the estimated yield when base [ CO 2 ] = 380 ppm, and Y c is the estimated yield under [ CO 2 ] in the target year, which is given by Representative Concentration Pathways RCP; described in 2.2 . To account for the uncertainty in the key parameter, we changed a' between lower and upper 75 confidence limits Table 1 to derive the range of yield estimates. Hereinafter, the improvement in the rice yield estimation described above is referred to as the combined effect of high temperature and elevated [ CO 2 ] "HTEC" .

Estimation of rice quality
We used PCG calculated by the statistical model developed by Nishimori et al. 2019 as an indicator of rice quality. In estimating PCG Nishimori et al., 2020 , the regression formulas corresponding to each percentile 95th, 75th, 50th, 25th, 5th are obtained by the polynomial quantile regression method to represent the uncertainties in the prediction of chalky grains. The mean value of PCG was calculated by weighted averaging of the values calculated for each percentile with the probability of each percentile.
In addition to high temperature during the reproductive period, elevated [ CO 2 ] increases PCG as found in the rice FACE studies in Japan and China Yang et al., 2007;Usui et al., 2014;Usui et al., 2016 . To account for the effect of elevated [ CO 2 ] on PCG, we compiled PCG data from three published studies using the FACE experimental facility at Tsukubamirai Usui et al., 2014;Usui et al., 2016;Hasegawa et al., 2019 . Collectively, we collected 33 datasets including five growing seasons,13 cultivars differing in heat tolerance, three levels of nitrogen fertilization, and two levels of water temperatures to derive the PCG dependence on [ CO 2 ] . PCG under ambient [ CO 2 ] ranged largely from 2 -52 depending on years, cultivars, nitrogen or water temperatures, but corresponding PCG under elevated [ CO 2 ] was consistently greater by on average 49 as indicated by the highly significant slope b in Fig. 1b.
Thus, PCG under elevated [ CO 2 ] , PCG mod , can be estimated as: where ∆CO 2 is the increase in [ CO 2 ] in the target year above the baseline [ CO 2 ] of 380 ppm, b is the enhancement ratio of PGC in response to ∆CO 2 = 200 ppm Table 1 , and PCG o is the percentage of chalky grain calculated by the HD_m26 heat-stress index upper limit was set to 70 C · days Nishimori et al., 2020 . Note that the upper limit of ∆CO 2 was set to 200 ppm, as in FACE experiments Hasegawa et al., 2013 , to avoid large uncertainties due to extrapolation. We also changed the b value between the lower and upper 75 confidence limits Table 1 to present uncertainty in the appearance quality projection.

Data
The H/H model requires three types of input data: daily meteorological data, yearly [ CO 2 ] , and cultivation management practices such as transplanting date and cultivar. We implemented it with a third-order mesh composed of nearly 380 000 grid cells with ~1-km ~1-km spatial resolution , and thus all data were transformed to match this framework. The study area is illustrated in Fig. 2.
The H/H model needs daily mean, maximum, and minimum temperatures, downward shortwave radiation, relative humidity, and wind speed as inputs. All inputs were taken from projected future climate scenarios in Nishimori et al. 2019 . This dataset was developed for use in agricultural impact assessments and covers almost the whole land area of Japan at third-order resolution ~1-km ~1-km , coordinated by the Japanese Geodetic Datum 2000 JGD2000 with GRS80 ellipsoid . It comprises 10 climate scenarios, composed of 5 GCMs MIROC5, MRI-CGCM3, GFDL-CM3, HadGEM2-ES, CSIRO-Mk3-6-0 2 RCPs RCP2.6 and RCP8.5 , obtained originally from the fifth phase of the Coupled Model Intercomparison Project Taylor et al., 2012 . Biases in each original GCM were corrected by linear  Table 1. The cultivar used in each cell was the cultivar with the largest area in the prefecture to which the grid cell belongs, as obtained from the statistical yearbooks available on the MAFF website http://www.maff.go.jp/j/tokei/kouhyou/syokuryo_nenkan/; last accessed 25 February 2021 . If that cultivar was not included in the 15 cultivars used in the H/H model, the second one was selected.
The N fertilizer per area of paddy fields in each cell was drawn from statistical data on the consumption of fertilizers from the MAFF website http://www.maff.go.jp/ j/ tokei/ kouhyou/ noukei/seisanhi_nousan/; last accessed 25 February 2021 and the Handy Directory of Fertilizer 2009 MAFF, 2010 , as described by Ishigooka et al. 2017 .
The transplanting date in each cell was obtained from statistical data collated by sub-administrative regions by using the subregion codes, as described by Ishigooka et al. 2011 .

Simulation conditions
The model was implemented to run from 1981 to 2100. Rice yield and quality were calculated in cells where paddy fields covered ≥1 of the total land area within the cell. Of the total 377 449 cells, this criterion selected 152 927 cells. To evaluate the effect of considering the combined effects of high temperature and elevated [ CO 2 ] , we assumed the cultivation management conditions such as cultivars, N input, and transplanting date in each cell and the distribution and areas of paddy fields as in 2006 for the entire simulation period.

Rice yield
To reveal the effect on the yield calculation by incorporating the combined effect of high temperature and elevated [ CO 2 ] "HTEC" into the previous model, we averaged yields calculated by the previous without HTEC model and the modified with HTEC model for the periods 1981-2000 "baseline period" , 2031-2050 "mid period" , and 2081-2100 "late period" in each cell. Table 2 shows the 20-year average nationwide total rice production estimated by both models in the mid and late periods under each GCM RCP, all expressed as relative to that of the baseline period = 100 . In all cases, the estimated production was remarkably lower with HTEC than without HTEC. The differences were small under RCP2.6, ranging from -4.2 to -5.5 average, -4.9 in the mid period and from -3.7 to -4.3 average, -4.0 in the late period, but larger under RCP8.5, ranging from -7.9 to -9.7 average, -8.9 in the mid period and from -17.4 to -28.2 average, -23.8 in the late period.

Fig. 2.
Study area in Japan. Thick black lines, regional boundaries; thin gray lines, prefectural boundaries. The black shadows in the Japanese land area represent the geographical feature the darker the color, the higher the elevation .   Figure 3 shows the temporal changes in the 20-year average of relative production estimated by both models from 1981 to 2100, presenting the range estimated by each GCM and the average. Under RCP2.6, production showed slight continuous increases throughout the period both without HTEC and with HTEC, smaller with HTEC Fig. 3a . The differences between production without HTEC and with HTEC increased until 2021-2040 and then stabilized. The ranges of production among GCMs increased until 2041-2060 and then decreased in both models. Under RCP8.5, production estimated without HTEC increased greatly until 2041-2060 and then decreased, while that with HTEC increased slightly until 2021-2040 and then decreased greatly Fig. 3b . From 2061-2080, most estimates of relative production were >100 without HTEC, but most were <100 with HTEC. The differences became larger with time. The ranges among GCMs increased continuously throughout the period in both models. The variation in the projected impacts on production due to the uncertainty in the key parameter for HTEC a' in Eq. 2 also grew as the years progressed, but the range was consistently smaller than those among GCMs or emission scenarios Fig. 3 .
Spatial distribution maps of 20-year averages of relative yield estimated without HTEC and with HTEC in all 10 climate scenarios in the mid and late periods are shown in the Supplemental Information Supp. Figs. S1-S4 . As an example, maps of relative yields estimated by both models forced by MIROC-5 + RCP8.5 in the mid and late periods are shown in Fig. 4. In the mid period Fig. 4a, b , yields were estimated to increase relative to the baseline period without HTEC in most areas except for limited areas in the Tokai, Kinki, and Kyushu regions, but to decrease with HTEC in the plains in the eastern to western regions. In the late period Fig. 4c, d , as in the mid period, yields were estimated to increase without HTEC in most areas and to decrease in some of the plains from the eastern to western regions, but to decrease with HTEC in most areas except Hokkaido and the eastern side of the Tohoku region and in inland high-altitude areas from the eastern to western regions. The yields with HTEC were significantly lower than those without HTEC.

Rice quality
Yearly PCGs were estimated with HD_m26, which is based on the heading date in the H/H model, and are summarized as 20-year averages in each climate scenario Table 3 . In the baseline period, the nationwide average of PCG was 4 of the total grain. Projected PGCs differed largely between RCPs Fig. 5, Table 3 . Under RCP8.5, PCG increased almost exponentially in the first half of this century and then linearly toward the end of the century, whereas the increase was modest and plateaued by the mid-century under RCP2.6 Fig. 5 . In the mid and late periods, comparisons between the results estimated only by temperature "TO" and those estimated by considering the effect of elevated [ CO 2 ] "T&C" show that the PCGs were large under climate scenarios with large temperature increases Table 2 , and increased with [ CO 2 ] . Differences between T&C and TO were small under RCP2.6, ranging from 0.7 to 2.5 average, 1.5 in the mid period and from 0.7 to 2.2 average, 1.3 in the late period, and larger under RCP 8.5, ranging from 1.8 to 5.3 average, 3.3 in the mid period and from 8.9 to 15.5 average, 12.5 in the late period. The variation of PCG due to changes in the key parameter b Eq. 4 was less than 3 even in the late century Fig. 5c , which was much smaller than the variation among GCMs and RCPs.
Spatial distribution maps of 20-year averages of PCGs estimated with TO and with T&C in all 10 climate scenarios are shown in the Supplemental Information Supp. Figs. S5-S8 . As an example, Fig. 6 shows the spatial distribution of PCGs forced by MIROC-5 + RCP8.5 in the mid and late periods. Because the methods used to calculate the PCG are based on the heat stress index HD_m26 Nishimori et al., 2020 , the spatial pattern of the PCG distribution corresponds to that of temperature in general. Areas with large values of PCG are shown in the plains or basins from the eastern to central regions. In the western regions, the estimated PCG is low despite the high temperatures, as the heading date is later than in the other regions.
In the mid period Fig. 6a, b , the PCG was estimated with TO to be <20 in most of the country, except for part of the plains in the central region ~8 of the country; Fig. 6a , but with T&C increased by 29 , expanding the area with PCG > 20 to 20 of the country and the area with PCG > 30 into the plains from the eastern to western regions Fig. 6b . In the late period Fig. 6c, d , the PCGs were estimated with T&C to increase by 50 Fig. 6d relative to those estimated with TO Fig. 6c . The area with PCG > 20 increased from 71 to 88 of the country, that with PCG > 30 from 12 to 71 , and that with PCG > 40 from 0 to 32 .

Discussion
The inclusion of HTEC in the estimation of yield showed that estimated rice yields could be much lower than estimates without HTEC. In the late period, most GCMs under RCP8.5 estimated the total rice production to increase without HTEC but to decrease with HTEC Fig. 3b . Estimated increases became decreases in most of the country except for parts of the north Fig. 4c, d . These results suggest that this improvement in methods may lead to opposite results in assessing the effect of climate change on rice yields In the estimation of rice yield by a process-based rice growth model under the warming conditions associated with elevated [ CO 2 ] projected in the future, the change in yield is determined by the balance among effects causing a yield decrease, such as the shortening of the growing period Yoshida, 1981;Krishnan et al., 2011 and increase in spikelet sterility due to high temperatures Kim et al., 1996;Matsui et al., 1997 , and effects causing a yield increase, such as the CO 2 fertilization effect Kimball et al., 2002;Hasegawa et al., 2013 and the reduction of cold damage risk in cooler regions Hayashi et al., 2001;Yokozawa et al., 2009 . Among these factors, the CO 2 fertilization effect has played a significant role in the prediction of rice yield under the projected increased temperature and [ CO 2 ] , as it is expected to offset yield decreases due to high temperature across the country.
In our previous study Ishigooka et al., 2017 , the total Japanese rice production was estimated to increase slightly throughout the 21st century in most climate change scenarios, but to decrease significantly only in areas where the large increase in spikelet sterility due to higher temperatures exceeded the CO 2 fertilization effect under climate scenarios with particularly large temperature increases such as GFDL-CM3 +6.9 C and HadGEM2-ES in RCP8.5 +5.7 C from baseline to late periods . This study incorporated the empirical relationship between the CO 2 fertilization effect on yield and air temperature obtained from the rice FACE experiments Fig. 1a in the H/H model. The uncertainty of the key parameter in the relationship a' as expressed by the lower and upper 75 confidence limit was much smaller than that associated with GCMs and emission scenarios Fig. 3 . These results suggest that previous predictions of rice production in Japan under the projected climate conditions may be too optimistic and thus might be revised downward in subsequent impact assessments. The negative effect of higher temperatures on PCG is well documented and has been quantitatively evaluated e.g., Masutomi et al., 2019 . However, the effects of rising [ CO 2 ] on PCG and its interaction with temperature have not been considered previously.
The 20-year average of PCG in the mid period estimated by the TO function proposed by Nishimori et al. 2019 agreed well with those in the same period reported by Masutomi et al. 2019 . The mechanism by which elevated [ CO 2 ] exacerbates the heat-induced increase in PCG is not fully understood, but it could be related to increased canopy temperatures or to reduced grain protein content, both of which are known to increase PCG Usui et al., 2016 . Since the increase in PCG due to elevated [ CO 2 ] is assumed to be proportional to the increase in [ CO 2 ] from baseline, the increase in PCGs with T&C was more remarkable in RCP8.5, with a higher [ CO 2 ] , than in RCP2.6 Table 2 . According to the relationship between PCG and the proportion of each grade of rice quality shown in Fig. 1 of Nishimori et al. 2020 , the percentage of first-grade rice is <40 when PCG > 10 , <20 when PCG > 20 , and nil when PCG > 30 . Using this relationship, the PCGs calculated here can be evaluated in conjunction with the proportion of each grade to enable a more practical assessment of the effect on the quality of Japanese rice.
To the best of our knowledge, this study is the first to consider the combined effects of high temperature and elevated [ CO 2 ] on rice yield and quality prediction on a nationwide scale. However, there are some uncertainties in the combined effects on yield and quality, since the relationships were obtained from experiments with limited sites and cultivars. The temperature dependence of the yield enhancement effect due to elevated [ CO 2 ] was obtained by an analysis using just one cultivar Akitakomachi . Some cultivars retain positive effects of elevated [ CO 2 ] on yield under high temperature during the grain-filling period Hasegawa et al., 2013 , so that our yield estimates based on the temperature dependence of the CO 2 fertilization of Akitakomachi, a typical japonica cultivar, could be highly sensitive to high temperature. Although this finding is certainly relevant e.g., Nakagawa and Horie 2000 , it is necessary to test the relationship with many other cultivars in order to validate the results of this study. On the other hand, we showed that the effect of elevated [ CO 2 ] on PCG was consistent across cultivars differing in heat tolerance. A ~50 enhancement in PGC with +200 ppm is highly significant with small uncertainty Fig. 5 . In the FACE experiments, however, the elevated [ CO 2 ] had been set to 200 ppm above ambient, and thus estimation will be less reliable at higher values of [ CO 2 ] projected after 2057 under RCP8.5. To avoid calculating unrealistic values by extrapolation, we set the upper limit of the increase in PCG due to elevated [ CO 2 ] at 50 . Therefore, the estimated PCGs in the late period under RCP8.5, in which [ CO 2 ] > 580 ppm Table 2; Fig. 6d , are possibly underestimated. Although there are still considerable uncertainties in the methodologies in this study, it is meaningful that we have shown that the modified impact assessment model incorporating new knowledge obtained from recent experiments could predict significantly different impacts from those of previous studies. To reduce uncertainties in the impact assessment of projected climate change on crop productivity, we are just continuing further studies to elucidate the physiological responses of crops under high temperature and elevated [ CO 2 ] . Also, this suggests the importance of adaptive management in climate change adaptation planning by establishing alternative measures to