Understanding Crop Response to Climate Variability with Complex Agroecosystem Models

Newlands, Nathaniel K.; Espino-Hernández, Gabriela; Erickson, R. Scott

doi:https://doi.org/10.1155/2012/756242

International Journal of Ecology

On this page

Abstract Introduction Materials and Methods Results Discussion Acknowledgments References Copyright Related Articles

Research Article | Open Access

Volume 2012 | Article ID 756242 | https://doi.org/10.1155/2012/756242

Understanding Crop Response to Climate Variability with Complex Agroecosystem Models

Nathaniel K. Newlands,¹Gabriela Espino-Hernández,²and R. Scott Erickson¹

Academic Editor: Pavlos Kassomenos

Received30 Sept 2011

Accepted03 Nov 2011

Published13 Feb 2012

Abstract

Agroecosystem modeling studies often rely on relatively short time-series historical records for training/tuning empirical parameters and to predict long-term variation in crop production associated with trends in climate and hydrological forcing. While ecosystem models may exhibit similar prediction skill in validation studies, their sensitivity to climate variability can differ significantly. Such discrepancy often arises due to the need to tradeoff model complexity with data availability. We examine the sensitivity in predicting spring wheat crop productivity across agricultural sites with differing soil and climate conditions where long-term agronomic and climate records are available. We report significant changes in the model sensitivity accompanying changing climatic regime. If not corrected for, this can lead to substantial predictive error when simulating across time and space. Our findings lend further support for a hierarchical (componentwise) approach for reducing model complexity and improving prediction skill.

1. Introduction

Understanding natural processes within terrestrial ecosystems can be made easier by using ecosystem models [1]. Such system-level models are used in research efforts to understand plant biomass production in both agricultural and nonagricultural ecosystems. Ecosystem models estimate production as a function of deterministic growth processes and stochastic environmental variables. They enable study of numerous physical, environmental, and other variables as well as the processes and interrelationships that connect them. Although the models have tended to be inclusive with respect to the variables considered, they can be difficult to parameterize correctly [2] and often become more complex than needed for certain applications. This may result in poor predictive performance due to the compounding of parameter uncertainties within the model [3]. While ecosystem models may exhibit similar prediction skill in validation studies, the sensitivity of their long-term predictions to climate variability can differ significantly [4–6]. Such discrepancies often relate to an outstanding dilemma in ecosystem science—how best to tradeoff model (i.e., parameter and structural) complexity with data availability. Modeling studies often rely on relatively short time-series historical records for training that involves tuning parameter values, to predict long-term variation in crop production associated with trends/variability in climate and underlying hydrological forcing. In particular, the robustness of ecosystem model predictions of net primary production (NPP) and soil water content (SWC) relies on ability to handle environmental variables that are strongly coupled and vary nonlinearly. Model sensitivity analysis can be carried out to quantify the effects of given input parameters on certain model outputs. Local sensitivity analysis varies the inputs one at a time while keeping the remaining parameters fixed [7], but is not appropriate in situations where the model is nonlinear [8]. Global sensitivity analysis allows inputs to vary simultaneously within given ranges [9] and has recently become more favourable due to advances in computational power [10]. Sensitivity analyses can identify relationships between model components, as well as the lack of such relationships and can therefore guide efforts in model simplification. A variety of methods has been devised, each with their own drawbacks and benefits. The range of available techniques include: Differential (e.g., UCODE [11]), nominal range, pearson/spearman correlation, response surface, mutual information index, classification and regression tree (CART), and variance-based methods like Sobol’s method and Fourier Amplitude Sensitivity Test (FAST) [12].

While the consideration of a broad set of explanatory variables is important and beneficial for calibrating and propagating error through system-level models, resulting overparameterization can lead to problems in obtaining robust solutions and reliable interpretation and insights of model simulation output. Findings from reported studies on the sensitivity of ecosystem models in both forestry and agricultural (i.e., cropland) applications reveal that the majority of total output variability is attributed to only a few leading or “relevant” parameters [13]. Also, by reducing uncertainty in leading input parameters, substantial gains in overall model precision can be attained. Confalonieri [14] identified that 3 leading parameters out of 15 for CropSyst, and 3 out of 24 for the WOFOST crop models explained >90% of variance in model simulation output. In a companion study, Confalonieri et al. [15] further highlight the crucial importance of exploring site-climate combinations/interactions as a prerequisite to evaluate novel crop-modeling approaches and/or the application of known modeling solutions not previously explored. Their study examined the sensitivity of the WARM rice model across lowland littoral, semicontinental, and high-continental climatic regimes in Europe, for which air temperature and site location (distance from coast) were considered the most important driving variables for rice grown under flooded conditions. For the WARM rice model, radiation-use efficiency (RUE), optimum temperature for growth (), initial leaf area index (LAI_init), partition coefficient to leaf at early stages (RipL0), and the extinction coefficient for solar radiation () were the most sensitive parameters to changes in site and climatic regimes, respectively. Lamboni et al. [16] devised a multivariate sensitivity analysis method with time-series input, comparing it to sequential-based assumptions. They explored the sensitivity of the CERES-EGC model simulating nitrous oxide (N₂O) emission from cropland, and a model of wheat-biomass (called WWDM), concluding overall that the estimation of parameters associated with highest sensitivity indices led to strong reduction in model prediction error. Congruent with Confalonieri et al.’s findings [15], their analysis identified that radiation-use efficiency, coefficient of LAI increase, and air temperature had the leading impact on winter wheat biomass under conventional management and a haplic Calcisol soil, respectively [16]. Three parameters, namely, the fraction of denitrified nitrogen, potential denitrification rate, and half-saturation constant of denitrification had a leading effect on N₂O emission predicted by CERES-EGC of a total of 15 model parameters. Richter et al. [17] devised a model to simulate water-limited crop growth to obtain several agroecological indicators. Their complex ecosystem model linked micrometeorology hydrology with crop development and growth for hilly terrain (via the consideration of spatial soil type variability and topography). They parameterised their model for durum wheat grown on two Mediterranean sites (i.e., a silty-clay-loam in Volturino, Italy versus a drought-prone sandy loam in Nabeul, Tunisia—both exhibiting a skewed yield distribution). The most sensitive process “domains” within both of the sites were phenological development, leaf area dynamics, and light interception, followed by morphological (allocation) and ecophysiological parameters. However, model sensitivity due to intersite variability in yield was best explained by temperature, rainfall and soil water content. Winter wheat yield variability has also been explored using the AZODYN dynamic agroecosystem model by Makowski et al. [18] using observational data from crop experiments in clay-loam soil in Grignon, northern France, conducted during 1991–2002. Mineral nitrogen and climate (warm temperature and precipitation) are considered not limiting for yield at this site. They focused on identifying the relative extent to which 13 genetic parameters (of 69 total model parameters) that track genotypic variation contributed to sensitivity in the model output. Maximum yield, ratio of intercepted to incident radiation, and ratio of leaf area index to critical nitrogen (only 3 parameters) consistently explained all the observed yield variation, while for grain protein content simulated output, ratio of total to above ground nitrogen, and fraction of remobilized nitrogen. Varella et al. [19] have applied an extended version of the FAST technique to explore the sensitivity of a highly complex model (STICS wheat) having >200 parameters under different soil, climatic, and crop conditions. They examined model sensitivity in Cambry, northern France, across a range of climatic regimes, soil depths and crops (i.e., sugar-beet and peas). Overall, their study identified that climate and soil depth had a leading effect on yield at harvest, leaf area index, and nitrogen absorbed, with initial soil water content most sensitive in a drier climate regime.

