Evaluation of parameterization strategies for rice modelling

In a previous study, the crop simulator CropSyst was evaluated against crop data collected on rice varieties grown in Northern Italy. The need of model re-parameterization became apparent after investigating new field data, as inconsistencies in the simulation of leaf area index emerged for Indica-type varieties. Key parameters (specific leaf area, stem-leaf partition, extinction coefficient, light-to-biomass conversion efficiency) derived from field measurements (respectively, 27 m2 kg-1, 3.6 m2 kg-1, 0.53, 3.2 g MJ-1) considerably differed from those previously obtained via calibration (39 m2 kg-1, 1.5 m2 kg-1, 0.50, 3.0 g MJ-1). Such new parameters are informative for suitable modelling of rice systems. The agreement between simulated and observed above ground biomasses was similar with both parameter sets: average general standard deviation = 25% (previous) and 26% (new); average modelling efficiency = 0.90 (previous) and 0.87 (new). Such comparisons demonstrate as the accumulation of aerial biomass in crop models can be depicted in different ways and reasonable estimations can be achieved by different pathways (not all acceptable). A check on parameters like the one performed here (field measurements versus calibrated parameters) is worth to give protection against spurious conclusions while indicating whether the parameterization is conceptually consistent and related to reality. Additional key words: above ground biomass, CropSyst, leaf area index, model evaluation, model parameterization, Oryza sativa L.


