Modeling the rice phenology and production in China with SIMRIW: sensitivity analysis and parameter estimation

Abstract

Crop models are robust tools for simulating the impact of climate change on rice development and production, but are usually designed for specific stations and varieties. This study focuses on a more adaptable model called Simulation Model for Rice-Weather Relations (SIMRIW). The model was calibrated and validated in major rice production regions over China, and the parameters that most affect the model’s output were determined in sensitivity analyses. These sensitive parameters were estimated in different ecological zones. The simulated results of single and double rice cropping systems in different ecological zones were then compared. The accuracy of SIMRIW was found to depend on a few crucial parameters. Using optimized parameters, SIMRIW properly simulated the rice phenology and yield in single and double cropping systems in different ecological zones. Some of the parameters were largely dependent on ecological zone and rice type, and may reflect the different climate conditions and rice varieties among ecological zones.

Keywords

rice / phenology / parameter optimization / SIMRIW / simulation

Introduction

Rice is the staple food of Asians, Africans and Latin Americans (Yoshida, 1981; Peng et al., 1995). Throughout the past several decades, a warming trend has been documented around the world. This global warming trend has presented special obstacles and challenges to rice production (Tao et al., 2003; Peng et al., 2004; Tao et al., 2006). The effect of climate change on rice production, which has raised considerable concern among scientists (Kropff et al., 1993; Horie et al., 1997; Matthews et al., 1997; Hayashi and Jung, 2000; Xiong et al., 2001; Lin et al., 2005), is frequently investigated and predicted by crop models. Since their introduction in the 1960s (Brouwer and De Wit, 1969; Duncan, 1971), crop models have benefitted from advances in agricultural science and computer technology. The Simulation Model for Rice-Weather Relations (SIMRIW) is a simplified process model that simulates the growth and yield of irrigated rice under different weather conditions. The model rationally simplifies the underlying physiological and physical processes of rice crop growth, and has been used to assess the impact of climate change on paddy rice yield (Horie et al., 1995).

The large number of parameters in typical models may introduce uncertainty in the simulated results (Anderson, 2010). To ensure consistency between model predictions and corresponding observations, we must calibrate the parameter values in crop models. Although several rice models have been developed and applied (Spitters, 1986; Gao et al., 1992; Yin and Kropff, 1996; Bouman and Van Laar, 2006; Wopereis et al., 2009), these have usually been limited to specific station and varieties, and do not account for the varying parameter values among ecological zones.

In this paper, we conduct a sensitivity analysis of the SIMRIW parameters and select the parameters that most greatly affect the output of the model. These sensitive parameters are optimized in different ecological zones. The simulated results of single and double rice cropping systems are then compared in different ecological zones.

Materials and methods

Description of model

SIMRIW is a simplified process model for simulating the growth and yield of irrigated rice under different weather conditions. This model assumes that the grain yield (Y_G, g·m^‒2) constitutes a specific proportion of the total dry matter production (W_total, g·m^‒2) of rice:

(1)Y_G=h{W}_{\mathrm{t}\mathrm{o}\mathrm{t}\mathrm{a}\mathrm{l}},

where h is the harvest index.

Phenological development of the crop

The rice development rate in SIMRIW is described by the development index (DVI). DVI is defined as 1.0 at heading, and 2.0 at maturity, i.e., 1<DVI<2. The DVI at day t after transplanting is calculated by summing the rice developmental rate (DVR) over time:

(2)DV{I}_t={\displaystyle \sum _{i=\mathrm{0}}^{i=t}DV{R}_i},

where DVR_i is the development rate at day i (day^‒1). The daily development rate is given by the daily temperature and day length, and depends on the crop stage as follows:

(3)DVR=\mathrm{1}/G_v\left\{\mathrm{1.0}+\exp \left[-A_T\left(T-T_h\right)\right]\right\}\mathrm{for}DVI\le DV{I}^{*},

DVR=\left\{\mathrm{1}-\exp \left[B\left(L-L_c\right)\right]\right\}/G_v\left\{\mathrm{1.0}+\exp \left[-A\left(T-T_h\right)\right]\right\}

(4)\mathrm{for}DVI>DV{I}^{*}\mathrm{and}L\le L_c,

(5)

(6)

where G_v is the minimum number of days required for heading (day), A_T is the sensitivity of the developmental rate to air temperature (dimensionless), T is the daily mean air temperature (°C), T_h is the air temperature at which DVR equals half of the maximum rate at the optimum temperature (°C), DVI* is the DVI at which the crop becomes sensitive to the photoperiod (day^‒1), L is the day length (h), L_c is the critical day length (h), and K_r and T_cr are empirical constants. G_r is the minimum number of days in the grain-filling period.

