Global sensitivity analysis of outputs over rice-growth process in ORYZA model

Environmental Modelling & Software 83 (2016) 36e46

Contents lists available at ScienceDirect

Environmental Modelling & Software journal homepage: www.elsevier.com/locate/envsoft

Global sensitivity analysis of outputs over rice-growth process in ORYZA model Junwei Tan, Yuanlai Cui, Yufeng Luo* State Key Laboratory of Water Resources and Hydropower Engineering Science, Wuhan University, Wuhan 430072, China

a r t i c l e i n f o

a b s t r a c t

Article history: Received 14 August 2015 Received in revised form 15 March 2016 Accepted 1 May 2016

Dynamic crop models usually have a complex structure and a large number of parameters. Those parameter values usually cannot be directly measured, and they vary with crop cultivars, environmental conditions and managements. Thus, parameter estimation and model calibration are always difficult issues for crop models. Therefore, the quantification of parameter sensitivity and the identification of influential parameters are very important and useful. In this work, late-season rice was simulated with meteorological data in Nanchang, China. Furthermore, we conducted a sensitivity analysis of 20 selected parameters in ORYZA_V3 using the Extended FAST method. We presented the sensitivity results for four model outputs (LAI, WAGT, WST and WSO) at four development stages and the results for yield. Meanwhile, we compared the differences among the sensitivity results for the model outputs simulated in cold, normal and hot years. The uncertainty of output variables derived from parameter variation and weather conditions were also quantified. We found that the development rates, RGRLMN and FLV0.5 had strong effects on all model outputs in all conditions, and parameters WGRMX and SPGF had relative high effects on yield in cold year. Only LAI was sensitive to ASLA. Those influential parameters had unequal effects on different outputs, and they had different effects at four development stages. With the interaction effects of parameter variation and different weather conditions, the uncertainty of model outputs varied significantly. However, the weather conditions had negligible effects on the identification of influential parameters, although they had slight effects on the ranks of the parameters' sensitivity for outputs in the panicle-formation phase and the grain-filling phase, including yield at maturity. The results suggested that the influential parameters should be recalibrated in priority and fine-tuned with higher accuracy during model calibration. © 2016 Elsevier Ltd. All rights reserved.

Keywords: Crop model Extended FAST Uncertainty Weather conditions Model calibration

1. Introduction Dynamic crop models are valuable tools for agro-environment research. They are very useful in yield prediction (Ma et al., 2013; Wu et al., 2013), climate change impact evaluation (Asseng et al., 2013; Wang et al., 2014), decision making on field management (Feng et al., 2007; Yadav et al., 2011; Fang et al., 2012), and crop genotype improvement (Boote et al., 2001; Li et al., 2013). These models are always highly complex and include many parameters. Inaccurate values of these parameters can lead to unreliable model predictions (Makowski et al., 2006; Lamboni et al., 2009; Confalonieri et al., 2010a). The estimation of certain parameter is

* Corresponding author. E-mail addresses: [email protected] (J. Tan), [email protected] (Y. Cui), yfl[email protected] (Y. Luo). http://dx.doi.org/10.1016/j.envsoft.2016.05.001 1364-8152/© 2016 Elsevier Ltd. All rights reserved.

based on specific measurement or extracted from literature sources (Wallach et al., 2001; Confalonieri et al., 2006), and many parameter values vary with crop cultivars, environmental conditions and managements (Wang et al., 2013; Vazquez-Cruz et al., 2014). Therefore, parameter estimation and model calibration are usually tough jobs due to the uncertainty of the parameter values (Gallagher and Doherty, 2007; Wallach et al., 2014). Generally, the model parameters cannot be calibrated against observed data simultaneously at one time (Tremblay and Wallach, 2004). In fact, only a small number of parameters is responsible for most of the variability in the model outputs; these parameters deserve an accurate determination, whereas others with a small influence can be fixed at a nominal value (Makowski et al., 2006; Varella, 2012; Wang et al., 2013). Sensitivity analysis (SA) is defined as the study of how the uncertainty in the output of a model can be apportioned to different sources of uncertainty in the model parameters (Saltelli et al., 2004). SA is used to identify the

J. Tan et al. / Environmental Modelling & Software 83 (2016) 36e46

sensitive parameters, to simplify the number of model parameters and to evaluate model robustness by setting insensitive parameters as default values (Cariboni et al., 2007; Vanuytrecht et al., 2014). Usually, a simple One-at-a-Time (OAT) sensitivity analysis method (local sensitivity analysis) is likely to be used, whereas it may be perfunctory for complicated nonlinear models (Saltelli and Annoni, 2010). Instead, variance-based global sensitivity analysis (GSA) methods are suggested, e.g., the Morris method and the Fourier Amplitude Sensitivity Test (FAST; Saltelli et al., 2008). Global sensitivity analysis can be used to analyze the model output uncertainties over the entire range of parameters. GSA quantifies not only the individual effect of a parameter (the main effect) but also potential interactions among selected parameters (the interaction effect). The SA results are influenced by many factors; the ranking of parameter sensitivity may depend on parameter variation range (Richter et al., 2010; Wang et al., 2013) or on the climate, site, weather conditions (Confalonieri et al., 2010a, 2010b; Shin et al., 2013) or agriculture management (Zhao et al., 2014). ORYZA (Bouman et al., 2001) is a crop-growth model that can be used to simulate the rice-growing process and predict yield under multiple conditions with interaction of genotype
environment
management. The model has been applied and validated in different regions on a global scale (Bouman and van Laar, 2006). However, in the past literature, there is few research with a systematical description of parameter sensitivity analysis for the ORYZA model. Soundharajan and Sudheer (2013) used Sobol's method to assess parameter sensitivity on rice yield, and based on the results of the sensitivity analysis, the author proposed an auto-calibration procedure for ORYZA. Nevertheless, the parameter sensitivity and uncertainty on output state variables over rice-growth process in the model remain unknown. Therefore, the model user lacks knowledge about key parameters in calibration and how much influence they have on model outputs. In this study, we used the Extended FAST method implemented in the professional sensitivity analysis tool SIMLAB software (http:// simlab.jrc.ec.europa.eu/) to pursue following objectives: (i) to distinguish influential parameters and uninfluential parameters, (ii) to explore the effects of weather conditions on the sensitivity results and the uncertainty of model outputs, (iii) to explore the temporal characteristics of parameter sensitivity over the ricegrowth period, and (iv) to quantify the uncertainty of model outputs at different development stages derived from parameter variation. The answers to these objectives are expected to be helpful to quickly and accurately calibrate a model by identifying influential parameters, and to evaluate model balance and model improvement by offering a good understanding of the effects of parameters on model outputs at different development stages. 2. Materials and methods In this work, ORYZA_V3, the latest version since ORYZA Version 2.13 released in 2009 (https://sites.google.com/a/irri.org/ oryza2000/home), was used for the parameter sensitivity analysis. In the model, the development stage (DVS) of rice crop is used to define its physiological age, and the life cycle of rice crop is divided into four phenological phases: the basic vegetative phase, the photoperiod-sensitive phase, the panicle-formation phase and the grain-filling phase (Bouman et al., 2001). At the end of each phase (DVS ¼ 0.4, DVS ¼ 0.65, DVS ¼ 1.0, DVS ¼ 2.0), the leaf area index (LAI, ha/ha), the dry weight of stems (WST, kg/ha), the total aboveground dry matter (WAGT, kg/ha) and the dry weight of storage organs (WSO, kg/ha) were chosen as output state variables. The Extended FAST algorithm was adopted to calculate the main effect index and the total sensitivity index of each considered parameter.