Introduction
Process-based simulation models for estimating crop growth and development were built and revised over the past decades (e.g.volume 18 of the European Journal of Agronomy: issues 1-2, 2002; issues 3-4, 2003).Typically, such models estimate crop growth and its development by mathematical representations of bio-physical processes, and differ from functionoriented models, where only the empirical correlations between inputs and outputs are identified (e.g.Kandiannan et al., 2002).Process-based crop models give estimates of yield and harvest time based upon soil characteristics and weather dynamics, under different management scenarios.For this reason, modelling approaches have been used in a wide variety of applications at field, farm and regional scale (see, for instance, proceedings from Farming Systems Design 2007 -International symposium on Methodologies for Integrated Analysis of Farm Production Systems, Catania, Italy, http://www.iemss.org/farmsys07).
The production of high-power computers and the increase of scientific knowledge are allowing simulation models to exponentially increase their complexity.Nowadays simulation models are constituted by hundreds of algorithms and, in many cases, this "growing" process is not without risks.Monteith (1996), in an article titled "The quest for balance in crop modelling", gave a provocative idea of to what such risks refer.According to the author, in fact, one of the most common dangers is that parameter values are given that do not really correspond to the morphological and physiological features of the crop species (or cultivar) under study.With crop models this may occur, for instance, when a set of parameters are assigned values (generally obtained by calibration against actual datasets) to match marketable yield only, not considering the time-course of basic variables (e.g.Banterng et al., 2003).Quantitative information regarding biomass accumulation and leaf area evolution remains therefore limited or inconclusive.
In the case of complex models, the vector of all parameters can be represented by high dimensional hyperspace, with the parameters placed on different axes (Acutis and Confalonieri, 2006).Theoretically, a crop is characterized (both morphologically and physiologically) by a single set of parameter values, whilst other combinations (also likely to give plausible estimates under ordinary conditions) must be considered as not logically acceptable (Sinclair and Seligman, 2000).Under prevailing conditions, experimenting with different combinations of parameter values may return the same output value (distance from the origin in the multi-dimensional space).When outlying conditions occur (for example, as a result of anomalous meteorological events) some combinations of parameters may turn into different, yet inconsistent outputs.These concepts, in appearance purely speculative and theoretical, describe well what happened when the crop model CropSyst (Stöckle et al., 2003), calibrated for three groups of rice (Oryza sativa L.) varieties grown in Northern Italy (Indica, early and mediumlate maturity Japonica) by Confalonieri and Bocchi (2005), was tested on a new Indica dataset (Confalonieri et al., 2006) including morphological and physiological measurements not present in previous datasets.Indica and Japonica types refer to the two main groups of rice cultivars, which differ for morpho-physiological aspects and market destination.The model provided reliable estimate of above ground biomass (AGB) but revealed some inconsistencies on simulated leaf area index (LAI) values.This called for more accurate modelling to support the analysis of rice systems.This is particularly important in the European Union, where models are largely used for the purposes of planning, management, and monitoring (http://mars.jrc.it/mars/Bulletins-Publications),especially in the Mediterranean area where Italy and Spain are the main producers (e.g.Confalonieri and Bocchi, 2005;Castejón-Muñoz et al., 2007).A further analysis on the new experimental data in relation to the value of some parameters suggested it was worth re-parameterizing the model using the new measured data wherever possible.
The objectives of this study were (i) re-parameterizing CropSyst for Indica varieties according to key model parameters determined using field data; (ii) compare the CropSyst performances obtained with the new and the old parameterizations.

Experimental data
The data used in this study are from experiments described in detail by Confalonieri and Bocchi (2005) and Confalonieri et al. (2006), referring to rice crop trials (Indica-type varieties) grown between 1999 and 2004 on loam soils at three sites across the Po Valley (Northern Italy): Velezzo Lomellina (lat.45º 9' N, long.8º 44' E) in 1999 (cv.Thaibonnet, sowing date April 1), Vignate (lat.45º 29' N, long.9º 22' E) in 2002 (cv.Sillaro, sowing date April 29), Opera (lat.45º 22', long.9º 12' E) in 2002 (cv.Thaibonnet, sowing date: April 29), and Opera in 2004 (cv.Gladio, sowing date: May 24).The Po valley is a large basin characterized by a transitional climate between the Mediterranean climate dominated by anti-cyclonic patterns and the Central European climate, dominated by the oceanic influence of westerly circulations.For the study-area, the rainfall variability is limited to 650-800 mm year -1 (with about 75 rainy days per year), in contrast to about 950-1100 mm year -1 of reference crop evapotranspiration.The mean annual air temperature is about 12.5-13.5°C, with minimum and maximum values registered in January-February and July-August respectively.
Datasets were selected from treatments where nitrogen was not a limiting factor.Water management was comparable in all the experiments with fields flooded before sowing.The fields kept flooded for the whole growing season, maintaining the water level at about 0.10 m.Management always allowed the prevention of water and nutrient stresses, and full control of weeds, pests, and diseases.Measured variables for all the datasets are described by Confalonieri and Bocchi (2005) and Confalonieri et al. (2006) and are not reproduced here.Moreover, during the experiment carried out in Opera in 2004, measurements of light transmittance into the canopy were carried out using a LAI-2000 (LI-COR, Lincoln, NE) (Stroppiana et al., 2006).

Model parameterization and evaluation
CropSyst (Stöckle et al., 2003) version 3.02.23 (January 8, 2002) was used.Model algorithms are described in detail in the CropSyst manual (http:// www.bsyse.wsu.edu/cropsyst) and a brief description of the main routines involved with growth and devel-opment can be found in Confalonieri and Bocchi (2005) and Bechini et al. (2006).
In this study, only potential production levels (van Ittersum and Rabbinge, 1997) were simulated.The dataset of Opera 2004, the one which revealed the inconsistencies in the parameters' set provided by Confalonieri and Bocchi (2005), was used for model parameterization while the validation of model outputs was performed on the datasets of Velezzo L. 1999, Opera 2002, and Vignate 2002.During the experiment carried out in Opera in 2004, the following parameters were determined: specific leaf area (SLA, m 2 kg -1 ), steam/leaf partition coefficient (SLP, m 2 kg -1 ), photosynthetically active radiation (PAR), light extinction coefficient (k), light to above ground biomass conversion (i.e.radiation use efficiency; LtBC, g MJ -1 ), and above ground biomass-transpiration coefficient (BTR, kPa kg m -3 ).Other parameters were assigned the same values used by Confalonieri and Bocchi (2005).In the latter, crop parameters were obtained using a trial-and-error calibration procedure.In this work, this method was applied to estimate BTR, while the other parameters were derived from measured data.
Differently from other crop simulation models (e.g.those from the Wageningen School, van Ittersum et al., 2003), CropSyst assumes a constant SLA for the whole growing season.The value of this parameter was set as the average of measurements carried out during the early growth stages (Stöckle C.O., Washington State University, personal communication).The obtained constant SLA and measured AGB and LAI values were used to determine SLP by re-arranging equation [1] (Stöckle et al., 2003) and averaging the values corresponding to all pairs LAI-AGB. [1] Parameter k (also in this case CropSyst uses a constant value) was determined by using LAI-2000 light transmittance measurements, following the procedure proposed by Dingkuhn et al. (1999) and averaging the values obtained in different moments during the crop cycle.LtBC was considered equal to the ratio between periodic AGB samplings and the corresponding cumulated absorbed PAR (Sinclair and Muchow, 1999) in the period between emergence and anthesis.Since LtBC corresponds to the radiation use efficiency under nonlimiting thermal conditions, an inverse of the CropSyst factor for temperature limitation (a linear function between base and optimal temperatures) was applied to the measured AGB to obtain AGB values referring to fully potential (although virtual) conditions.The obtained value is 12% higher than the one obtained without considering the temperature effect.Further details on the methods used for parameter estimation are provided by Boschetti et al. (2006).
The agreement between measured and simulated values was evaluated by using IRENE_DLL (Fila et al., 2003), a dedicated tool for evaluation of model performance.In particular, the general standard deviation (GSD, %; Jørgensen et al., 1986; 0 ÷ +∞, optimum = 0%) was used for evaluating the average error of model estimates and the modelling efficiency (EF, Greenwood et al., 1985; -∞ ÷ 1; optimum = 1; if positive, it indicates that the model is a better estimator than the average of measured values) to evaluate the agreement between measured and simulated trends.In order to compare the results obtained from using the two sets ("old" and "new") of model parameters, a fuzzy logic based multiple-statistics assessment system (Bellocchi et al., 2002) was used to aggregate GSD and EF in order to obtain a synthetic indicator of model performance.For the computation of the fuzzy-based aggregated indicator (FAI), the favourable and unfavourable thresholds for GSD were set, respectively, to 20% and 40%.The corresponding thresholds for EF were set to 0.9 and 0.4.The same relative importance was assigned to the two indices (expert weight equal to 0.5).

Results and discussion
Modified parameter values compared to those used by Confalonieri and Bocchi (2005) are in Table 1.The new values of SLA, k, and LtBC shown in the table were derived through field data measured by Boschetti et al. (2006).The new value of SLA (27 m 2 kg -1 ), set to the average of the first three measurements (22, 33, and 37 days after sowing), falls within the range of values reported by other authors (e.g.Dingkuhn et al., 1998).The value is lower than that used by Confalonieri and Bocchi (2005) for Indica-type varieties (39 m 2 kg -1 ) and approximates the value used by the same authors for simulation of Japonica-type varieties (~28 m 2 kg -1 ).Also the new value of SLP (3.60 m -2 kg -1 ) is similar to that used in the original parameterization for Japonica-type varieties (3.75 m -2 kg -1 ).The new combination of parameters addresses the main differences between Indica and Japonica-types varieties to the interception and use of radiation.Values of both k and LtBC reflect the tendency of Indica-type varieties recently selected (semi-dwarf, high yielding) to improve the efficiency in the use of solar radiation.Therefore, the value of k used in the original parameterization (0.5) can be considered suitable for Japonica-type but a higher value (0.53) is now used for Indica-type varieties.The new value of k is consistent with those measured in the first part of the crop cycle by Kiniry et al. (2001) for rice cvs.Cocodrie and Cypress.These cultivars are close derivations from L202 (renamed Thaibonnet in Europe) and are therefore related to cv.Gladio used in this study (which is also related to L202 itself; Spada et al., 2004).High values of k for these types of cultivars were also reported by Casanova et al. (2000) and Campbell et al. (2001).These authors indicated, respectively, average k values of 0.54 and 0.68.Our estimation of LtBC (3.20 g AGB MJ -1 intercepted PAR) is higher than those commonly found in the literature (e.g.2.39 g MJ -1 , Kiniry et al., 1989) probably because our procedure led to obtain a potential value which is not affected by sub-optimal thermal conditions.Eliminating the air temperature effect from the procedure for LtBC estimation, the value computed using our data is 2.8 g MJ -1 , close to that reported by other authors (e.g.Horie and Sakuratani, 1985).The value used to run new simulations (Table 1) is anyway lower than the value of 3.52 g MJ -1 estimated by Campbell et al. (2001).
The results obtained from the new and the old simulations were compared, as shown in Fig. 1 and Confalonieri and Bocchi (2005) in the column "Old value"; measured (or derived from measurements) parameters generated in this study in the column "New value" 2. Systematic model under-estimation after heading stage (Fig. 1) already affected the results obtained with Indica-type varieties in the simulation carried out by Confalonieri and Bocchi (2005).Such behaviour is probably due to high crop response to nitrogen fertilization at panicle initiation.The consequent, rapid increase in the concentration of plant nitrogen (PNC), which is mainly transferred to leaves in case of high uptake rates, could lead to an increased synthesis of proteins directly involved in photosynthesis (e.g.RuBisCO, Confalonieri et al., 2006).A dynamic value of LtBC during the period after panicle initiation as a function of PNC should be tested to solve under-estimation problems.Simulation of LAI was tested against the dataset collected at Opera in 2004, which included LAI measurements.By using the new parameterization, close agreement was observed between simulated and measured values of LAI (GSD=39%, EF=0.72), while the combination of parameters proposed for Indica-type varieties by Confalonieri and Bocchi (2005) led to high over-estimations (GSD almost 100%, EF negative).This could be considered as an argumentation suggesting that, whereas a calibration work has to be carried out because crop parameter values are not available, more reliable results are expected to come out by adjusting the model fitting against more than one state variable.
Fig. 1 (internal graph) shows the comparison between the AGB values simulated with the two sets of model parameters.It is worth noting the good agreement between them (R 2 >0.99): the two parameters sets (which strongly differ from each other for important parameters) led to almost the same results for AGB accumulation.
GSD and EF are similar with both parameterizations (Table 2).For the datasets of Opera 2004, Vignate 2002 and Opera 2002, the old combination of crop parameters performed slightly better, while for the dataset of Velezzo L. 1999, the new parameterization shows a better agreement between simulated values and observations.A possible explanation is related to the different tendency of the model in producing photosynthetic material at early stages with the two parameterizations.The old parameterization was characterized by a higher production of LAI units immediately after sowing (higher SLA), while the new one reaches high final production levels mainly because of a higher photosynthetic efficiency.In the Velezzo L. 1999 experiment, rice was sown earlier than in the other cases.In this condition, the old parameterization led to overestimate the AGB measured in the first samplings, since the model compensated the low photosynthetic efficiency due to the sub-optimal thermal and irradiative conditions (which characterize rice growth in early spring) with higher interception.In all cases, GSD and EF gave satisfactory values, both being always close to their optima.Fuzzybased aggregated index showed good discrimination ability, amplifying the relative differences between the  Confalonieri and Bocchi (2005) and the new parameterization.GSD: general standard deviation (Jørgensen et al., 1986), EF: modelling efficiency (Greenwood et al., 1985), FAI: fuzzy-based aggregated indicator (see text for details).Opera 2004 dataset was used for calibration, the others for validation results obtained with the two combinations of parameters.The values of GSD and EF were arithmetically averaged for the four datasets, in order to obtain global GSD and EF values independent of the number of observations that compose the different datasets.The FAI calculated aggregating these indices were 0.07 and 0.10 for the old and the new combination of parameters, respectively.This allows considering the overall performances of the two parameterizations for AGB simulation as very similar in the explored conditions, although the new one fits better the real morpho-physiological plant features whereas the old one was reproducing plant behaviour by means of compensating errors in the parameters' values.

Conclusions
Plant measurements present a series of challenges due to the subjectivity of an observer or to a definition that is not unambiguously applied in the field.This emphasizes the importance of meta-data associated with original observations and development of model parameters.During the period of this research, there was a great deal of change in our fundamental understanding involving CropSyst-based modelling of rice growth.The understanding that finally emerged has been presented in this paper, which contrasts the parameterization used in a previous simulation-based paper.This contribution draws on the application of the CropSyst crop parameters calibrated by Confalonieri and Bocchi (2005) for Indica-type rice varieties to experimental data collected in Northern Italy.The main model parameters (specific leaf area, stem-leaf partition, extinction coefficient and light-to-biomass conversion coefficient), directly measured over a new dataset, were different from those calibrated in the previous work.However, despite these differences, above ground biomass values simulated with both combinations of parameters accurately fitted observations and resulted in model estimates that where close to each other, but leaf area index was better reproduced by the new parameter set.Baseline values were provided in this study for parameters of rice crops grown under well-supplied conditions.This work brings new awareness that model parameterization should be mostly based on measured data (rather than using calibration), so preventing crop modelling from degrading towards a mere fitting exercise not reproducing the actual system's behaviour.From a practical point of view, these results represent a sort of rec-tification to the parameter values proposed in a previous study.Although the database used is not extensive and a single dataset was used for model parameterization, the satisfying results obtained for the three datasets used for validation lead us feeling confident about the proposed set of parameters.Potential CropSyst users (and virtually all crop modellers pointing their attention on rice) are now provided with a more reliable set of parameters for modelling rice.However, these data are presented and discussed here not only to make available the "correct" set of CropSyst parameters for Indica-type varieties.Parameterization is not a one-off event or a "once-andfor-all" activity, and the paper lays emphasis on the importance of continually revisiting data collected until the researcher is persuaded that the parameters used to describe the system under study are a truthful and accurate reflection of the data.This is recommended for any modelling study to ensure a progressive evolution of the models by time.With the new set of parameters, a better knowledge and a more reliable modelling of rice systems are indeed obtained for use in Italy and other rice systems in the Mediterranean area.The authors hope this study will draw once more the attention of the modellers' community on the caution necessary when complex tools such as simulation models are used.There are strong limits on the role of re-parameterization of existing models for reducing estimation uncertainty, but work towards this goal can be stimulated by identifying key controls of agro-system processes and their sensitivities to driving inputs.

)Figure 1 .
Figure 1.Comparison between simulated and measured above ground biomass (AGB) values.Circles refer to the Opera 2004 dataset (calibration), triangles, squares and diamonds, to the datasets of Vignate 2002, Velezzo L. 1999 and Opera 2002 (validation), respectively.White and black symbols correspond to the results obtained with the new and the old parameterization, respectively.Regression lines are shown, and compared with the 1:1 line.Internal graph: comparison between above ground biomass values simulated by the new and old sets of parameters (continuous line is the regression line, dotted line is the 1:1 line).

Table 1 .
in Table Crop model parameters for Indica-type varieties.Parameters calibrated by

Table 2 .
Performance indices between measured above ground biomass data and those simulated using both the parameterization from