The Biome-BioGeoChemical Cycles (BGCs) model simulates the mass-balanced dynamics of carbon, water, nitrogen, and energy in forest and nonforest terrestrial ecosystems [20–22]. The model operates on a daily time step and is a generalization of the Forest-BGC computational model [23, 24]. Sensitivity analyses have been performed on the Biome-BGC model in nonagricultural ecosystems. For example, White et al. [25] conducted a local sensitivity analysis of net primary production (NPP) of forest, grassland, and shrub biomes to variations in ecophysiological variables such as turnover, C : N ratio, allocation of plant material, leaf morphology, and others. Their study revealed that NPP is strongly affected by variations in leaf and fine root C : N in all biomes. In woody biomes, NPP was strongly affected by percent leaf nitrogen in rubisco, maximum stomatal conductance, and specific leaf area. In nonwoody biomes, NPP was sensitive to fire mortality and litter quality. None of the parameters examined demonstrated strong interaction effects. Biome-BGC has been adapted for managed forests in central Europe and performed a local sensitivity analysis to determine the effects of site and ecophysiological parameters on NPP and carbon pools in biomass, litter and soil [26]. Their analysis revealed that all the tested variables were highly sensitive to site parameters such as total precipitation, ambient CO₂, and nitrogen input, as well as ecophysiological parameters such as leaf and fine root C : N ratio, specific leaf area, maximum stomatal conductance, fire mortality, PLNR, and others. Furthermore, such findings indicate that NPP can respond with greater sensitivity to a 20% increase in precipitation or a doubling of atmospheric CO₂ than to a 2°C temperature increase [26].

In summary, previous ecosystem modeling studies indicate that agro- and forest-ecosystem model predictions vary across a range of complexity, structural, empirical, and nonlinearity assumptions. Within agroecosystems, response sensitivity has been explored across a range of different crops (i.e., rice, sugar-beet, peas, wheat) revealing consistent findings that within-site variation is most sensitive to (in order of strength): light-interception, leaf-area dynamics, and phenological development/staging. Such ranking is reasonably consistent also across different wheat genetic varieties and changes in phenological development (i.e., fall-planted hard/soft red winter wheat, spring-planted durum, and spring wheat types). In contrast, between sites, air temperature, rainfall, solar radiation (via radiation use-efficiency/light-interception), and soil depth are all consistently reported to drive overall sensitivity in crop response. Yet despite these general insights, considerable questions still exist regarding statistically how the sensitivity of environmental variables varies in time and space. Only with sufficiently long time-series of climate data/record, and sufficient knowledge of how model parameter sensitivity varies under changing soil and climate conditions can the best corresponding sets of environmental variables and model parameters be identified in predicting crop response.

Wheat is a major crop grown across North America and around the world. Canada and USA are major world producers of wheat, with 2008 production estimates exceeding 28 and 68 M tonnes, respectively (United States Department of Agriculture, USDA). The majority of wheat grown in North America is of the spring-seeded type. Recent questions regarding climate change have raised concerns about the security of this and other important food supplies. It is presumed that changes in global and regional climates may cause production of food crops to fluctuate, with the potential to adversely affect world supplies and thereby cause societal problems. In this paper, we report findings from a global sensitivity analysis of spring wheat crop response to annual climate variability (i.e., across seasons) using the Biome-BGC complex model. Although previous studies have analyzed local sensitivity of NPP to certain parameters, no global sensitivity analysis of the Biome-BGC model has been performed to date. In addition, Biome-BGC has been applied to agricultural/cropland ecosystems in only a few instances [27, 28]. We chose to apply the Sobol variance-based method, as it has been widely applied, thereby enabling our findings to be compared with previous work. Also, the method has been coded, verified, and made available as a software package (Simlab) that facilitated its linkage with ecosystem model source code. We examine the predictions of NPP and SWC from a complex ecosystem model across agricultural sites with differing soil and climate conditions where long-term agronomic and climate records are available. Soil water content (SWC) and net primary productivity (NPP) are defined as the quantity of water on the soil, and the net amount of available carbon for crop allocation to leaves, stem, roots, and growth per unit of space and time, respectively. We parameterize an ecosystem model for cropland, linking it with a global sensitivity analysis routines that can be also further used and applied by other agricultural scientists and ecologists. The long-term goal of our work is to better inform how to obtain better prediction skill in modeling the response of crops to climate variability in time and space.