37

2.1. ORYZA model The ORYZA model is an eco-physiological simulation model of the ‘School of De Wit’ (Bouman et al., 1996). The model can be used to simulate the growth, development, and water balance of rice under potential production, water-limited and N-limited situations (Bouman et al., 2001). For all situations, the crop is assumed to be well protected against diseases, pests, and weeds. A detailed description of ORYZA is given by Bouman et al. (2001). The ORYZA model simulates the daily rate of dry matter production in plant organs and the daily rate of phenological development. By integrating these rates over time, dry matter production and the development stage are simulated throughout the growing season. The phenological development rate is tracked as a function of the daily mean temperature and the photoperiod. The daily assimilation rate is calculated by integrating instantaneous rates of leaf CO2 assimilation over the day and the depth of the canopy. The integration is based on an assumed sinusoidal time course of radiation during the day and on an exponential extinction of radiation within the canopy. Photosynthesis of single leaves depends on leaf N content (on area basis, g/m2), radiation, and temperature. The daily dry matter increment is obtained by subtracting maintenance and respiration requirements from total assimilation rate. Then, the daily dry matter increment is partitioned among roots, leaves, stems and panicles using experimentally derived partitioning factors as a function of the development stage. The spikelet density at flowering is derived from the total biomass accumulated from panicle initiation to first flowering. The spikelets turn into grains filled with assimilates until their maximum grain weight is reached. Additionally, spikelet sterility is taken into consideration in cases of temperatures that are too high or too low. The leaf area growth includes two phase: exponential growth phase and linear growth phase. In the early phase of growth, the leaf area grows exponentially as a function of the temperature sum and the relative leaf-growth rate. When the canopy overlaps, e.g., LAI is larger than 1, the leaf area growth is limited by the amount of carbohydrates available for leaf growth. In the linear phase, the increase in leaf area is obtained from the increase in leaf weight using the specific leaf area. When the rice is transplanted, LAI and biomass become reduced due to the dilution of plant density. Crop growth resumes only after an end of “transplanting shock”, which is derived from the accumulated temperature in the seedbed (Bouman and van Laar, 2006; Feng et al., 2007; Amiri and Rezaei, 2010). 2.2. Model parameters In the ORYZA model, all parameter values are listed in external data files and can be changed by the model user. About 10% of crop parameters are expected to be variety specific and need empirical derivation (Bouman and van Laar, 2006). In this work, 20 parameters were considered. The description and change variation of these parameters are presented in Table 1. For most parameters, the base values are obtained from default values in model file. In particular, the partitioning factors, leafdeath rates and the development rates are set as an average of those values calculated with data (unpublished) from experiments conducted from 2011 to 2013 at the Jiangxi Irrigation Experiment Station in Nanchang, Jiangxi Province, China. The parameter change variation was set as ±30% perturbation of its base value. 2.3. Study site and data The Jiangxi Irrigation Experiment Station is located in Poyang

38

J. Tan et al. / Environmental Modelling & Software 83 (2016) 36e46

Table 1 The parameter description and variance in the ORYZA model. Parameter

DVRJ DVRI DVRP DVRR RGRLMX RGRLMN ASLA BSLA CSLA DSLA SLAMAX FLV0.5 FLV0.75 FST1.0 DRLV1.0 DRLV1.6 DRLV2.1 FSTR SPGF WGRMX

Description

Unit

i¼1

Fraction of carbohydrates allocated to stems stored as reserve Spikelet growth factor Maximum individual grain weight

Vi þ

X

Vij þ / þ V12…n

(1)

1i < jn

whereVðYÞ denotes the total variance of model output Y induced by the ORYZA model, Vi ¼ V½EðY=xi Þ denotes the variance allocated to each parameterxi, EðY=xi Þ denotes the expectation of output variableY conditional on a fixed value of each parameterxi, and Vij ; …; V12…n denotes the variance allocated to interactions among parameters. For each parameter, two sensitivity indices are calculated. The first-order (main effect) index (Si ) measures the contribution of single parameterxi to the uncertainty of output; Si is defined by Eq. (2):

Si ¼

Vi VðYÞ



Fraction shoot dry matter partitioned to the stems at DVS ¼ 1.0 Leaf death coefficient as a function of development stage

The Extended FAST is a variance-based global sensitivity analysis algorithm, which was developed by integrating the high computing efficiency of the Fourier Amplitude Sensitivity Test (FAST) algorithm and Sobol's capability to calculate high-order interactions among parameters (Saltelli et al., 1999). We usedY to represent the output of the ORYZA model. The total variance VðYÞ of model output can be decomposed by Eq. (1) as follows (Makowski et al., 2006): n X