Dry matter production

In SIMRIW, the amount of radiation absorbed by the canopy (S_s) is a function of the leaf area index (LAI) (F) (m²·m^‒2).

(7)S_s=S_o\left\{\mathrm{1}-r-\left(\mathrm{1}-r_o\right)\exp \left[-\left(\mathrm{1}-m\right)kF\right]\right\},

where S_o is the daily incident solar radiation (MJ·m^‒2·day^‒1),r is the reflectance of the canopy, r_o is the canopy of bare soil, k is the extinction coefficient of the canopy to daily short-wave radiation, and m is the scattering coefficient. The canopy reflectance (r) is given by the following equation:

(8)r=r_f-\left(r_f-r_o\right)\exp \left(-{\displaystyle \frac{F}{\mathrm{2}}}\right),

where, r_f is the reflectance of a surface completely covered by vegetation. We assigned the following values to the parameters in the above functions: k = 0.6, m = 0.25, r_f = 0.22, and r_o = 0.1.

The conversion efficiency of absorbed short-wave radiation to rice crop biomass C_s is influenced by the atmospheric CO₂ concentration:

(9)C_s=C_o\left[\mathrm{1}+{\displaystyle \frac{R_m\left(C_a-\mathrm{330}\right)}{\left(C_a-\mathrm{330}\right)+K_c}}\right],

where C_o is the radiation conversion efficiency at 330 ppm<FootNote>

1 ppm=1×10^-6

</FootNote> CO₂ (g·MJ^‒1). R_m is the asymptotic limit of the relative response to CO₂, and K_c is an empirical constant (ppm).

The relative growth rate of leaf area index (F) is related to the daily mean temperature (T) in the period before heading:

(10){\displaystyle \frac{\mathrm{1}}{F}}\times {\displaystyle \frac{\mathrm{d}F}{\mathrm{d}t}=A\left\{\mathrm{1}-\exp \left[-K_r\left(T-T_{cf}\right)\right]\right\}\left[\mathrm{1}-{\left(F/F_{as}\right)}^h\right],}

where A is the maximum relative growth rate of LAI (m²·m^‒2), obtained under optimum conditions of non-limiting temperature, solar radiation, nutrients, pests and diseases, T_cf is the minimum temperature for LAI growth (°C), F_as is the asymptotic leaf area index when the temperature is non-limiting (m²·m^‒2), and K_r and h are empirical constants.

Yield formation

The harvest index h is the lower value of the low-temperature stress h_c and high-temperature stress h_h:

(11)h=\min \left(h{}_c,h_h\right).

The low-temperature harvest index is formulated as:

(12)h_c=h_m\mathrm{\gamma }\left(\mathrm{1}-{\mathrm{\gamma }}_c\right),

(13)\mathrm{\gamma }=\mathrm{1}-\exp \left[-K_h\left(DV{I}_t-\mathrm{1.22}\right)\right],

where h_m is the maximum harvest index in a given prefecture, obtained under optimum climatic conditions and cultivation practices. K_h is an empirical constant and γ_c is a sterile spikelet:

(14){\mathrm{\gamma }}_c={\mathrm{\gamma }}_o-K_q{Q}^{C_{\mathrm{c}\mathrm{o}\mathrm{o}\mathrm{l}}}.

In Eq. (14), γ_o and K_q are empirical constants, and C_cool is the curvature factor of the low-temperature spikelet sterility. Q is calculated as follows:

(15)Q={\displaystyle \sum \left(T^{*}-T\right)},

where T* is the base temperature, and the summation is performed over the period \mathrm{0.75}\le DVI\le \mathrm{1.20}.

The high-temperature harvest index is formulated as:

(16)h_h=h_m\left(\mathrm{1}-\mathrm{0.95}{\mathrm{\gamma }}_h\right),