2. Materials and Methods

The Biome-BioGeoChemical Cycles (BGCs) model simulates the mass-balanced dynamics of carbon, water, nitrogen, and energy in forest and non-forest terrestrial ecosystems. The model operates on a daily time step and is a generalization of the Forest-BGC computational model. Biome-BGC considers not only the simple effects of input parameters, but also interactions among them. Benchmarked, public-released code for the model (Version 4.1.1, 2005) was utilized in this study, focusing on its 10 key site-specific, and 53 vegetation-linked input parameters (Table 1). The total number of parameters in the latest model version 4.2 is 67 to simulate the transformation of carbon and nitrogen between four different litter and soil pools. In particular, the ecophysiological parameters in the model are used to distinguish different biomes and are described in detail elsewhere [24, 25]. Initial runs of the model are typically conducted to bring soil organic matter into dynamic equilibrium prior to conducting simulation runs-termed “spin-up.” While guidance on the specification of default parameter values is available, small uncertainties propagate in the model to generate large net variability in simulation output. For this reason, fine-tuning of parameter values for various applications is required. Furthermore, what portion of such variability is attributable to the model parameter values as opposed to model structural assumptions is often a very confounded question.

The modeling framework is detailed below and is adapted from the Gaussian-process metamodeling methodology of Marrel et al. [29] for the specific agricultural application in this study and its extension to the following more realistic and less restrictive statistical assumptions: (i) Gaussian distributed temperature and non-Gaussian, Weibull-distributed precipitation inputs, with consideration of heavy-tailed distribution extremes; (ii) adaptation of stationary response function (one-degree polynomial function) to nonstationary, multivariate response surface; (iii) exponential-type local spatial correlation function to non-local extended spatial correlation. The study considers a multivariate crop response function, comprising separable terms each with their specific functional description and association, as a simplified representation of a complex ecosystem model. Statistical information and probability distributions associated with these contributions to crop response to climate variability were obtained via global sensitivity analysis of the ecosystem model and statistical information contained in measured individual, total, and interaction indices related to variances and their relative strength. A simplified statistical crop response function that approximates the crop response simulated by a complex agroecosystem model is given by The first is a deterministic term related to intrinsic plant/crop growth in response to water, carbon, and nitrogen concentrations under different levels of availability and use-efficiency (i.e., as complex field model-measurement versions of standardized empirical plant/crop growth curves obtained from laboratory measurement). This term is associated with the internal logic and interaction rules of the ecosystem model on plant/crop growth. The second and third terms are also deterministic and are associated with external climate trends forcing crop growth responses at the seasonal and inter-annual scale, respectively. Two additional terms considered are associated with first-order interactions between intrinsic crop growth and extrinsic climate variables related to specified inputs that are measurable or observable. The model also comprises hidden variables that are internally specified based on intervariable and interparameter relationships and internal logic (i.e., various ecological thresholds and constraints), so an additional crop response-climate interaction term considered to be Gaussian-distributed with known mean and variance was considered. Crop response from second-order and higher-variable interactions were aggregated in a final stochastic term considered to be driven primarily by environmental fluctuations and assumed to be normally distributed white noise, that is, N(0,1).

The statistical modeling framework used for this study applies the Sobol variance-based method for global sensitivity analysis to a complex ecosystem simulation model, calibrated for agricultural cropland. Given that most complex ecosystem models contain some degree of nonlinearity, it is not advised to apply local sensitivity analysis methods, as they can give misleading results. Instead, global sensitivity analysis should be applied. Here, the method used for the global sensitivity analysis was Saltelli’s [30–32] modification of Sobol’s [33] method. In Sobol’s method, the uncertainty associated with any important input parameter is propagated through the model resulting in a large contribution to the output variability [32, 33]. Main (first-order) effect and total effect indices are calculated for each combination of inputs, and a significant difference between the two indices signals an important interaction in the input(s) concerned. Decomposition of the output variance is given by: When the inputs are partitioned into different sets, the variance decomposition is The law of total variance states that: The main effect index is given by normalizing as follows: The main effects (first-order) sensitivity index indicates the relative importance of an individual input X_i in driving the output variance in a model. However, when one studies nonadditive models one must identify also noninfluential inputs (no significant effect on the output), and it is necessary to analyze higher-order terms. Being this the case, it would be necessary to estimate 2p-1 sensitivity indices which frequently involves a high computational cost—known as the “curse of dimensionality.” Saltelli et al. [31] showed that a necessary and sufficient condition for identifying noninfluential parameters consists in analyzing a more compact measure, the so-called total effect index. Let us assume that the set in (3) above contains only the input vector , and vector the remaining model inputs (i.e., ), then decomposition of the output variance can be rewritten as: where is the total contribution to the output variance due to all input parameters but X_i. Based on (4) (law of variance), can be written as as the expected amount of variance remaining unexplained if only input varies over its uncertainty range. In this way, if the difference is negligible, then input is noninfluential.