C/day C/day C/day  C/day  C/day  C/day e e e e ha/kg e e e e e e e no./kg kg/grain 

Maximum value of SLA Fraction of shoot dry matter partitioned to the leaves at DVS ¼ 0.5 and DVS ¼ 0.75

2.4. Extended FAST method

VðYÞ ¼



Development rate in juvenile phase Development rate in photoperiod-sensitive phase Development rate in panicle development Development rate in reproductive phase Maximum relative growth rate of leaf area Minimum relative growth rate of leaf area Parameters of a function to calculate specific leaf area (SLA, ha/kg)

Lake Basin, Jiangxi Province, China. This region is characterized by a warm and humid subtropical monsoon climate with a mean annual temperature of 18.1  C and a mean annual precipitation of 1636 mm. Late-season rice was the major focus of this study, and the phenology investigation was carried out among farmers around the station. The ORYZA_V3 model was used for late-season rice-growth simulation under the potential production condition. The meteorological data were obtained from an automated weather station (115 550 E, 28 360 N) located in the study region, including daily values of minimum temperature, maximum temperature, sunshine hours, vapor pressure and mean wind speed 2 m above the ground surface.

Base value

Min

Max

0.0007 0.000525 0.000595 0.0014 0.00595 0.0028 0.00168 0.00175 5.85 0.098 0.00315 0.42 0.21 0.28 0.014 0.021 0.035 0.175 45,430 0.0000175

0.0013 0.000975 0.001105 0.0026 0.01105 0.0052 0.00312 0.00325 3.15 0.182 0.00585 0.78 0.39 0.52 0.026 0.039 0.065 0.325 84,370 0.0000325

0.001 0.00075 0.00085 0.002 0.0085 0.004 0.0024 0.0025 4.5 0.14 0.0045 0.6 0.3 0.4 0.02 0.03 0.05 0.25 64,900 0.000025

The total sensitivity index (STi ) measures the sum of all contributions involving parameterxi to the uncertainty of output; STi is defined by Eq. (3):

STi ¼

X

Si þ

X

Sij þ

jsi

X

Sijm þ / þ S1…i…n

isjsm

¼

VðYÞ  Vi VðYÞ

(3)

whereSij denotes the second-order sensitivity index for the couple of parameterxi and any other parameterxj, Sijm denotes the thirdorder sensitivity index for the combination of parameterxi and any other two parameters, and so on, and Vi denotes the sum of all the contributions that do not include parameterxi to the uncertainty of output. Si and STi are both in the range of (0, 1). Low values indicate negligible effects, and high values close to 1 indicate more important effects (Vanuytrecht et al., 2014).

2.5. Sensitivity analysis In this study, we set the parameter sample size to N ¼ 11,000 for the sensitivity analysis based on the knowledge of sensitivity indices' convergence. To assess the influence of weather conditions on the sensitivity results, the late-season rice was simulated in three years (1980, 1985, 2008) of contrasting weather types, which were chosen based on the annual accumulated temperature in the late-season rice-growing duration with a 95% (cold year), 50% (normal year), and 5% (hot year) probability of exceedance, according to the Pearson III distribution over 60 years (1954e2013) in Nanchang, China. Several statistical indicators were also used to quantify the uncertainty of output variables as follows from Eq. (4) to Eq. (6):

MðYÞ ¼ y ¼ (2)

Variance

N X i¼1

, yi

N

(4)

J. Tan et al. / Environmental Modelling & Software 83 (2016) 36e46

vffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi , ffi u n uX 2 SDðYÞ ¼ t ðyi  yÞ N

(5)

i¼1

CVðYÞ ¼ SDðYÞ=MðYÞ

(6)

whereMðYÞ andy represent the average values of the model output Y, yi represents the value of the model output obtained from model execution with a combination of parameter samples,N represents the sample size of model parameters, SDðYÞ represents the standard deviation of the model outputY, andCVðYÞ represents the coefficient of the variation of the model outputY.

3. Results 3.1. Sensitivity indices of parameters for model outputs at different development stages 3.1.1. Sensitivity indices of parameters for model outputs at basic vegetative phase Fig. 1 shows the sensitivity indices of the 20 parameters for model outputs at the end of the basic vegetative phase (DVS ¼ 0.4) simulated in cold, normal and hot years, respectively. Different weather conditions have almost no influence on the sensitivity results of the selected parameters for LAI, WAGT and WST at this stage; e.g., the rank of those influential parameter sensitivities is the same, though some sensitivity indices change slightly. The development rate in the juvenile phase (DVRJ) and the maximum relative growth rate of the leaf area (RGRLMX) are the most influential parameters for all outputs. The fraction of shoot dry matter partitioned to the leaves at DVS ¼ 0.5 (FLV0.5) shows essential effects on LAI and WST and negligible effects on WAGT. It is worthy to note that in the model, the partition factors for organs at different development stages are linearly interpolated between the adjacent values of two development stages. Before panicle-formation, the fraction of dry matter partitioned to stems can be calculated by subtracting the fraction of leaves from 1. Moreover, according to Fig. 1, the interaction effects among parameters on all model outputs are notable. For outputs LAI, WAGT and WST, the interaction indices (the white part of the bars in the figure) of DVRJ contribute to 41%, 32%, and 34%, respectively, of the total sensitivity index value at the average level under three weather conditions, while that of RGRLMX contributes to approximately 32%, 28%, and 30%, respectively. The interaction effects of parameter FLV0.5 contribute to approximately 31% for the variance of WST and approximately 62% for LAI.

3.1.2. Sensitivity indices of parameters for model outputs at photoperiod-sensitive phase Fig. 2 shows the sensitivity indices of the 20 parameters for model outputs at the end of the photoperiod-sensitive phase (DVS ¼ 0.65) simulated in cold, normal and hot years, respectively. At this stage, the effects of single parameters on all model outputs are significant, whereas the interaction effects are non-significant. Additionally, the weather type make no differences in the sensitivity results. All model outputs are sensitive to the development rate in the photoperiod sensitive phase (DVRI), DVRJ and RGRLMX. Parameter FLV0.5 is the second most sensitive parameter for LAI, and it is the fourth most sensitive for WAGT, whereas it is insensitive for WST. Parameter ASLA (one of parameters to calculate specific leaf area) has slight effects on WAGT and WST, whereas it plays an important part in the variation of LAI.