(17)\{\begin{array}{l}\begin{array}{lll}{\mathrm{\gamma }}_h={\left[\left({\displaystyle \frac{T_{\max }-T_b}{T_o-T_b}}\right){\left({\displaystyle \frac{T_c-T_{\max }}{T_c-T_o}}\right)}^{\frac{T_c-T_o}{T_o-T_b}}\right]}^{C_{\mathrm{h}\mathrm{o}\mathrm{t}}}\hfill & \mathrm{f}\mathrm{o}\mathrm{r}\hfill & T_{\max }>T_o,\hfill \end{array}\hfill \\ \\ \begin{array}{lll}{\mathrm{\gamma }}_h=\mathrm{0}\hfill & \mathrm{f}\mathrm{o}\mathrm{r}\hfill & T_{\max }\le T_o.\hfill \end{array}\hfill \end{array}

Rice phenology data and climate data

Data on rice phenology and yield were obtained from the agro-meteorological experimental stations of the China Meteorological Administration (CMA) from 1981 to 2009. Early-rice, late-rice and single-rice data were collected from 120, 117, and 160 stations, respectively.

Following Mei et al. (1988), rice cultivation areas in China were classified into six ecological zones according to: northeastern China (zone I), North China Plain (zone II), north of the Yangtze River (zone III), south of the Yangtze River (zone IV), southwestern China (zone V) and South Costal Areas (zone VI).

The mean, maximum, and minimum temperatures from 1981 to 2009 were also obtained from the CMA. Day length is a function of latitude and day of year (Spitters, 1986).

Solar radiation (Rad) was estimated from sunshine duration observations and the Angstrom–Prescott equation (Angström, 1924; Prescott, 1940):

(18)R_s=\left(a_s+b_s\cdot {\displaystyle \frac{n}{N}}\right)R_a,

where R_s denotes solar or shortwave radiation (MJ·m^-2·day^-1), n is the actual sunshine duration, R_a is the solar radiation at the top of the atmosphere, and N is the maximum possible duration of sunshine or daylight hours. On cloudless days, the actual sunshine duration equals the number of daylight hours (n = N) and their ratio is one. a_s is a regression constant expressing the fraction of extraterrestrial radiation reaching the earth on overcast days (when n= 0), and a_s + b_s denotes the fraction of extraterrestrial radiation reaching the earth on clear days (n = N).

Sensitivity analysis

Sensitivity analysis measures the reaction of models to changes in parameters. The purpose of this analysis is to select the parameters that largely affect the model output. The accuracy of these sensitive parameters determines the reasonableness of the simulation results.

SIMRIW exports the heading date, maturity date and yield from the observational dataset at each station. These variables are used for calibration and validation. Thus, the sensitivity analysis was conducted on the phenology (heading and maturity dates) and yield.

The sensitivity of each parameter in SIMRIW (Table 1) is determined by varying that parameter while retaining the other parameters at their default values.

Parameter estimation and model validation

Based on the results of the sensitivity analysis, a SIMRIW parameter estimation was undertaken by the SCE–UA method (Duan et al., 1994).

The phenological development and yield calculation parameters were calibrated using phenological data (heading date and maturity date) collected from 1990 to 2000, and yield data collected from 2000 to 2009.

The accuracy of the models was evaluated by the root mean square error (RMSE) between the observed and simulated values:

(19)RMSE=\sqrt{{\displaystyle \frac{{\sum _{i=\mathrm{1}}^n\left(Simulate{d}_i-Observe{d}_i\right)}^{\mathrm{2}}}{n}}},

where n is the number of comparisons.

Results

Sensitivity analysis

The sensitivity analysis identified six (B, L_c, G_v, A_T, T_h, DVI*), three (K_r, T_cr, G_r) and four (h_m, T_cf, F_as, K_r) parameters relevant to the heading date, maturity date and yield calculation, respectively.

The phenological development of rice during the DVI<1 period (growth stage before heading) involves six SIMRIW parameters. The sensitivity of each parameter is shown in Fig. 1. The heading date is largely affected by the parameter G_v, consistent with the inverse relationship between DVR and G_v shown in Eq. (3). The simulation result also significantly depended on L_c, A_T and T_h (Fig. 1(a)), although these parameters exerted less effect on heading date during the 1<DVI<2 period (Fig. 1(b)).

The yield formulation (contain the dry matter production) in SIMRIW is chiefly affected by h_m, and is insensitive to T_cf and F_as. The parameters are optimized in the parameter estimation (Fig. 1(c)).

Parameter optimization

Along with the default values, Table 1 lists the ranges of optimal values of the estimated parameters in SIMRIW.

Among the four important SIMRIW parameters governing the phenological development of rice, G_v (the minimum number of days required for heading) varies much more widely among ecological zones than the other parameters during the 0<DVI<1 period. A_T, denoting the sensitivity of the development rate (DVR) to air temperature, increases in the order late-rice>single-rice>early-rice, although late-rice and early-rice developments are reasonably consistent among ecological zones. The A_T of single-rice is lower in zone I than in the other zones. T_h, the air temperature at which DVR is half the maximum rate at the optimum temperature (°C), is lower for late-rice than for early-rice and single-rice (Table 2).

The 1<DVI<2 period is most greatly affected by G_r; the remaining parameters (K_r and T_cr) little vary among the ecological zones (Table 2).

Concerning the three parameters in the yield formulation, h_m, denoting the maximum harvest index of a given prefecture under optimum climatic conditions and cultivation practices, is relatively invariant among ecological zones. F_as, the asymptotic value of the leaf area index under non-limiting temperature conditions, is high for single-rice and lower for late-rice, whereas T_cf (minimum temperature for LAI growth) is lower for early-rice than single-rice and late-rice (Table 2).

Model validation

The primary purpose of rice models is to simulate the rice phenology. The accuracy of the phenology simulation is critical in yield simulations. Since the rice phenology in each ecological zone depended on the climatic conditions, we optimized the parameter sets in each simulated growing season in each ecological zone. The simulated heading and maturity dates were validated in each zone.

The mean simulated heading date for single-rice was DOY 218.71, 231.21, 224.66, 210.60, and 220.27 in zones I, II, III, IV, and V, respectively. The RMSE was minimized in zones I and III. The mean heading date for early-rice was 176.91 in zone III (RMSE= 3.82 days), 171.81 in zone IV, and 160.25 in zone VI. The mean heading date of late-rice was 256.00, 256.24, and 274.50 in zones III, IV, and VI, respectively. The SIMRIW model more accurately simulated the heading date for single-rice than for early-rice and late-rice. It also output more accurate single-rice heading dates in zone III than in other zones. Conversely, the heading date of early-rice was most accurately simulated in zone IV (Table 3).

The simulated maturity date of single-rice in zone I was 263.12, whereas that of early-rice, single-rice and late-rice in zone II was 199.49, 257.91 and 292.87, respectively. In zone III, the mean simulated maturity date of early-rice, single-rice and late-rice was 198.10, 235.70 and 296.53, respectively. Similarly to the heading date, the SIMRIW-simulated maturity date was more accurate for single-rice than for early-rice and late-rice (Table 3).

Table 3 also lists the SIMRIW-simulated rice yields. Once calibrated, the SIMRIW model properly simulated the rice yield in major rice production areas in China. Comparing the observed and simulated yields in Table 3, we observe that single-rice yield was most accurately simulated in zone V, and that single-rice yields were more accurately simulated than early-rice and late-rice yields.

Discussion

Only a few of the many parameters in SIMRIW exert a large impact on the model output. Once these parameters were optimized for a given ecological zone and rice type, SIMRIW accurately simulated the rice phenology and yield in single and double cropping systems in different ecological zones. Some of the parameters largely vary among ecological zones and rice types, probably reflecting the different climate conditions and rice varieties among ecological zones.

The SIMRIW-estimated parameters clearly depend on crop region. The minimum number of days required for heading, G_v exerts a major impact on the heading date calculation, mainly through the inverse relationship between DVR and G_y shown in Eq. (3). During the 1<DVI<2 period, heading date is most sensitive to the parameter G_r. Collectively, G_v and G_r determine the base length of the rice growth duration. In most ecological zones, the duration of single-rice exceeds that of early-rice and late-rice. The optimized A_T was relatively high for late-rice, and low for early-rice. This indicates that late-rice is more sensitive to temperature than single-rice and early-rice. The A_T of early-rice and late-rice was relatively invariant among ecological zones, but was lower in zone I than in other zones for single-rice. Therefore, single-rice is less sensitive to temperature in zone I than in other zones.

Conclusions

This study has identified the sensitive parameters in the SIMRIW model. These sensitive parameters varied widely among the investigated ecological zones. For single-rice, the heading date was more accurately simulated than the maturity date. The RMSE of the heading date simulation ranged from 1.25 to 4.54 days, whereas that of the maturity date simulation ranged from 2.28 to 5.11 days. The simulated heading and maturity dates were more accurate for single-rice than for early-rice and late-rice. The yield RMSEs ranged from 0.78 t/ha to 0.89 t/ha, with the best results obtained in zone V. We conclude that once the parameters were optimized, SIMRIW properly simulated the rice phenology and yield in single and double cropping systems in various ecological zones.

References

[1] Anderson T R (2010). Progress in marine ecosystem modelling and the “unreasonable effectiveness of mathematics”. J Mar Syst, 81(1–2): 4–11

[2] Angström A (1924). Solar and terrestrial radiation. Report to the international commision for solar research on actinometric investigations of solar and atmospheric radiation. Q J R Meteorol Soc, 50(210): 121–126

[3] Bouman B A M, Van Laar H H (2006). Description and evaluation of the rice growth model ORYZA2000 under nitrogen-limited conditions. Agric Syst, 87(3): 249–273

[4] Brouwer R, De Wit C T (1969). A simulation model of plant growth with special attention to root growth and its consequences. In: Proceedings of the 15th Easter School in Agri Science. University of Nottingham, 224–242

[5] Duan Q, Sorooshian S, Gupta V K (1994). Optimal use of the SCE-UA global optimization method for calibrating watershed models. J Hydrol (Amst), 158(3–4): 265–284

[6] Duncan W G (1971). Simulation of growth and yield in cotton: A computer analysis of the nutritional theory. Proc Beltwide Cotton Prod Res Conf, (78): 45–61

[7] Gao L, Jin Z, Huang Y, Zhang L (1992). Rice clock model— A computer model to simulate rice development. Agric Meteorol, 60(1–2): 1–16

[8] Hayashi Y, Jung Y S (2000). Paddy rice production under possible temperature fluctuation in East Asia. Global Environmental Research, 3: 129–137

[9] Horie T, Centeno G S, Nakagawa H, Matsui T (1997). Effect of elevated carbon dioxide and climate change on rice production in East and Southeast Asia. In: Proceedings of the International Scientific Symposium on Asian Paddy Fields. Saskatchewan: University of Saskatchewan, 913–967

[10] Horie T, Nakagawa H, Ceneno H G S, Kropff M (1995). Rice production in Japan under current and future climates. In: Matthews R B, Kropff M J, Bachelet D, eds. Modeling the Impact of Climate Change on Rice in Asia. Wallingford: CAB International, 143–164

[11] Kropff M J, Centeno G S, Bachelet D, Lee M H, Mohan Dass S, Horie T, De Feng S, Singh S, Penning De Vries F W T (1993). Predicting the impact of CO₂ and temperature on rice production. In: IRRI Seminar Series on Climate Change and Rice. Los Banos: International Rice Research Institute

[12] Lin E, Xiong W, Ju H, Xu Y, Li Y, Bai L, Xie L (2005). Climate change impacts on crop yield and quality with CO₂ fertilization in China. Philos Trans R Soc Lond B Biol Sci, 360(1463): 2149–2154

[13] Matthews R B, Kropff M J, Horie T, Bachelet D (1997). Simulating the impact of climate change on rice production in Asia and evaluating options for adaptation. Agric Syst, 54(3): 399–425

[14] Mei F, Wu X, Yao C, Li L, Wang L, Chen Q (1988). Rice cropping regionalization in China. Chinese Journal Rice Science, 2(3): 97–110 (in Chinese)

[15] Peng S, Huang J, Sheehy J E, Laza R C, Visperas R M, Zhong X, Centeno G S, Khush G S, Cassman K G (2004). Rice yields decline with higher night temperature from global warming. Proc Natl Acad Sci USA, 101(27): 9971–9975

[16] Peng S, Ingram K T, Neue H U, Ziska L H (1995). Climate Change and Rice. Berlin, Heidelberg: International Rice Research Institute, Springer-verlag

[17] Prescott J A (1940). Evaporation from a water surface in relation to solar radiation. Transactions of the Royal Society of South Australia, 64: 114–125

[18] Spitters C J T (1986). Separating the diffuse and direct component of global radiation and its implications for modeling canopy photosynthesis Part II. Calculation of canopy photosynthesis. Agric Meteorol, 38(1–3): 231–242

[19] Tao F, Yokozawa M, Hayashi Y, Lin E (2003). Future climate change, the agricultural water cycle, and agricultural production in China. Agric Ecosyst Environ, 97(1–3): 361

[20] Tao F, Yokozawa M, Xu Y, Hayashi Y, Zhang Z (2006). Climate changes and trends in phenology and yields of field crops in China, 1981‒2000. Agric Meteorol, 138(1–4): 82–92

[21] Wopereis C S, Defoer T, Idinoba P, Diack S, Dugué M J (2009). Effects of temperature on rice. In: Wopereis M C S, Defoer T, Idinoba P, Diack S, Dugue M J, eds. Curriculum for Participatory Learning and Action Research (PLAR) for Integrated Rice Management (IRM) in Inland Valleys of Sub-Saharan Africa: Technical Manual. Cotonou: Africa Rice Center, 39–42

[22] Xiong W, Tao F, Xu Y, Lin E (2001). Simulation of rice yields under future climate scenarios in China. Agric Meteorol, 22: 9–14

[23] Yin X, Kropff M J (1996). Use of the Beta function to quantify effects of photoperiod on flowering and leaf number in rice. Agric Meteorol, 81(3–4): 217–228

[24] Yoshida S (1981). Fundamentals of Rice Crop Science. Los Banos: International Rice Research Institute, 269

RIGHTS & PERMISSIONS

Higher Education Press and Springer-Verlag Berlin Heidelberg