The total effect index for a given is then Accordingly, the total effect results from the addition of the main effect and the interactions between inputs. Estimating the set of all main effect and total effect indices provides a good description of the model sensitivities in an efficient way. Nominal values and ranges of variation for the four sites are presented in Table 2. Specific values for spring wheat seeding dates were established based on the research of Chan et al. [34] and for other parameters using the knowledge of the authors. Ranges of variation were defined either as a deviation of 75% from the nominal value, using information from White et al. [25], or as scaling values.

Currently, for Sobol’s method there is no readily available method for estimating confidence intervals for sensitivity indices, but standard error in the total and first/main indices can be estimated by resampling of model output via bootstrap numerical simulation. Sobol’s method was selected because previous bootstrap simulation results show that it is very robust in quantifying sensitivities and ranking parameters despite a large number of model executions [35]. For Sobol’s method, the total number of model evaluations , where is the base sample size and is the number of model parameters involved in the sensitivity analysis. Given that the accuracy of Sobol’s sensitivity method depends on the base sample size used, the sensitivity code was first run on the total number of version 4.1.1 Biome-BGC variables to identify all those variables that displayed negligible/no influence on NPP and SWC output. We then analyzed the model code/equations and their inter-connectivity to further verify that the variables displaying negligible/no response were irrelevant/noninfluential on the output NPP and SW (as per (1)). These variables were subsequently eliminated before performing the sensitivity runs based on the 10,000 model executions, thereby enabling us to focus our attention to a relevant or the influential subset and further reduced subset of model variables. As in our approach here, Confalonieri et al. (2010) [15], based on their sensitivity analysis of the rice model WARM in Europe, also first identify a subset of relevant/influential parameters, focusing their analysis on variables and performed 12,228 model executions. For an increasing number of ecological model executions (up to 49,152 runs), they report that 384 model executions were sufficient to identify variables with the highest influence on the output. Also, sensitivity and parameter ranking become constant with number of model executions greater than 6,144 runs (base sample of 439). A total of independent realizations/model evaluations of the ecosystem model were performed and response was evaluated at four long-term agricultural experimental sites situated in distinct soil and climate zones. For 10,000 model executions , and variables correspond to base sample sizes of and , respectively. By comparison, the number of base samples in our sensitivity analysis corresponds closely to the reported limit where sensitivity becomes constant with number of model executions (i.e., base sample 439, 6,144 executions for variables). Our study relies both on previous reported results focusing on these aspects for Sobol’s method and considers two different base sample sizes, measuring changes in sensitivity between the subset of to a reduced number of parameters/variables—whereby one approaches the reported limit where an increased number of model executions does not improve sensitivity estimates nor add increased accuracy/validity.

Implementation of Sobol’s method was performed using version 3.2.6 of SimLab software, which is freely distributed and can be downloaded from http://simlab.jrc.ec.europa.eu/. Custom software was written in C# and C++ to serve as a user interface and as an intermediary between the Biome-BGC model and SimLab software packages. The interface allowed the user to choose which input parameters would be analyzed, the type of distribution (and related parameters) to be used for sampling each input, and the number of times the model should be run. Orthogonal and uniform distribution of input parameters was assumed. This information was provided to the SimLab framework, and SimLab generated the run samples using Sobol’s method. Each run sample defined the value of each input parameter for a single execution of the model. The custom software then used the samples to automate the model runs. For each run, the parameter text files required by the model were written, the model was executed, and the selected output values were parsed. Once completed, the output values for all of the models runs were passed to SimLab, which performed the sensitivity analysis and generated, for each of the input parameters, the first order and total effect order sensitivity indices for each selected output.

Our sensitivity analysis considers realizations of the biogeochemical model, Biome-BGC at locations/sites. For each model simulation run, , the set of input parameters comprise a d-dimensional, real-valued array, denoted as . For model simulations, the full set of model inputs is then expressed as the array, . The corresponding model realization output from all simulation runs is the set where each simulation run output . With the model evaluated at each of s locations/sites, implies , where and . The set of model realizations for all sites can then be expressed as the matrix, For the agricultural application, two output variables (NPP and SWC) were selected as primary indicator variables for crop yield/production. Model input parameters included site-dependent soil/landscape, and site-independent crop growth/ecophysiological parameters.

The sites selected for this study were in Colorado, USA; and three sites in Canada: Harrow, Melfort, and Swift Current (Figure 1). The Colorado site (40°19′48′′ N, 104°42′ W) is located in the southern part of Weld County, Colorado, USA at an elevation of 1430 m. Its soil is classified as loam. The Harrow site (42°12′ N, 82°42′ W) is located in southern Ontario, Canada at an elevation of 185 m. Its soil is classified as clay loam. Melfort and Swift Current are located on the Canadian prairies in the province of Saskatchewan, Canada. The Melfort site (52°49′ N, 104°36′ W) is at an elevation of 482 m, and the Swift Current site (50°16′ N, 107°44′ W) is at an elevation of 825 m. The Melfort site is located in the Black soil zone (Udic Boroll), and its soil is classified as silty clay with slight evidence of wind erosion. The Swift Current site is in the Brown soil zone (Aridic Haploboroll) and its soil is classified as loam to silt loam [36, 37]. According to the Koeppen Climate classification system [38], the climate from Colorado is classified as cold semi-arid (Bsk), Harrow is classified as humid-continental (Dfa), and Swift Current and Melfort are classified as continental (Dfb).