39

3.1.3. Sensitivity indices of parameters for model outputs at panicleformation phase Fig. 3 presents the sensitivity indices of the 20 parameters for model outputs at the end of the panicle-formation phase (DVS ¼ 1.0) simulated in cold, normal and hot years, respectively. At this stage, the weather type slightly changes the ranks of some influential parameters' sensitivity, whereas it has no effects on the identification of influential parameters. However, the difference in parameter sensitivity ranks caused by weather type exists only for those parameters with equivalent effects on output variation, and it shows a poor change rule. RGRLMX, DVRJ and DVRI still have a significant influence on outputs LAI, WAGT and WST, whereas they have negligible effects on WSO. Parameter ALSA is the most sensitive for LAI, whereas it has minor effects on the others. Parameter FLV0.5 has an important effect on LAI and WAGT but little effect on WST and WSO. The development rate in the panicle-formation phase (DVRP) has the strongest effects on WSO, and the maineffect indices of that contribute to approximately 60% of the variation of WSO at an average level under different weather conditions. The fraction of dry matter partitioned to the storage organ at DVS ¼ 0.1 (FST1.0) ranks second for sensitivity, and the main effect indices of that contribute to approximately 20% of the variation of WSO. Furthermore, the interaction effects among those parameters are negligible for the outputs in this stage, particularly for WAGT and WST. 3.1.4. Sensitivity indices of parameters for model outputs at grainfilling phase Fig. 4 shows the sensitivity indices of the 20 parameters for model outputs at the end of the grain-filling phase (DVS ¼ 2.0) simulated in cold, normal and hot years, respectively. The weather type has more effects on the sensitivity results for WST and WSO at this stage than at previous stages; e.g., the WST at DVS ¼ 2.0 is slightly sensitive to RGRLMX in hot years but not in cold or normal years. Similarly, the WST at DVS ¼ 2.0 is slightly sensitive to the maximum individual grain weight (WGRMX) in cold years but not in other years. Additionally, the development rate in the grainfilling phase (DVRR) has more minor effects on WSO in cold year than in normal and hot years. Nevertheless, the sensitivity results for other parameters have no differences under contrasting weather conditions. Parameter DVRP and RGRLMX have strong effects on WAGT and WSO but minor effects on LAI, WST. The spikelet growth factor (SPGF) and parameter WGRMX have minor effects on WST and WSO, whereas they have negligible effects on LAI and WAGT. For WST, the interaction effects between DVRR and other parameters is particularly high, contributing to approximately 66% of the total sensitivity index value at an average level under different weather conditions. 3.2. Coefficient of variance for state variables over rice-growth period Table 2 shows the statistics results of output variables at different development stages simulated in three weather types. According to Table 2, at DVS ¼ 0.4, the mean values of LAI, WAGT and WST simulated in cold year are the highest, whereas those in hot year are the lowest. However, the coefficients of variation (CV) of LAI and WAGT are lower in normal year than in cold and hot years. In particular, the CV of WST is lower in hot year than in normal year, whereas the CV of WST is higher in cold year than in normal year. In addition, the common point under the three weather conditions is that the CV of LAI is higher than that of both WAGT and WST. At DVS ¼ 0.65, the mean values of LAI, WAGT and WST in normal year are the highest, whereas those in hot year are the lowest.

40

J. Tan et al. / Environmental Modelling & Software 83 (2016) 36e46