The climate data used for the study included daily measurements of maximum and minimum air temperature, and precipitation. The climate data for the Colorado site ranged from 1980–2003 and were obtained as part of the DAYCENT model software package [39], from the Natural Resource Ecology Laboratory, Colorado State University, Fort Collins, Colorado, USA. Two climate data sets were used for the Harrow site: data ranging from 1980–2002 were provided by W. N. Smith, Eastern Cereal and Oilseed Research Centre, Agriculture and Agri-Food Canada, Ottawa, ON, Canada; and a second set ranging from 2003–2006 and was obtained through Environment Canada’s Meteorological Data Service, applicable to the climate station “Harrow Cda Auto Ontario” (Climate ID 6133362 and WMO ID 71298). For Melfort and Swift Current sites, daily values of climate variables were obtained from the AAFC-AAC Real Time Weather Network for the years 1980–2006. Other climate data required by Biome-BGC included daily values for average daytime temperature vapour pressure deficit, solar radiation, and day length. These values were generated independently for each site and time period using MT-CLIM software version 4.3 (http://www.ntsg.umt.edu/bioclimatology/mtclim/).

3. Results

Monthly average series of SWC and NPP were obtained for the three different (cold semiarid, humid-continental, and continental) climate distributions. For both SWC and NPP, the series from the humid-continental climate were larger than those from the cold semiarid and continental climates. The kernel density estimates for the monthly SWC and NPP averages are shown in Figures 2(c) and 2(d), respectively. As expected, the distribution of SWC from the continental and humid-continental climates had a thicker right tail than the distribution from the cold semiarid climate. For the humid-continental climate, the higher precipitation resulted in higher likelihood of larger SWC values, compared to the other two climates. A similar but less drastic trend was observed for the NPP distributions.

(a)

(b)

(c)

(d)

The sensitivity of NPP was explained by a large number of parameters. Sobol’s indices for main effects, obtained by considering the entire set of input parameters for each site and season, are given in Figure 3(a). Most of the input parameters had a main-effects index smaller than 0.2, except for the fraction of leaf nitrogen in rubisco (NR), which was the input that individually explained most of the variability in NPP for most cases. For the Melfort and Swift Current sites, the main indices for spring season were about 0.23, and for the Swift Current site, the index for the parameter day of the year to start new growth (NewG) during winter was 0.22. The main indices obtained from the sensitivity analysis performed on the reduced set of input parameters are presented in Figure 3(b). In this case, although there were more indices with values greater than 0.2, the number of important parameters was still small. The NewG parameter for winter was the most important input for all the sites, with values between 0.26 and 0.32. The parameter NR had main index values around 0.22 for Melfort and Swift Current during the spring. For the Harrow site, the input parameter current growth proportion (CG) had a main index of 0.2 for spring. The contribution of interactions to the NPP output variance is represented by the difference between the total effects and main-effects indices. The total effects (Sobol) index (TSI) is the addition of main effects and interaction index (TSI for all set of model inputs is computed as the summation of index values shown separately within insets (a) (main effects index) +(c) (interaction index), for a given parameter), whereby TSI > 0.8 (very important), 0.5 < TSI < 0.8 (important), 0.3 < TSI < 0.5 (unimportant), TSI < 0.3 irrelevant. These values are given in Figures 3(a) and 3(c) for the case of all influential input parameters, and Figures 3(b) and 3(d) for the reduced set of inputs. The parameter, NewG, is shown to play an important role explaining the variability of NPP across seasons, even during the winter season, with TSI ranging within 0.6-0.7, across experimental sites. But, for a reduced set of model inputs, TSI ranges within 0.25–0.35, becoming unimportant, mainly attributable to reduced interaction with other parameters. Canopy average specific leaf area (SLA) also had a significant interaction effect on NPP variability in both the full and partial sets of parameters, with interaction effect indices ranging from 0.38 to 0.58 for the summer season at all sites. These parameters besides NewG and SLA were: fraction of leaf nitrogen in rubisco (NR), canopy light extinction coefficient (), current growth proportion (CG), new fine root carbon to new leaf carbon ratio (FRC : LC), and maximum stomatal conductance (). Reducing the total number of parameters (reduced model) and without including interactions, sensitivity response for the leading parameter, , does not appear (see Figure 3).

(a)

(b)

(c)

(d)

Figure 3

Sensitivity indices for NPP based on (a) main effects using all inputs; (b) main effects using a reduced set of inputs; (c) interaction effects using all inputs; (d) interaction effects using a reduced set of inputs. The total effects (Sobol) index (TSI) is the addition of main effects and interaction index (TSI for all set of model inputs is computed as the summation of index values shown separately within insets (a) (main effects index) +(c) (interaction index), for a given parameter. Likewise, TSI for the reduced set of model inputs is the summation of indices shown within insets (b) (main effects index) +(d) (interaction index), for a given parameter, and TSI significance/importance is classified/ranked according to: TSI > 0.8 (very important), 0.5 < TSI < 0.8 (important), 0.3 < TSI < 0.5 (unimportant), and TSI < 0.3 irrelevant.

Main-effects and interaction indices across all model parameters are shown in Figures 4(a) and 4(b), respectively. The sensitivity of SWC was much less than for NPP and is mostly attributable to climate variability mediated by soil depth and texture. Sensitivity in model response to climate variability is evident via several additional parameters, namely: , CG, FRC : LC, , , NR, and SLA (refer to Table 1). These parameters were also involved in the sensitivity response of NPP. For all situations (total and seasonal) and in both sets of parameters for all the sites, the main effect indices for soil depth were between 0.60 (Colorado) and 0.92 (Harrow). As a result, the interaction effect indices for this parameter were small. During the summer season, the soil depth parameter had the smallest main-effect indices for all the sites, coinciding with the season in which the highest temperatures were present. For the SWC output distributions, Kolmogorov-Smirnov (KS) statistics resulted in -values of less than , except for the comparison between cold semiarid and continental climates, where the -value was . In all cases the results indicated that the compared distributions were significantly different. For the NPP distributions, the KS statistics resulted in -values of 0.0005 for the comparison of cold semi-arid and humid-continental climates, 0.0284 for the comparison of cold semiarid and continental climates, and 0.0086 for the comparison of humid-continental and continental climates. Therefore, the NPP distributions of the three different climates were significantly different.

(a)

(b)

(c)

(d)

4. Discussion

This study demonstrates that the input parameters in a complex terrestrial ecosystem model vary in importance with respect to their contribution to model outputs such as SWC and NPP. The parameter with the greatest influence on NPP in this study was NR, which finding is consistent with the studies of Tatarinov and Cienciala [26] and White et al. [25]. The parameter with the greatest contribution to variability in SWC in this study was soil depth. The development of any soil profile occurs over a long period of time, and generally reflects the accumulation of the effects of variations in short-term weather and long-term climate. Soils with greater depth often reflect a longer history of development under conditions of higher precipitation. Therefore, a strong connection between soil depth and SWC is not surprising. The finding that some model input parameters are more important than others has significant implications for many subdisciplines in the broader field of ecosystems modeling. For example, the development of agricultural production monitoring systems over large areas is currently being undertaken worldwide. These systems are intended to generate crop yield estimates at multiple scales using ecosystems models, but they face strict limitations in available computing resources. In order to develop an efficient system that can deliver estimates in a timely fashion, it is highly desirable to limit the number of variables within the system. Sensitivity analyses such as the current study and others will help to guide the developers of agricultural monitoring systems in their choice of parameters to include and exclude in such efforts. Our findings reveal large seasonal differences in the sensitivity of Biome-BGC ecosystem model parameterized for spring wheat. Significant changes in the variance and skew of climate variable distributions accompanies changing climatic regime, in addition to changes in the mean, and can lead to overestimation in complex models, if not corrected for. There are significant opportunities to increase the precision of ecosystem models to predict wheat yield and its variability due to climate (radiation, air temperature, precipitation) by staging these models in time; such that parameters within a given phenological stage are most sensitive to leading climate covariates within each stage. Also, the best opportunity to reduce complexity in complex models is under drier moisture-limited environments, where there is significant soil variability (semiarid/arid), where leaf-area dynamics and stomatal conductance are key mediators in development, growth and survival. Our findings lend further support (in the case of agricultural cropland ecosystems) for the need to apply a hierarchical, componentwise approach in reducing complex ecosystem models and predicting crop response, as recently applied to forest ecosystems by Wang et al. [40]. For a typical cropland ecosystem of spring wheat, we find the same “domain” order holds for sets of input parameters ranked according to their relative sensitivity on model output, as previously reported by Wang et al. [40] for Biome-BGC applied to a modeled forest ecosystem, namely, phenology, leaf-area dynamics, light-interception, morphological, ecophysiological. This indicates promise for further testing and adapting Biome-BGC to agroecosystems, characterized by rapid rates of carbon, nitrogen and water cycling due to tillage, fertilization, and irrigation disturbances, within realistic input parameter ranges. We find that NPP becomes far more sensitive to input parameters when parameter interactions are considered. Specifically, reducing the total number of parameters (reduced model) and without including interactions, sensitivity response for the leading parameter, , (maximum stomatal conductance) does not appear. Given the importance of for controlling plant/crop water-use efficiency, interactions should be considered, so that output is significantly (and more realistically) sensitive to its variation.

The significant seasonal changes in sensitivity across phenology, leaf area dynamics and light-interception parameters, indicates that one should not attempt to simply seasonally adjust only phenology parameters alone—as this is not sufficient. Nonetheless, our findings also show the major control that phenology parameters have on NPP and SWC predictions. By ranking input parameters into domains of similar relative sensitivity, hierarchical reduction of complex ecosystem models can be performed more objectively. This approach offers potential for improving predictions of NPP and SWC via phenology and other sets of parameters linked with annual climate variability (i.e., seasonal). This would enable a better determination of storage and cycling of water and nutrients at the start and end of the growing season. Improved prediction linked with sensitivity ranking of parameters would also support more objective re-structuring of complex models and gaining more reliable, clearer insights associated with input parameters and interactions within these models on simulation output. This would enhance agricultural decision-making and improved system resilience when informed using complex ecosystem models.

Acknowledgments

This research was funded by the Canadian Federal Government Program for Sustainable Agriculture Environmental Systems (SAGES) Program. Climate data for the Harrow site was provided by W. N. Smith and B. B. Grant (Agriculture and Agri-Food Canada). The authors also acknowledge Mr. J. Lays (University of Lethbridge, Co-op Internship Program) for computer programming assistance and Dr. Tracy A. Porcelli for helpful comments on earlier drafts of this manuscript. The authors thank P. Thornton of the National Center for Atmospheric Research (NCAR), Boulder, CO and the Numerical Terradynamic Simulation Group (NTSG) at the University of Montana, Missoula, MT, USA, for providing Biome-BGC version 4.1.2 computer software. NCAR is sponsored by the National Science Foundation.

References

C. D. Canham, J. J. Cole, and W. K. Lauenroth, “Models in Ecosystem Science,” in Models in Ecosystem Science, C. D. Canham, J. J. Cole, and W. K. Lauenroth, Eds., pp. 1–10, Princeton University Press, Princeton, NJ, USA, 2003.
View at: Google Scholar
L. R. Ahuja and L. Ma, “Parameterization of agricultural system models: Current approaches and future needs,” in Agricultural System Models in Field Research and Technology Transfer, L. R. Ahuja, L. Ma, and T. A. Howell, Eds., pp. 273–316, Lewis Publishers, Boca Raton, Fla, USA, 2002.
View at: Google Scholar
G. M. Cox, J. M. Gibbons, A. T. A. Wood, J. Craigon, S. J. Ramsden, and N. M. J. Crout, “Towards the systematic simplification of mechanistic models,” Ecological Modelling, vol. 198, no. 1-2, pp. 240–246, 2006.
View at: Publisher Site | Google Scholar
M. Fuentes, T. G. F. Kittel, and D. Nychka, “Sensitivity of ecological models to their climate drivers: statistical ensembles for forcing,” Ecological Applications, vol. 16, no. 1, pp. 99–116, 2006.
View at: Publisher Site | Google Scholar
J. Cariboni, D. Gatelli, R. Liska, and A. Saltelli, “The role of sensitivity analysis in ecological modelling,” Ecological Modelling, vol. 203, no. 1-2, pp. 167–182, 2007.
View at: Publisher Site | Google Scholar
J. W. Hansen and J. W. Jones, “Scaling-up crop models for climate variability applications,” Agricultural Systems, vol. 65, no. 1, pp. 43–72, 2000.
View at: Publisher Site | Google Scholar
M. P. R. Haaker and P. J. T. Verheijen, “Local and global sensitivity analysis for a reactor design with parameter uncertainty,” Chemical Engineering Research and Design, vol. 82, no. 5, pp. 591–598, 2004.
View at: Publisher Site | Google Scholar
S. Wagner, “Global Sensitivity Analysis of Predictor Models in Software Engineering,” IEEE Computer Society Note 12690513, 2007.
View at: Google Scholar
A. Saltelli, S. Tarantola, and K. P. S. Chan, “A quantitative model-independent method for global sensitivity analysis of model output,” Technometrics, vol. 41, no. 1, pp. 39–56, 1999.
View at: Google Scholar
D. M. King and B. J. C. Perera, “Sensitivity analysis for evaluating importance of variables used in an urban water supply planning model,” in Proceedings of the International Congress on Modelling and Simulation (MODSIM '07), L. Oxley and D. Kulasiri, Eds., pp. 2768–2774, Canberra: Modelling and Simulation Society of Australia and New Zealand, 2007.
View at: Google Scholar
E. P. Poeter, M. C. Hill, E. R. Banta, E. R. Mehl, and S. Chistensen, UCODE_2005 and Six Other Computer Codes for Universal Sensitivity Analysis, Calibration, and Uncertainty Evaluation: Constructed, using the JUPITER API (Joint Universal Parameter Identification and Evaluation of Reliability), U.S. Department of the Interior and U.S. Geological Survey (USGS), 2005.
A. Saltelli, M. Ratto, S. Tarantola, and F. Campolongo, “Sensitivity analysis for chemical models,” Chemical Reviews, vol. 105, no. 7, pp. 2811–2827, 2005.
View at: Publisher Site | Google Scholar
P. Young, S. Parkinson, and M. Lees, “Simplicity out of complexity in environmental modelling: Occam's razor revisited,” Journal of Applied Statistics, vol. 23, no. 2-3, pp. 165–210, 1996.
View at: Google Scholar
R. Confalonieri, “Monte Carlo based sensitivity analysis of two crop simulators and considerations on model balance,” European Journal of Agronomy, vol. 33, no. 2, pp. 89–93, 2010.
View at: Publisher Site | Google Scholar
R. Confalonieri, G. Bellocchi, S. Tarantola, M. Acutis, M. Donatelli, and G. Genovese, “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,” Environmental Modelling and Software, vol. 25, no. 4, pp. 479–488, 2010.
View at: Publisher Site | Google Scholar
M. Lamboni, D. Makowski, S. Lehuger, B. Gabrielle, and H. Monod, “Multivariate global sensitivity analysis for dynamic crop models,” Field Crops Research, vol. 113, no. 3, pp. 312–320, 2009.
View at: Publisher Site | Google Scholar
G. M. Richter, M. Acutis, P. Trevisiol, K. Latiri, and R. Confalonieri, “Sensitivity analysis for a complex crop model applied to Durum wheat in the Mediterranean,” European Journal of Agronomy, vol. 32, no. 2, pp. 127–136, 2010.
View at: Publisher Site | Google Scholar
D. Makowski, C. Naud, M. H. Jeuffroy, A. Barbottin, and H. Monod, “Global sensitivity analysis for calculating the contribution of genetic parameters to the variance of crop model prediction,” Reliability Engineering and System Safety, vol. 91, no. 10-11, pp. 1142–1147, 2006.
View at: Publisher Site | Google Scholar
H. Varella, M. Guérif, and S. Buis, “Global sensitivity analysis measures the quality of parameter estimation: the case of soil parameters and a crop model,” Environmental Modelling and Software, vol. 25, no. 3, pp. 310–319, 2010.
View at: Publisher Site | Google Scholar
S. W. Running and E. R. Hunt Jr., “Generalization of a forest ecosystem process model for other biomes, BIOME-BGC and an application for global-scale models,” in Scaling Physiological Processes: Leaf to Globe, J. R. Ehleringer and C. B. Field, Eds., pp. 141–158, Academic Press, New York, NY, USA, 1993.
View at: Google Scholar
P. E. Thornton, Description of a numerical simulation model for predicting the dynamics of energy, water, carbon and nitrogen in a terrestrial ecosystem, Ph.D. thesis, University of Montana, Missoula, Mont, USA, 1998.
P. E. Thornton, B. E. Law, H. L. Gholz et al., “Modeling and measuring the effects of disturbance history and climate on carbon and water budgets in evergreen needleleaf forests,” Agricultural and Forest Meteorology, vol. 113, no. 1–4, pp. 185–222, 2002.
View at: Publisher Site | Google Scholar
S. W. Running and J. C. Coughlan, “A general model of forest ecosystem processes for regional applications I. Hydrologic balance, canopy gas exchange and primary production processes,” Ecological Modelling, vol. 42, no. 2, pp. 125–154, 1988.
View at: Google Scholar
S. W. Running and S. T. Gower, “FOREST-BGC, A general model of forest ecosystem processes for regional applications. II. Dynamic carbon allocation and nitrogen budgets,” Tree Physiology, vol. 9, pp. 147–160, 1991.
View at: Google Scholar
M. A. White, P. E. Thornton, S. W. Running, and R. R. Nemani, “Parameterization and sensitivity analysis of the BIOME-BGC terrestrial ecosystem model: net primary production controls,” Earth Interactions, vol. 4, pp. 1–85, 2000.
View at: Google Scholar
F. A. Tatarinov and E. Cienciala, “Application of BIOME-BGC model to managed forests. 1. Sensitivity analysis,” Forest Ecology and Management, vol. 237, no. 1–3, pp. 267–279, 2006.
View at: Publisher Site | Google Scholar
H. X. Su and W. G. Sang, “Simulations and analysis of net primary productivity in Quercus liaotungensis forest of Donglingshan mountain range in response to different climate change scenarios,” Journal of Integrative Plant Biology, vol. 46, pp. 1281–1291, 2004.
View at: Google Scholar
Q. Wang, M. Watanabe, and Z. Ouyang, “Simulation of water and carbon fluxes using BIOME-BGC model over crops in China,” Agricultural and Forest Meteorology, vol. 131, no. 3-4, pp. 209–224, 2005.
View at: Publisher Site | Google Scholar
A. Marrel, B. Iooss, M. Jullien, B. Laurent, and E. Volkova, “Global sensitivity analysis for models with spatially dependent outputs,” Environmetrics, vol. 22, no. 3, pp. 383–397, 2011.
View at: Google Scholar
A. Saltelli, “Making best use of model evaluations to compute sensitivity indices,” Computer Physics Communications, vol. 145, no. 2, pp. 280–297, 2002.
View at: Publisher Site | Google Scholar
A. Saltelli, “Sensitivity analysis for importance assessment,” Risk Analysis, vol. 22, no. 3, pp. 579–590, 2002.
View at: Publisher Site | Google Scholar
A. Saltelli, S. Tarantola, F. Campolongo, and M. Ratto, Sensitivity Analysis in Practice: A Guide to Assessing Scientific Models, John Wiley & Sons, New York, NY, USA, 2004.
I. M. Sobol, “Sensitivity estimates for nonlinear mathematical models,” Mathematical Modelling and Computational Experiments, vol. 1, pp. 407–414, 1993.
View at: Google Scholar
K. Chan, A. Saltelli, and S. Tarantola, “Sensitivity analysis of model output: variance-based methods make the difference,” in Proceedings of the 29th Conference on Winter Simulation, IEEE Computer Society, 1997.
View at: Google Scholar
J. Yang, “Convergence and uncertainty analyses in Monte-Carlo based sensitivity analysis,” Environmental Modelling and Software, vol. 26, no. 4, pp. 444–457, 2011.
View at: Publisher Site | Google Scholar
C. A. Campbell and R. P. Zentner, “Crop Production and Soil Organic Matter in Long-Term Crop Rotations in the Semi-Arid Northern Great Plains of Canada,” in Soil Organic Matter in Temperate Agroecosystems: Long-Term Experiments in North America, E. A. Paul, K. Paustian, E. T. Elliot, and C. V. Cole, Eds., pp. 317–334, CRC Press, Boca Raton, Fla, USA, 1997.
View at: Google Scholar
C. A. Campbell, G. P. Lafond, A. P. Moulin, L. Townley-Smith, and R. P. Zentner, “Crop production and soil organic matter in long-term crop rotations in the sub-humid Northern Great Plains of Canada,” in Soil Organic Matter in Temperate Agroecosystems: Long-Term Experiments in North America, E. A. Paul, K. Paustian, E. T. Elliot, and C. V. Cole, Eds., pp. 297–315, CRC Press, Boca Raton, Fla, USA, 1997.
View at: Google Scholar
M. Kottek, J. Grieser, C. Beck, B. Rudolf, and F. Rubel, “World map of the Köppen-Geiger climate classification updated,” Meteorologische Zeitschrift, vol. 15, no. 3, pp. 259–263, 2006.
View at: Publisher Site | Google Scholar
E. Stehfest, M. Heistermann, J. A. Priess, D. S. Ojima, and J. Alcamo, “Simulation of global crop production with the ecosystem model DayCent,” Ecological Modelling, vol. 209, no. 2–4, pp. 203–219, 2007.
View at: Publisher Site | Google Scholar
W. Wang, K. Ichii, H. Hashimoto et al., “A hierarchical analysis of terrestrial ecosystem model Biome-BGC: equilibrium analysis and model calibration,” Ecological Modelling, vol. 220, no. 17, pp. 2009–2023, 2009.
View at: Publisher Site | Google Scholar

Copyright

Copyright © 2012 Nathaniel K. Newlands et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

PDF Download Citation

Download other formats

Order printed copies

Views

2051

Downloads

1629

Citations