Fig. 1. Sensitivity indices of 20 parameters in ORYZA_V3 model calculated by Extended FAST method for rice-growth state variables at DVS ¼ 0.40: (a) outputs for cold year (1980, (b) outputs for normal year (1985), and (c) outputs for hot year (2008). First order indices are in black and interactions are in white.

However, the CV of WAGT and WST is higher in normal year than in hot and cold years. Compared to DVS ¼ 0.4, the mean value of LAI increases from 2.09 to 4.30 at DVS ¼ 0.65 in normal yeardincreasing by 1.06 timesdwhereas it increases by 0.80 and 0.97 times in cold and hot years, respectively. At an average level, the mean values of WAGT and WST from DVS ¼ 0.4 to DVS ¼ 0.65 increase by 1.87 and 2.32 times, respectively. In addition, the CV of all output variables at DVS ¼ 0.65 is higher than those at DVS ¼ 0.40. At DVS ¼ 1.0, the mean values of LAI, WAGT, WST and WSO in normal year are the highest among those simulated with different weather types. The mean values of LAI, WAGT and WST in cold year are the lowest, whereas the mean values of WSO in hot year is the lowest. Even so, the CV of LAI, WAGT, WST and WSO in cold year is higher than that in normal and hot years, whereas those in hot year are the lowest. On the contrary, the relative increments of LAI, WAGT, WST and WSO from DVS ¼ 0.4 to DVS ¼ 0.65 in hot year are the highest among those in three weather types, whereas those in cold year are the lowest. Compared to DVS ¼ 0.65, the mean values

of LAI, WAGT and WST increase by 0.31, 2.14, and 2.73 times at an average level, respectively; however, the CV of LAI, WAGT and WST decreases by 53%, 48%, and 49% at the average level, respectively. At DVS ¼ 2.0, the mean values of LAI, WAGT, WST and WSO in normal year are the highest among those simulated with different weather types, whereas the mean value of LAI in cold year is the lowest, and the mean values of WAGT, WST and WSO in hot year are the lowest. The CV of LAI and WSO is lower in normal year than in hot and cold year. The CV of WAGT, and WST is lower in normal year than in cold year, and it is higher than that in hot years. Compared to DVS ¼ 1.0, at an average level, the mean value of LAI greatly decreases by approximately 61%, and the mean value of WST slightly decreases by approximately 10%. On the contrary, the mean value of WSO rapidly increases by 5.16 times, and the mean value of WAGT modestly increases by 0.64 times. Moreover, the CV of LAI and WST at DVS ¼ 2.0 increases significantly compared to that at DVS ¼ 1.0, whereas the CV of WAGT and WSO decreases slightly, except for the CV of WSO in hot year. The data in Table 2 indicate that for late-season rice, when the

J. Tan et al. / Environmental Modelling & Software 83 (2016) 36e46

41

Fig. 2. Sensitivity indices of 20 parameters in ORYZA_V3 model calculated by Extended FAST method for rice-growth state variables at DVS ¼ 0.65: (a) outputs for cold year (1980), (b) outputs for normal year (1985), and (c) outputs for hot year (2008). First order indices are in black and interactions are in white.

selected parameters change in space, the CV of LAI gradually increases from the basic vegetation phase until the maximum value in the photo-sensitive phase. Then, it greatly decreases to the lowest point in the panicle-formation phase. Finally, in the grainfilling phase, it increases to a new high again. The variation of the CV of WST is the same to that of LAI. Meanwhile, the CV of WAGT first increases to a maximum value in the photo-sensitivity phase. Then, it greatly decreases in the panicle-formation phase, and finally, it slows down to the lowest value in the grain-filling phase. In addition, the CV of WSO remains relatively steady from flowering (DVS ¼ 1.0) to maturity (DVS ¼ 2.0). The mean value of LAI has an impressive expansion in the basic vegetative phase and the photo-sensitive phase and a modest increase in the panicle-formation phase. Then, it decreases after flowering. The mean value of WST maintains stable growth in the photo-sensitive phase before increasing greatly in the panicleformation phase. Then, it slows down and even has a modest decrease in the grain-filling phase. The mean value of WAGT greatly increases in the panicle-formation phase. Then, it slows down and obtains a maximum increment in the grain-filling phase.

Meanwhile, the mean value of WSO rapidly increases in the grainfilling phase.

3.3. Sensitivity indices of parameters for yield at maturity Normally, crop yield remains the greatest concern in crop models. Accordingly, we calculated the main effect indices and the total sensitivity indices of the 20 parameters on WRR14 (the yield at maturity with 14% moisture) simulated with the three weather types, as shown in Fig. 5. The figure shows that the result is similar to WSO at DVS ¼ 2.0 (see Fig. 4); however, the difference caused by weather conditions is more significant. It is worth noting that the WRR14 in the ORYZA_V3 model is calculated by the growth rate of the number of grains; it is not simply calculated from WSO at maturity. For WRR14, the most influential parameter in cold year is DVRP, whereas it is DVRR in normal and hot year. Parameter SPGF has stronger effects on WRR14 in cold and normal years than in hot year, and parameter WGRMX has essential effects on WRR14 only in cold year, even though both parameters have slight main effects

42

J. Tan et al. / Environmental Modelling & Software 83 (2016) 36e46

Fig. 3. Sensitivity indices of 20 parameters in ORYZA_V3 model calculated by Extended FAST method for rice-growth state variables at DVS ¼ 1.0: (a) outputs for cold year (1980), (b) outputs for normal year (1985), and (c) outputs for hot year (2008). First order indices are in black and interactions are in white.

on WRR14 in all conditions. After all, the development rates in different stages (DVRJ, DVRI, DVRP, DVRR), RGRLMX and FLV0.5 are the most influential parameters. Parameters ASLA, SPGF and WGRMX have negligible main effects and low total sensitivity indices. The other parameters are non-influential.

in the grain-filling phase. Parameter RGRLMX stays stably influential over the rice-growth period. Additionally, parameter ASLA has the strongest effects on LAI in the panicle-formation phase. Meanwhile, parameter FLV0.5 has the strongest effects on LAI in the photo-sensitive phase, and it has the strongest effects on WST in the basic vegetative phase. Parameters SPGF and WGRMX have slight effects only in the grain-filling phase.

3.4. Sensitivity change of key parameters over rice growth period According to Figs. 1e4, we selected 10 influential parameters (key parameters) in the ORYZA_V3 model; then, we analyzed the temporal characteristics of the total sensitivity indices of these key parameters over the rice-growth period, as presented in Fig. 6. The figure shows that there is no significant difference in the results obtained with the three weather types. Parameter DVRJ has the strongest effects on all output variables at DVS ¼ 0.4. Then, the effects decrease, but there remain relatively high total sensitivity indices (more than 0.2) in following rice-growing stage. Similarly, parameter DVRI has relatively high effects on model outputs in the photo-sensitive phase, and DVRP has relatively high effects in the panicle-formation phase, whereas DVRR has relatively high effects

4. Discussion 4.1. Sensitivity analysis of selected parameters for outputs Figs. 1e4 show that the development rates (DVRJ, DVRI, DVRP, DVRR) have a strong influence on model outputs over the ricegrowth period simulated in the ORYZA_V3 model, although parameter DVRP has negligible effects on LAI at all times. Additionally, parameter RGRLMX has strong effects on all output variables in different stages, and it is one of the most influential parameters in the model (Soundharajan and Sudheer, 2013). FLV0.5 has important effects on yield, whereas FLV0.75 is non-influential

J. Tan et al. / Environmental Modelling & Software 83 (2016) 36e46

43

Fig. 4. Sensitivity indices of 20 parameters in ORYZA_V3 model calculated by Extended FAST method for rice-growth state variables at DVS ¼ 2.0: (a) outputs for cold year (1980), (b) outputs for normal year (1985), and (c) outputs for hot year (2008). First order indices are in black and interactions are in white.

at all times in rice growth; it shows a redundancy in those partition factor tables for different development stages (Richter et al., 2010). However, the leaf death coefficients (DRLV) have no effects on model outputs in this work; this finding is inconsistent with the results drawn by Soundharajan and Sudheer (2013). This discrepancy is perhaps due to the parameter range space that is not wide enough to make much of a difference in the model outputs in this study. Parameter ASLA has strong effects on LAI over the rice-growth period. Considering that the specific leaf area over rice growth varies greatly across seasons and different nitrogen-application levels (Bouman and van Laar, 2006), parameter ASLA may should not be fixed in the model. Instead, it should be supplied with different empirical values or a function for different development stages and fertilizer levels. Figs. 1e4 also show that the sensitivity indices of each parameter are distinct for different output variables and vary across development stages over the rice-growth period (see Fig. 6). This finding indicates that for complex crop models, the parameter

sensitivity analysis for other output variables is as important as that for yield (Varella, 2012; Wang et al., 2013), and some sensitive parameters would be missed if single output variable or objective function was used in the sensitivity analysis (Shin et al., 2013). 4.2. Effects of weather conditions on sensitivity results For late-season rice, the weather types have slight effects on the ranks of those influential parameters' sensitivity for outputs only in the grain-filling phase (see Figs. 4e5). In the ORYZA_V3 model, the number of spikelets is calculated as a function of accumulated dry matter in the panicle-formation phase, and parameter SPGF, and rice spikelet sterility are calculated because of cold temperatures and high temperatures over the flowering period. In 2008 (characterized as a hot year in Nanchang, China), the maximum average temperature over the flowering period did not exceed 35  C, which is the threshold temperature that would result in spikelet sterility (Bouman et al., 2001). The cold temperature over the rice-flowering period in 1980 (cold year) and 1985 (normal year) had some effects

44

J. Tan et al. / Environmental Modelling & Software 83 (2016) 36e46

Table 2 Statistic characters of outputs at different development stages over rice-growth process of late-season rice. Year

Phase

Statistic

LAI

WAGT

WST

WSO

Cold year (1980)

DVS ¼ 0.4

Mean SD CV Mean SD CV Mean SD CV Mean SD CV Mean SD CV Mean SD CV Mean SD CV Mean SD CV Mean SD CV Mean SD CV Mean SD CV Mean SD CV

2.23 1.34 0.602 4.02 2.75 0.684 5.19 1.73 0.333 1.92 0.75 0.389 2.09 1.20 0.574 4.30 3.02 0.703 5.68 1.82 0.321 2.23 0.84 0.376 2.01 1.23 0.613 3.97 2.82 0.709 5.30 1.77 0.333 2.23 0.85 0.380

1114.12 528.22 0.474 2981.33 1696.99 0.569 9181.41 2839.20 0.309 15547.79 4186.06 0.269 1040.25 472.53 0.454 3196.70 1920.25 0.601 10101.43 3084.39 0.305 16500.87 4278.13 0.259 1019.02 468.22 0.459 2917.09 1728.65 0.593 9320.34 2777.72 0.298 14930.08 3832.92 0.257

439.38 229.81 0.523 1353.00 715.06 0.529 4935.00 1406.01 0.285 4524.66 1508.73 0.333 410.91 208.49 0.507 1464.66 826.16 0.564 5473.06 1519.38 0.278 4846.61 1481.21 0.306 401.06 200.97 0.501 1324.14 742.34 0.561 5054.53 1382.99 0.274 4483.62 1269.90 0.283

e e e e e e 1269.29 326.53 0.257 8045.40 1980.46 0.246 e e e e e e 1360.30 340.17 0.250 8385.51 1908.25 0.228 e e e e e e 1240.86 285.10 0.230 7420.97 1744.16 0.235

DVS ¼ 0.65

DVS ¼ 1.0

DVS ¼ 2.0

Normal year (1985)

DVS ¼ 0.4

DVS ¼ 0.65

DVS ¼ 1.0

DVS ¼ 2.0

Hot year (2008)

DVS ¼ 0.4

DVS ¼ 0.65

DVS ¼ 1.0

DVS ¼ 2.0

Fig. 5. Sensitivity indices of 20 parameters in ORYZA_V3 model for WRR14 at maturity: (a) outputs for cold year (1980), (b) outputs for normal year (1985), and (c) outputs for hot year (2008). First order indices are in black and interactions are in white.

on spikelet sterility. Therefore, parameters WGRMX and SPGF had relatively little effect on yield in the hot year for the late-season rice in Nanchang, China, as shown in Fig. 5.

4.3. Uncertainty analysis of outputs In the basic vegetative phase (DVS ¼ 0e0.4) and the photosensitive phase (DVS ¼ 0.4e0.65), the high CV of output variables

J. Tan et al. / Environmental Modelling & Software 83 (2016) 36e46

45

Fig. 6. Total sensitivity indices of key parameters in ORYZA_V3 model varied in four development stages: (a) outputs for cold year (1980), (b) outputs for normal year (1985), and (c) outputs for hot year (2008).

is likely caused by two main aspects: On one hand, at DVS ¼ 0.4, the seedlings may either have been transplanted or still be in the seedbed, due to the effects of parameter DVRJ. For those transplanted, the values of all output variables would be diluted and reinitialized on the transplanting day in the model. Consequently, the CV has a higher value. On the other hand, when the canopy is not closed, the plants grow exponentially as a function of the temperature sum and the relative leaf-area growth rate (RGRL); this would have strong effects on the plant-growth rate in the photosensitive phase. Accordingly, the CV of the model outputs in the photo-sensitive phase are the highest over the rice-growth period. The CV of output variables varies significantly across the three weather types. Compared to normal year, the cold weather exaggerates the uncertainty of outputs in the ORYZA_V3 model derived from parameter variation, except for that at DVS ¼ 0.65. The hot weather reduced the uncertainty of outputs, although it exaggerated it at some stages. This finding may be because first, the weather conditions of some stages in the normal year may have negative effects on rice growth, which exaggerates the uncertainty of model outputs. Second, hot weather may increase stomatal closure and limit photosynthesis; to some extent, it may reduce the uncertainty of model outputs derived from parameter variation. Additionally, in the grain-filling phase, the CV of LAI increases because of leaf senescence, whereas the CV of WST increases due to the transportation of nutrients from stems to storage organs. The results suggested that influential parameters should be recalibrated with priority and higher accuracy in model calibration (Ma et al., 2011), particularly when the model is used for years characterized by colder or hotter weather than normal years. Therefore, the development rates should be accurately calculated with observed phenological data or fine-tuned for different settlements. Otherwise, inaccurate values of development rates may

result in bad simulations. The other non-influential parameters in the model, such as the minimum relative growth rate of the leaf area (RGRLMN) and those parameters in the function of specific leaf area (BSLB, CSLA, DSLA), have no effects over the rice-growth period. Therefore, by conducting parameter SA in all possible conditions, we can fix those non-influential parameters with constant values to simplify the number of model parameters without reducing the quality of the model simulation (Confalonieri et al., 2006; Richter et al., 2010). 4.4. Limitations and further study We primarily considered 20 crop parameters under potential production situation in the ORYZA_V3 model for a sensitivity analysis, with no consideration of soil parameters. Considering that the great variation of soil parameters and the inaccurate values of soil parameters resulted simulation errors (Varella, 2012), further research on sensitivity and uncertainty analysis for soil parameters would improve the understanding of the ORYZA_V3 model. 5. Conclusions In the ORYZA_V3 model, only a few parameters have important effects on model output. The development rates (DVRJ, DVRI, DVRP, DVRR), RGRLMN and FLV0.5 had strong effects on all model outputs over the rice-growth period. Only the LAI was sensitive to ASLA, and only the WSO at DVS ¼ 1.0 was sensitive to FLV1.0. Parameters WGRMX and SPGF had relatively high effects on yield but negligible effects on the other model outputs. Furthermore, those influential parameters in the ORYZA_V3 model had unequal effects on different outputs, and they had different effects at four development stages.

46

J. Tan et al. / Environmental Modelling & Software 83 (2016) 36e46

In general, the uncertainty of model outputs derived from parameter variation would be exaggerated by the cold weather. On the contrary, with the interaction effects of parameter variation and hot weather condition, the uncertainty of model outputs would be reduced to some extent. Furthermore, there is high uncertainty of outputs derived from parameter variation in the basic vegetative phase and the photo-sensitive phase under all conditions. Nevertheless, the weather conditions have negligible effects on the identification of influential parameters in the ORZYA_V3 model. Finally, the results of the sensitivity analysis calculated by the Extended FAST algorithm in this study are expected to be valuable for model calibration and model improvement. Acknowledgments This work was financially supported by the National Natural Science Foundation of China (NSFC 51579184). We are grateful for the help of Dr. Tao Li and Bouman B. A. M. of International Rice Research Institute (IRRI) in relation to the ORYZA2000 model. The valuable comments and suggestions provided by the two anonymous reviewers are also gratefully acknowledged. References Amiri, E., Rezaei, M., 2010. Evaluation of waterenitrogen schemes for rice in Iran, using ORYZA2000 model. Commun. Soil Sci. Plant Anal. 41 (20), 2459e2477. Asseng, S., Ewert, F., Rosenzweig, C., Jones, J.W., Hatfield, J.L., Ruane, A.C., Boote, K.J., Thorburn, P.J., Rotter, R.P., Cammarano, D., Brisson, N., Basso, B., Martre, P., Aggarwal, P.K., Angulo, C., Bertuzzi, P., Biernath, C., Challinor, A.J., Doltra, J., Gayler, S., Goldberg, R., Grant, R., Heng, L., Hooker, J., Hunt, L.A., Ingwersen, J., Izaurralde, R.C., Kersebaum, K.C., Muller, C., Kumar, S.N., Nendel, C., O'Leary, G., Olesen, J.E., Osborne, T.M., Palosuo, T., Priesack, E., Ripoche, D., Semenov, M.A., Shcherbak, I., Steduto, P., Stockle, C., Stratonovitch, P., Streck, T., Supit, I., Tao, F., Travasso, M., Waha, K., Wallach, D., White, J.W., Williams, J.R., Wolf, J., 2013. Uncertainty in simulating wheat yields under climate change. Nat. Clim. Change 3 (9), 827e832. Boote, K.J., Kropff, M.J., Bindraban, P.S., 2001. Physiology and modelling of traits in crop plants: implications for genetic improvement. Agric. Syst. 70 (2e3), 395e420. Bouman, B.A.M., Keulen, H.V., Laar, H.H.V., Rabbinge, R., 1996. The school of de wit crop growth simulation models a pedigree and historical overview. Agric. Syst. 52, 171e198. Bouman, B.A.M., kropff, M.J., Tuong, T.P., Wopereis, M.C.S., Berge, H.F.M.T., Laar, H.H.V., 2001. ORYZA2000: Modelling Lowland Rice. International Rice Reaserch Institute, Wageningen University and Research Centre, Los Banos, Philippines, Wageningen, Netherlands, p. 235. 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), 249e273. Cariboni, J., Gatelli, D., Liska, R., Saltelli, A., 2007. The role of sensitivity analysis in ecological modelling. Ecol. Model. 203 (1e2), 167e182. Confalonieri, R., Acutis, M., Bellocchi, G., Cerrani, I., Tarantola, S., Donatelli, M., Genovese, G., 2006. Exploratory sensitivity analysis of CropSyst WARM and WOFOST:a case study with rice biomass simulations. Italian J. Agrometeorol. 3, 17e25. Confalonieri, R., Bellocchi, G., Bregaglio, S., Donatelli, M., Acutis, M., 2010a. Comparison of sensitivity analysis techniques: a case study with the rice model WARM. Ecol. Model. 221 (16), 1897e1906. Confalonieri, R., Bellocchi, G., Tarantola, S., Acutis, M., Donatelli, M., Genovese, G., 2010b. Sensitivity analysis of the rice model WARM in Europe: exploring the effects of different locations, climates and methods of analysis on model sensitivity to crop parameters. Environ. Model. Softw. 25 (4), 479e488. Fang, Q.X., Malone, R.W., Ma, L., Jaynes, D.B., Thorp, K.R., Green, T.R., Ahuja, L.R., 2012. Modeling the effects of controlled drainage, N rate and weather on nitrate loss to subsurface drainage. Agric. Water Manag. 103, 150e161. Feng, L., Bouman, B.A.M., Tuong, T.P., Cabangon, R.J., Li, Y., Lu, G., Feng, Y., 2007.

Exploring options to grow rice using less water in northern China using a modelling approach. Agric. Water Manag. 88 (1e3), 1e13. Gallagher, M., Doherty, J., 2007. Parameter estimation and uncertainty analysis for a watershed model. Environ. Model. Softw. 22 (7), 1000e1020. Lamboni, M., Makowski, D., Lehuger, S., Gabrielle, B., Monod, H., 2009. Multivariate global sensitivity analysis for dynamic crop models. Field Crops Res. 113 (3), 312e320. Li, T., Raman, A.K., Iii, M.Marcaida, Kumar, A., Angeles, O., Radanielson, A.M., 2013. Simulation of genotype performances across a larger number of environments for rice breeding using ORYZA2000. Field Crops Res. 149 (0), 312e321. Ma, L., Ahuja, L.R., Saseendran, S.A., Malone, R.W., Green, T.R., Nolan, B.T., Bartling, P.N.S., Flerchinger, G.N., Boote, K.J., Hoogenboom, G., 2011. In: Ahuja, L.R., Ma, L. (Eds.), A Protocol for Parameterization and Calibration of RZWQM2 in Field Research. Methods of Introducing System Models into Agricultural Research. American Society of Agronomy, Crop Science Society of America, Soil Science Society of America, Madison, WI, pp. 1e64. Ma, G., Huang, J., Wu, W., Fan, J., Zou, J., Wu, S., 2013. Assimilation of MODIS-LAI into the WOFOST model for forecasting regional winter wheat yield. Math. Comput. Model. 58 (3e4), 634e643. Makowski, D., Naud, C., Jeuffroy, M.-H., Barbottin, A., Monod, H., 2006. Global sensitivity analysis for calculating the contribution of genetic parameters to the variance of crop model prediction. Reliab. Eng. Syst. Saf. 91 (10e11), 1142e1147. Richter, G.M., Acutis, M., Trevisiol, P., Latiri, K., Confalonieri, R., 2010. Sensitivity analysis for a complex crop model applied to Durum wheat in the Mediterranean. Eur. J. Agron. 32 (2), 127e136. Saltelli, A., Annoni, P., 2010. How to avoid a perfunctory sensitivity analysis. Environ. Model. Softw. 25 (12), 1508e1517. Saltelli, A., Ratto, M., Andres, T., Campolongo, F., Cariboni, J., Gatelli, D., Saisana, M., Tarantola, S., 2008. Variance-Based Methods. Global Sensitivity Analysis. In: The Primer. John Wiley & Sons, Ltd, pp. 155e182. Saltelli, A., Tarantol, S., Campolongo, F., Ratto, M., 2004. Sensitivity Analysis in Practice. A Guide to Assessing Scientific Models. John Wiley & Sons. Saltelli, A., Tarantol, S., K.P.S, C., 1999. A quantitative model-independent method for global sensitivity analysis of model output. Technometrics 41 (1), 39e56. Shin, M.-J., Guillaume, J.H.A., Croke, B.F.W., Jakeman, A.J., 2013. Addressing ten questions about conceptual rainfallerunoff models with global sensitivity analyses in R. J. Hydrol. 503 (0), 135e152. Soundharajan, B., Sudheer, K.P., 2013. Sensitivity analysis and auto-calibration of ORYZA2000 using simulation-optimization framework. Paddy Water Environ. 11 (1e4), 59e71. Tremblay, M., Wallach, D., 2004. Comparison of parameter estimation methods for crop models. Agronomie 24 (6e7), 351e365. Vanuytrecht, E., Raes, D., Willems, P., 2014. Global sensitivity analysis of yield output from the water productivity model. Environ. Model. Softw. 51, 323e332. Varella, H., 2012. Global sensitivity analysis for choosing the main soil parameters of a crop model to be determined. Agric. Sci. 03 (07), 949e961. Vazquez-Cruz, M.A., Guzman-Cruz, R., Lopez-Cruz, I.L., Cornejo-Perez, O., TorresPacheco, I., Guevara-Gonzalez, R.G., 2014. Global sensitivity analysis by means of EFAST and Sobol’ methods and calibration of reduced state-variable TOMGRO model using genetic algorithms. Comput. Electron. Agric. 100, 1e12. Wallach, D., Goffinet, B., Bergez, J.-E., Debaeke, P., Leenhardt, D., Aubertot, J.-N., 2001. Parameter estimation for crop models. Agron. J. 93 (4), 757e766. Wallach, D., Makowski, D., Jones, J.W., Brun, F., 2014. In: Brun, D.W.M.W.J. (Ed.), Chapter 6-Parameter Estimation with Classical Methods (Model Calibration). Working with Dynamic Crop Models, second ed. Academic Press, San Diego, pp. 205e276. Wang, J., Li, X., Lu, L., Fang, F., 2013. Parameter sensitivity analysis of crop growth models based on the extended Fourier Amplitude Sensitivity Test method. Environ. Model. Softw. 48, 171e182. Wang, W., Yu, Z., Zhang, W., Shao, Q., Zhang, Y., Luo, Y., Jiao, X., Xu, J., 2014. Responses of rice yield, irrigation water requirement and water use efficiency to climate change in China: historical simulation and future projections. Agric. Water Manag. 146 (0), 249e261. Wu, L., Liu, X., Wang, P., Zhou, B., Liu, M., Li, X., 2013. The assimilation of spectral sensing and the WOFOST model for the dynamic simulation of cadmium accumulation in rice tissues. Int. J. Appl. Earth Obs. Geoinform. 25, 66e75. Yadav, S., Li, T., Humphreys, E., Gill, G., Kukal, S.S., 2011. Evaluation and application of ORYZA2000 for irrigation scheduling of puddled transplanted rice in north west India. Field Crops Res. 122 (2), 104e117. Zhao, G., Bryan, B.A., Song, X., 2014. Sensitivity and uncertainty analysis of the APSIM-wheat model: interactions between cultivar, environmental, and management parameters. Ecol. Model. 279, 1e11.