Articles | Volume 26, issue 13
Research article
07 Jul 2022
Research article |  | 07 Jul 2022

Revisiting parameter sensitivities in the variable infiltration capacity model across a hydroclimatic gradient

Ulises M. Sepúlveda, Pablo A. Mendoza, Naoki Mizukami, and Andrew J. Newman

Despite the Variable Infiltration Capacity (VIC) model being used for decades in the hydrology community, there are still model parameters whose sensitivities remain unknown. Additionally, understanding the factors that control spatial variations in parameter sensitivities is crucial given the increasing interest in obtaining spatially coherent parameter fields over large domains. In this study, we investigate the sensitivities of 43 soil, vegetation and snow parameters in the VIC model for 101 catchments spanning the diverse hydroclimates of continental Chile. We implement a hybrid local–global sensitivity analysis approach, using eight model evaluation metrics to quantify sensitivities, with four of them formulated from runoff time series, two characterizing snow processes, and the remaining two based on evaporation processes. Our results confirm an overparameterization for the processes analyzed here, with only 12 (i.e., 28 %) parameters found to be sensitive, distributed among soil (7), vegetation (2) and snow (3) model components. Correlation analyses show that climate variables – in particular, mean annual precipitation and the aridity index – are the main controls on parameter sensitivities. Additionally, our results highlight the influence of the leaf area index on simulated hydrologic processes – regardless of the dominant climate types – and the relevance of hard-coded snow parameters. Based on correlation results and the interpretation of spatial sensitivity patterns, we provide guidance on the most relevant parameters for model calibration according to the target processes and the prevailing climate type. Overall, the results presented here contribute to an improved understanding of model behavior across watersheds with diverse physical characteristics that encompass a wide hydroclimatic gradient from hyperarid to humid systems.

1 Introduction

Over the past 4 decades, the increasing demand for more realistic spatial representations of water storages and fluxes across large domains has motivated the development of more complex physics-based models (e.g., Niu et al., 2011; Clark et al., 2015; Lawrence et al., 2019). The progress in this field has been partly facilitated by new observational datasets (e.g., McCabe et al., 2017; Berg et al., 2018) and advances in computing (see discussion on tradeoffs in Clark et al., 2017), enabling hydrological characterizations of national (e.g., Tian et al., 2017; Zink et al., 2017), continental (e.g., Xia et al., 2012; Abbaspour et al., 2015) and global (e.g., Schmied et al., 2014; Arheimer et al., 2020) domains.

Although spatially constant parameters can be used for large-domain applications, improving model realism requires the specification of parameter values that reflect spatial heterogeneities in landscape properties. Because increasing model complexity is often associated with a larger number of parameters, many models rely on lookup tables to assign soil thermal and hydraulic parameters as well as vegetation optical and physiological parameters for each modeling unit (e.g., Mitchell et al., 2004; Yang et al., 2011). However, parameter uncertainties may be considerable (Rosero et al., 2010; Hou et al., 2012), and such a problem is exacerbated by the existence of parameters that, despite being “adjustable” (e.g., runoff generation parameters), remain fixed and hard-coded (Mendoza et al., 2015a; Cuntz et al., 2016).

To address overparameterization problems that are typical in environmental models, sensitivity analysis has become a key tool that provides information on which parameter values are the most influential for the dynamics of specific model responses (Razavi and Gupta, 2015). The outcomes of sensitivity analysis not only help to improve the understanding of model functioning but also to inform decisions regarding parameter estimation problems. The literature provides many examples of sensitivity analysis studies with process-based hydrological models, including the Biosphere–Atmosphere Transfer Scheme, BATS (Bastidas et al., 1999); TOPKAPI (Foglia et al., 2009); PRMS (Mendoza et al., 2015b); the Noah land surface model (Bastidas et al., 2006; Rosero et al., 2010); Noah-MP (Mendoza et al., 2015a; Cuntz et al., 2016); the Simple Biosphere model, SiB3 (Prihodko et al., 2008; Rosolem et al., 2012); the MESH modeling system (Razavi and Gupta, 2016); the Community Land Model, CLM (Göhler et al., 2013; Massoud et al., 2019); and the Variable Infiltration Capacity (VIC) model (e.g., Demaria et al., 2007; Melsen et al., 2016), which is one of most popular modeling platforms in the hydrology community (Addor and Melsen, 2019).

The VIC model (Liang et al., 1994; Hamman et al., 2018) has been used for myriad applications all over the world, including snow modeling (Andreadis et al., 2009; Chen et al., 2014), streamflow forecasting (Wood et al., 2005; DeChant and Moradkhani, 2014), water balance studies (Mizukami et al., 2016; Vásquez et al., 2021), extreme event characterization (Melsen et al., 2019), land use change impacts (Chawla and Mujumdar, 2015) and climate change impact assessments (e.g., Vano and Lettenmaier, 2014; Chegwidden et al., 2019). Despite the large number of parameters contained in VIC – either “free” (e.g., the infiltration shape parameter “INFILT”, the exponent in baseflow curve) or “observable” (e.g., leaf area index) – many studies have relied on the calibration of only two or three soil-related parameters (Huang and Liang, 2006; Chawla and Mujumdar, 2015). Conversely, other authors have advocated for characterizing parameter sensitivities by using different approaches and sensitivity metrics as well as including different parameters. For example, Liang and Guo (2003) assessed the sensitivity of annual runoff, annual evapotranspiration (ET), annual mean soil moisture and annual mean sensible heat flux to variations in five soil and vegetation parameters at three experimental locations (i.e., point scale), finding that sensitivities varied with climatic and physiographic site characteristics. Demaria et al. (2007) examined sensitivities in simulated catchment-scale runoff responses using lumped VIC configurations, a Monte Carlo method and five objective functions computed for four basins with varying hydroclimates. They concluded that (i) 3 (out of 10) soil parameters dominated the simulated runoff response, (ii) the INFILT parameter and the drainage parameter (Expti) depended strongly on local hydroclimatology, and (iii) the baseflow model formulation is overparameterized.

Subsequent studies aiming at calibrating the VIC model to simulate observed catchment-scale responses have revisited its parametric sensitivity. Mendoza et al. (2015b) applied the Distributed Evaluation of Local Sensitivity Analysis method (DELSA; Rakovec et al., 2014) to find the parameters that provided the largest sensitivities in root-mean-square errors (RMSE) between simulated and observed streamflow; they showed that 9 (out of 34) parameters provided the largest sensitivities in three subcatchments from the Upper Colorado River basin. Melsen et al. (2016) also used the DELSA method to find influential parameters in three catchments located in Switzerland, identifying 4 (out of 28) very sensitive parameters for three calibration metrics. Wi et al. (2017) applied the method of Morris (1991) to quantify parameter sensitivities on the Nash–Sutcliffe efficiency (NSE; Nash and Sutcliffe, 1970) with daily flows, finding 6 soil parameters and 2 temperature threshold parameters as the most influential (out of 15). Gou et al. (2020) characterized the sensitivities provided by 13 soil parameters across 14 catchments in China, finding that INFILT, Depth1 and Depth2 dominated streamflow responses. Lilhare et al. (2020) applied the Variogram Analysis of Response Surfaces (VARS) method (Razavi and Gupta, 2016) to examine the sensitivities of three streamflow performance metrics to variations in six soil parameters across 10 catchments in Canada, finding that the INFILT and Depth2 parameters dominated streamflow responses. Yeste et al. (2020) quantified relative sensitivities provided by five soil parameters to water balance components across 31 basins on the Iberian Peninsula, concluding that INFILT and Depth2 control runoff components and evapotranspiration (ET). Finally, Melsen and Guse (2021) characterized VIC parameter sensitivities for a historical (1985–2008) and future (2070–2093) period in 605 catchments of the conterminous United States, finding that, in the historical period, Rmin, Depth2 and Expt2 controlled average streamflow, while Ds, DsMAX and many more parameters influenced streamflow timing. Melsen and Guse (2021) also projected increased (decreased) sensitivities to Depth2 (snow parameters) for the future period when examining average streamflow and increased (decreased) future sensitivities to deep soil (snow) parameters when looking at discharge timing.

Table 1A summary of the sensitivity analysis (SA) studies conducted with VIC that incorporate at least five parametersa.

a The studies are listed in order of publication date. The present study has been added for completeness. b We exclude two routing parameters that were found to be sensitive but were not used in the other studies. The abbreviations/acronyms used in the table are as follows: CDF – cumulative distribution function; KGE – Kling–Gupta efficiency (Gupta et al., 2009); PBIAS – percent bias in flow volumes; LAI – leaf area index; SWE – snow water equivalent; FFA – fractional factorial analysis (Montgomery, 1991); DELSA – Distributed Evaluation of Local Sensitivity Analysis (Rakovec et al., 2014); DT – Delta test (Pi and Peterson, 1994); SOT – sum-of-trees model (Chipman et al., 2010); MARS – multivariate adaptive regression splines (Friedman, 1991); VARS – Variogram Analysis of Response Surfaces (Razavi and Gupta, 2016); SRC – standardized regression coefficients (Saltelli et al., 2008); VISCOUS – Variance-based sensitivity analysis using copulas (Sheikholeslami et al., 2021).

Download XLSX

Table 1 summarizes the main characteristics of parameter sensitivity studies with VIC. Note that we have excluded a recent study conducted by Sheikholeslami et al. (2021), who quantified parameter sensitivities in a modified version of VIC – specifically, using a slow linear reservoir (Gharari et al., 2019) model instead of the traditional ARNO formulation. Most studies listed in Table 1 focused on streamflow responses, attributing the largest sensitivities to a few soil parameters (Demaria et al., 2007; Gou et al., 2020; Lilhare et al., 2020). Only two studies – also characterizing streamflow responses – have included a large number of soil-, vegetation- and snow-related parameters (Melsen et al., 2016; Mendoza et al., 2015b). Additionally, only two studies (Chaney et al., 2015; Bennett et al., 2018) aimed to characterize sensitivities across model grid cells. Chaney et al. (2015) quantified the effects of nine parameters on annual flow biases, runoff seasonality and daily flow extremes at 1 resolution grid cells across the globe. Bennett et al. (2018) examined the sensitivities of projected changes in water balance components to variations in 46 VIC parameters, across a suite of ∼7 km grid cells in the Colorado River basin. Sensitivity analyses at higher resolutions are particularly relevant for the hydrology community, considering recent developments in meteorological datasets (Tang et al., 2021) and global, gridded runoff datasets (e.g., Do et al., 2018; Ghiggi et al., 2019) as well as the increasing interest in improving the calibration density – i.e., using high-resolution data in the calibration process – in distributed hydrology and land surface models (Yang et al., 2019).

In this paper, we quantify VIC parameter sensitivities across 5574 grid cells (0.05×0.05) covering 101 catchments located in continental Chile, including a suite of 43 standard and hard-coded parameters as well as a set of metrics that span different runoff, ET and snow processes. The results presented here contribute to an improved understanding of model behavior across watersheds with diverse physical characteristics, spanning from hyperarid to humid hydroclimates. With this, we seek to answer the following research questions:

  1. Are there other vegetation and snow parameters, either standard or hard-coded, affecting simulated runoff responses in VIC?

  2. What are the effects of standard and hard-coded parameters on other simulated processes?

  3. How do parameter sensitivities relate to local climatic and physiographic characteristics?

Figure 1Spatial distribution of climatic and physiographic attributes across all grid cells: (a) mean annual precipitation (period 1979–2020), (b) mean annual temperature (period 1979–2020), (c) mean elevation, (d) aridity index and (e) bare-soil fraction. In each panel, the black thick lines represent the boundaries of six basins representative of the hydroclimatic diversity within the study domain; from north to south, the basins are as follows: (a) Loa River upstream of Lequena reservoir, (b) Pulido River at Vertedero, (c) Colorado River before the junction with the Maipo River, (d) Palos River at the junction with the Colorado River, (e) Biobío River at Rucalhue and (f) Cautín River at Cajón.

2 Study domain and data

In this work, we select 101 catchments with near-natural hydrological regimes from the CAMELS-CL dataset (Alvarez-Garreton et al., 2018). The selected basins span a total area of 139 350 km2 (i.e., 19 % of the territory of continental Chile) and meet the following criteria: (i) a maximum threshold value of 5 % for the relationship between the annual volume of water assigned as permanent consumptive rights and the average annual flow (Table 3 in Alvarez-Garreton et al., 2018) and (ii) the absence of large reservoirs within each catchment. The location, hydroclimatic and land cover characteristics across the domain are shown in Fig. 1, and the descriptors used to characterize the grid cells are listed in Table 2. These catchments cover a wide range of physiographic attributes, with drainage areas spanning from 100 to 7500 km2, mean elevations ranging from 119 to 4824 m a.s.l. (above sea level), mean slopes varying from 52 to 306 m km−1 and markedly different land cover types (ranging from completely covered by native forest or grassland to fully covered by impermeable land). Moreover, the selected basins represent the diversity in hydroclimatic conditions across the country. Figure 2 illustrates this by showing the aridity indices and seasonal cycles of rainfall, snowfall, runoff and temperature for a sample of six basins with very different hydroclimatic regimes. The hydrology of catchment 2101001 (Loa River before Lequena dam, Fig. 2a) is influenced by arid conditions between March and November, and altiplanic winter events triggering runoff increases between December and March; towards the south, there is a transition from arid to semiarid conditions (see progression in Fig. 2b and c), with precipitation events occurring mostly during fall and winter (especially May–August), favoring the accumulation of snow in the headwaters of Andean catchments, and thus snowmelt-driven regimes. Catchments 7115001 (Palos River at the junction with the Colorado River; Fig. 2d) and 8317001 (Biobío River at Rucalhue; Fig. 2e) reflect the transition towards mixed regimes, with larger contributions of winter rainfall events to runoff. Finally, catchment 9129002 (Cautín River at Cajón; Fig. 2f) has a rainfall-dominated hydrological regime, with the largest runoff volumes during the winter season (i.e., June–September).

Figure 2Seasonal cycles of catchment-averaged rainfall, snowfall, runoff and temperature (period 1979–2020) for six basins representative of the hydroclimatic diversity within the study domain. From north to south, the basins are as follows: (a) Loa River upstream of Lequena reservoir, (b) Pulido River at Vertedero, (c) Colorado River before the junction with the Maipo River, (d) Palos River at the junction with the Colorado River, (e) Biobío River at Rucalhue and (f) Cautín River at Cajón. The location of these catchments is shown in Fig. 1.


Table 2List of physiographic and hydroclimatic attributes used to characterize model grid cells.

Download Print Version | Download XLSX

In this study, meteorological forcing data are obtained from various sources. Time series of daily precipitation and maximum, average and minimum daily temperature are obtained from the CR2MET meteorological dataset (Boisier et al., 2018), introduced in DGA (2017), which provides data for continental Chile at a horizontal resolution of 0.05×0.05 (∼5 km) for the period from 1979 to 2016. The precipitation product builds upon a statistical post-processing technique that uses topographic descriptors and simulated meteorological variables from ERA-Interim (Dee et al., 2011) and ERA5 (C3S and Copernicus Climate Change Service, 2017) as predictors and uses daily precipitation records as the predictand. A similar approach is used to generate time series of daily maximum and minimum temperatures, including additional predictors from MODIS land surface products. Daily precipitation and temperature variables are disaggregated into 3-hourly time steps using the sub-daily distribution provided by ERA-Interim. Finally, relative humidity and wind speed were obtained from a blend between ERA-Interim and ERA5, which was subsequently rescaled to the CR2MET horizontal grid through spatial interpolation. It is important to note that this product combination was created because ERA5 was not available during the entire study period (1985–2015) at the time of data acquisition (early 2018, when only 2010–2016 data were available). However, the updated reanalysis information, despite the short time coverage, was included due to various developmental improvements.

3 Methods

3.1 Hydrological model

The Variable Infiltration Capacity (VIC; Liang et al., 1994) model is a semi-distributed, physically based hydrological model that simulates snow accumulation and melt, evapotranspiration (ET), canopy interception, surface runoff, baseflow, and other hydrological processes at daily or sub-daily time steps. While the model was originally designed as a land surface scheme for coupled simulations within Earth system models (Liang et al., 1994), most applications have involved uncoupled simulations for hydrological characterizations; accordingly, the literature reports many attempts to improve process representations (e.g., Liang et al., 1996, 1999; Cherkauer et al., 2003; Andreadis et al., 2009). VIC is predominantly used in the USA (Addor and Melsen, 2019), with many studies focused on water and energy balances (e.g., Andreadis and Lettenmaier, 2006; Cayan et al., 2010); however, its use has expanded to other geographical domains, including China (e.g., Zhao et al., 2013; Gou et al., 2021), Chile (e.g., Vásquez et al., 2021; Vicuña et al., 2021), Europe (e.g., Lohmann et al., 1998; Roudier et al., 2016) and globally (e.g., Shukla et al., 2013; Yang et al., 2021). Ongoing community efforts using the VIC model include the NASA Land Information System (LIS;, last access: 25 January 2022), the NASA Land Data Assimilation System (LDAS;, last access: 25 January 2022) and the Regional Arctic System Model (RASM;, last access: 25 January 2022). This study uses VIC version 4.1.2.g, which can be downloaded from (last access: 13 June 2022), along with other versions.

In VIC, the domain of interest can be spatially discretized into grid cells. Sub-grid land use type variability is accounted for by providing vegetation tiles and the fractional areas, for which water and energy balance equations are solved separately; model states and fluxes are then spatially averaged to provide results at the pixel scale. In VIC, each grid cell can have up to three soil layers: the two top layers represent the interaction between moisture and vegetation, and the bottom soil layer is used to simulate baseflow processes. VIC does not incorporate an aquifer at the bottom of the soil column nor lateral exchange of fluxes between grid cells. Finally, snowpack dynamics are simulated by a two-layer mass and energy balance model (Cherkauer et al., 2003; Andreadis et al., 2009), where the surface layer solves the energy exchange between the snowpack and the atmosphere, and the lower layer stores the excess snow mass from the upper surface layer.

3.2 Parameters considered for sensitivity analysis

We considered a suite of 43 parameters (Table 3) to incorporate most of the soil, vegetation and snow processes simulated by VIC. It should be noted that three of the snow parameters are not exposed to model users (NEW ALB, ALB AA and ALB THA), although the associated relationships and default values were proposed decades ago (USACE, 1956). For parameters with monthly variations, we examined sensitivities using regularization “superparameters” (Tonkin and Doherty, 2005), also called multipliers (Pokhrel and Gupta, 2010), which are uniformly applied over all monthly values. Hence, multipliers are used for the leaf area index (LAI), vegetation albedo (ALB), vegetation roughness length (ROU) and vegetation displacement (DIS). Despite the fact that some of these parameters are considered observable, a non-negligible degree of uncertainty may be involved in their determination; an example of this is the LAI parameter (Tian et al., 2002; Fang et al., 2012, 2013), whose implementation in many hydrology and land surface models is simplified through static monthly values. Assumptions like this motivate us to explore the sensitivity of these types of parameters, which may have the potential to be included in the calibration process.

Table 3Parameters of the VIC model considered in this study.

* This parameter is temporally distributed (monthly variations); therefore, its sensitivity is analyzed based on multipliers. Although the description and units refer to actual parameters of VIC, parameter values in bold represent the multiplier values (instead of actual parameters).

Download XLSX

Despite the aim to include the largest possible number of parameters, some of them were discarded for different reasons. For example, a few soil parameters (e.g., soil bubbling pressure) are not active unless the frozen soil algorithm is turned on. We also excluded the “trunk ratio” parameter – i.e., the ratio of total tree height that is trunk (no branches) – because it is activated only in those grid cells with forest (i.e., vegetation class with overstory, spanning 22 % of our study domain) as the land cover type. Finally, we found five mutually related soil parameters that do not allow independent variations: soil bulk density (bd), soil particle density (sd), fractional soil moisture content at the critical point (θcr), fractional soil moisture content at the wilting point (θwp) and residual soil moisture (θr). These parameters are related as follows:

(1) θ cr θ wp θ r 1 - bd sd .

From these five parameters, we only include θr and bd, as (1) perturbing θr and bd values did not affect numerical solutions, and (2) Bennett et al. (2018) showed that θcr, θwp, bd and sd did not have substantial effects on model simulations. Finally, those parameters that showed little or no sensitivity in the initial phases of the study were purposely discarded.

3.3 Sensitivity analysis approach

We used the Distributed Evaluation of Local Sensitivity Analysis (DELSA; Rakovec et al., 2014) method, which is a derivative-based, hybrid local–global approach. DELSA combines elements from the method of Morris (1991), the Sobol' method (Sobol', 2001) and regional sensitivity analysis (Hornberger and Spear, 1981), and it provides robust results with fewer model simulations compared with variance-based global methods such as the Sobol' technique. Although DELSA only examines first-order sensitivities, it has unexplored potential to be expanded with the aim of quantifying parameter interactions (Zegers et al., 2020), which could be achieved by including additional terms in the local total variance (Sobol' and Kucherenko, 2010).

Consider a transformation f and a vector θ with k parameters, which provides a metric Ψ describing model output:

(2) Ψ = f ( θ ) , f : R k R .

Given a sample point θ* in the parameter space, the gradient for metric Ψ and parameter θj around this point – i.e., Ψθj|θ* – is considered a measure of local sensitivity. In this work, we follow Rakovec et al. (2014) and compute such a gradient using a forward, finite difference approach with a 1 % change in the parameter value:

(3) Ψ θ j | θ * = Ψ θ j * + 0.01 θ j * - Ψ θ j * 0.01 θ j * .

In Eq. (3), Ψ(θ*) is a signature measure of hydrologic behavior, formulated by contrasting model output at the point θ* with that obtained from a reference parameter set θref in the grid cell of interest (see Sect. 3.3.2 for details). The first-order sensitivity measure for the jth parameter is calculated at each sample point as

(4) S j L = Ψ θ j | θ * 2 s j 2 V L θ * .

Here, sj2 is the a priori parameter variance of the jth parameter, and VL(θ*) is the linearized local variance:

(5) V L θ * = j = 1 k Ψ θ j | θ * 2 s j 2 .

Finally, sj2 is obtained from the variance of a uniform distribution (Rakovec et al., 2014), which is 112θj,max-θj,min2.

The first-order sensitivity indices vary between zero and one, and the sum of first-order sensitivities from all parameters at each sampling point is equal to one. Local sensitivities can be examined through their cumulative frequency distribution across the parameter space or by computing a specific statistical property. Here, we quantify the relative contribution of a specific parameter using the area above the curve of the full frequency distribution:

(6) IS j L = 1 - 0 1 F S j L d S j L .

3.4 Performance metrics

We use eight model evaluation metrics to quantify the sensitivity of simulated hydrological processes to variations in model parameters. The notation, brief description, mathematical formulation and physical process associated with each metric are detailed in Table 4. These metrics are computed by contrasting model output from sampling points produced for DELSA with a reference national-scale dataset with simulated states and fluxes obtained from the national water balance database (DGA, 2018, 2019a, b) for the historical period from 1985 to 2015. This dataset was developed by running the VIC model at the same grid discretization as employed here (i.e., 0.05×0.05) and using a combination of CR2MET version 2.0, ERA-Interim and ERA5 output as meteorological forcings. The spatially distributed parameter fields for our reference simulation were developed via parameter regionalization, based on the similarity between possible donor catchments – whose parameters were calibrated individually (Vásquez et al., 2021) – and each grid cell across the domain, following Beck et al. (2016). The reader is referred to Vásquez et al. (2021) and DGA (2018, 2019a, b; in Spanish) for more details on individual model calibration and parameter regionalization procedures used to generate the reference simulation.

Table 4Parameter sensitivity metrics used in this study.

The abbreviations used in the table are as follows: N – number of time steps; Qt – flow for time step t; Zt – flow transformed for time step t; Qm1m1 percentile flow of simulated flow duration curve, where m1=70; Qm2m2 percentile flow of simulated flow duration curve, where m2=30; R – grid-averaged mean annual runoff; P – grid-averaged mean annual precipitation; SWEt – snow water equivalent for time step t; Subl – grid-averaged mean annual sublimation; Transp – grid-averaged mean annual transpiration.

Download Print Version | Download XLSX

Four evaluation metrics are formulated from runoff time series. The first objective function is the root-mean-square error (RMSE), which is a standard metric that emphasizes high flows. The second metric selected is the transformed-root-mean-square error (TRMSE), for which the simulated and observed runoff time series are transformed using a Box–Cox transformation to emphasize low flows (Misirli et al., 2003). The third objective function is the flow duration curve (FDC) mid-segment slope difference (FMS), which represents the variability (or flashiness) of the flow magnitudes; thus, it measures how well a model captures the distribution of the mid-level flows. A steep slope of the FDC indicates flashiness of the streamflow response, whereas a flatter curve indicates a relatively damped response and a higher storage (Yadav et al., 2007; Casper et al., 2012). The fourth evaluation metric is the runoff ratio difference (RR), considered as a measure of the general water balance and, therefore, as a signature of the evapotranspiration model component (Mendoza et al., 2015b).

We use two metrics to characterize snow cover processes: the difference in the long-term simulated peak SWE (PeakSWE), which is an integrated measure of processes occurring during the snowfall season, and the difference in snow cover duration, quantified by the number of days with snow on the ground (SnowLength; Mizukami et al., 2014). Finally, we include two metrics based on evaporation fluxes: the sublimation difference (SUBL), which emphasizes the net sublimation from the snowpack surface, and the plant transpiration difference (TRANSP).

3.5 Experimental setup

We apply the DELSA method in 5574 grid cells across continental Chile, which are contained within the 101 catchments described in Sect. 2. In each grid cell, hydrologic model simulations are conducted at 3-hourly time steps for a 12-year period (April 1999–March 2011), with the first 2 years used to initialize model states (i.e., spin-up period). The model is run in full energy balance mode, which means that both energy and water balances are solved, and 3-hourly outputs are aggregated to daily time steps for subsequent analyses.

In this study, we use the Latin hypercube sampling (LHS) method to obtain 200 sample points across the parameter space, for which first-order sensitivity indices are computed. LHS is an efficient simulation technique, especially suitable for statistical and sensitivity calculations (Vořechovský, 2015). To stratify the parameter space, we sample uniformly in a 43-dimensional hypercube and then map onto the parameter space using the inverse cumulative distribution function of each parameter's prior. As DELSA is used here, all parameter distribution functions are assumed to be uniform. The computational cost of applying DELSA at each cell is Nl(k+1)=8800 model runs, where Nl is the number of sample points (200), and k is the number of parameters (43); hence, the total number of models runs required for this study is 49 051 200.

In this paper, a parameter is considered redundant or insensitive when the median value of the integrated first-order sensitivity index (i.e., median ISjL) across all grid cells is smaller than 0.05 for at least seven of the eight evaluation metrics listed in Table 4. Parameter sensitivity results are also examined per metric (Sect. 3.4) and grid cell climate type based on the aridity index (Table 5; UNEP, 1997; Verbist et al., 2010).

Table 5Climate classification used to group model grid cells.

Download Print Version | Download XLSX

Finally, we use the Spearman rank correlation coefficient (rs) to measure the degree of association between the parameter sensitivities and physiographic/hydroclimatic characteristics listed in Table 2. These attributes were chosen due to their ability to improve the prediction of hydrological signatures (Addor et al., 2018) and because they are relatively easy to obtain.

4 Results and discussion

4.1 Intra- and inter-basin variability in parameter sensitivities

Figure 3 shows the cumulative distribution functions (CDFs) of sensitivity indices SjL for evaluation metric–parameter combinations across six hydroclimatically different catchments – also displayed in Figs. 1 and 2. Each panel located in the first six rows includes the CDFs of all the grid cells contained in a specific basin, and the last row in Fig. 3 comprises the CDFs of all grid cells contained in the six basins. The results reveal high parametric sensitivities for RMSE–INFILT in basins located in northern Chile (arid regime) and lower sensitivities along central and southern catchments, whereas the opposite behavior is observed for TRMSE–DsMAX, FMS–Ds and RR–LAI (i.e., increasing sensitivities towards the south). Such dependence between hydroclimatic basin characteristics and parameter sensitivities was also reported by Demaria et al. (2007). Moreover, Gou et al. (2020) found that sensitivities were strongly related to environmental characteristics, including climate, vegetation, soil and topographic features. Figure 3 also enables the comparison between intra-catchment (top six rows) and inter-catchment (bottom row) variability in parameter sensitivities. For the sample of basins included here, one can note that inter-basin variability in sensitivities is larger than intra-basin variability in runoff-related metrics. Nevertheless, for some metric–parameter combinations, intra-catchment variability is comparable to inter-basin variations in parametric sensitivities. For example, the spread in the CDFs displayed for SUBL–z0SNOW at basin 8317001 (Biobío River at Rucalhue) – characterized by a wet hydroclimate – is comparable to the spread arising from all basins (see same column, last row).

Figure 3Comparison of cumulative frequency distributions of first-order DELSA indices (SjL) across six hydroclimatically different basins (displayed in different rows, sorted by latitude). The location and seasonal cycles for these basins are displayed in Figs. 1 and 2, respectively. Results are displayed for the most sensitive parameter associated with each evaluation metric (displayed in different columns); thus, each panel (except those in the last row) comprises the CDFs of all grid cells contained in a specific basin, for a particular metric–parameter combination. The most sensitive parameter was determined based on the median sensitivity index ISjL from all the grid cells contained in the study domain (see text for details). The number next to each basin code at the top of this figure is the catchment-scale aridity index.


To further illustrate intra-basin differences in parameter sensitivities, Figs. 4 and 5 show the spatial distribution of the ISjL indices for the leaf area index (LAI) and the snow albedo parameter ALB THA, respectively, over a cluster of subhumid and humid basins located in southern Chile. For the LAI parameter (Fig. 4), a west–east gradient in ISjL is observed for RMSE (high flows), TRMSE (low flows) and FMS (flashiness of runoff), with increasing sensitivity to LAI variations towards the coast, whereas an inverse pattern is observed for the same metrics and ALB THA (i.e., larger sensitivities towards the Andes, Fig. 5). For PeakSWE, SUBL and – to a smaller degree – SnowLength, LAI yields larger sensitivities in vegetated areas, where snow accumulates during winter (Fig. 4), matching those locations where forest is the dominant land cover type, which is also the only vegetation class with an overstory (e.g., trees). Notably, very large variations in PeakSWE sensitivities to LAI are observed over relatively short distances due to differences among grid cells with respect to the fraction of land cover defined as forest. Such dependence among SWE sensitivities, LAI and canopy fractions was also reported by Bennett et al. (2018). Figure 4 also shows that LAI does not yield a clear sensitivity pattern in RR and TRANSP throughout this subdomain, although ISjL values are higher for RR. For this metric, there are spatial singularities where the sensitivity is minimal or null, as the fraction of ground cover defined as bare soil in these areas increases considerably, reaching up to 100 % bare soil (LAI  0) in some grid cells.

Figure 4Spatial distribution of integrated first-order DELSA sensitivity indices (ISjL) for the leaf area index (LAI) across a humid subdomain located in southern Chile. Results are displayed for eight sensitivity metrics: (a) RMSE, (b) TRMSE, (c) FMS, (d) RR, (e) PeakSWE, (f) SnowLength, (g) SUBL and (h) TRANSP.

Figure 5Same as in Fig. 4 but for the base number in the snow albedo function ALB THA (melt).

The results presented in Fig. 5 reinforce the idea that hard-coded parameters should be exposed to users (Mendoza et al., 2015a; Cuntz et al., 2016). In particular, Fig. 5 shows the large effects of ALB THA variations on SnowLength (with a very pronounced east–west gradient) and, to a smaller degree, on PeakSWE and SUBL. ALB THA also affects runoff-based metrics along the Andes, especially simulated high (RMSE) and low (TRMSE) flows.

4.2 Identification of the most sensitive parameters

Figure 6 displays box plots comprising ISjL results from all grid cells in the study domain, for each parameter and evaluation metric (displayed in different panels). The results show that 72 % of the parameters analyzed (i.e., 31) yielded small sensitivity indices for the metrics examined here. Conversely, a suite of 12 sensitive parameters are associated with soil (INFILT, Ds, DsMAX, Ws, Expt2, Depth2 and Depth3), snow (NEW ALB, ALB THA and ALB AA) and vegetation (Rmin and LAI) processes. Figure 7 shows the spatial variability in ISjL for the 12 parameters identified as the most sensitive across the 101 basins of continental Chile.

Figure 6Box plots comprising integrated first-order DELSA sensitivity indices (ISjL) from all modeling units (5574 grid cells). Results are displayed for all parameters (x axis) and sensitivity metrics, which are presented in different panels: (a) RMSE, (b) TRMSE, (c) FMS, (d) RR, (e) PeakSWE, (f) SnowLength, (g) SUBL and (h) TRANSP.


Figure 7Integrated first-order DELSA sensitivity indices for all grid cells within our study basins. The results are displayed only for the 12 most sensitive parameters as well as their associated most impacted metrics.

For the case of high flows (RMSE), low flows (TRMSE) and flashiness of runoff (FMS), the parameters identified as sensitive are INFILT, Ds, DsMAX, Ws, Expt2, Depth2 and Depth3 (see top three panels in Fig. 6). The INFILT parameter controls the shape of the variable infiltration capacity curve (Zhao et al., 1980; Wood et al., 1992) and, thus, the partitioning of rainfall or snowmelt into infiltration and surface runoff. A higher INFILT value yields less infiltration and higher surface runoff. The RMSE and TRMSE metrics are particularly sensitive to INFILT, indicating a key role in the generation and timing of high and low flows. This parameter has been identified as sensitive in all of the studies listed in Table 1. DsMAX is the maximum velocity of baseflow, while Ds and Ws are the fraction of DsMAX and the fraction of the maximum soil moisture content in the third layer, respectively, where non-linear baseflow occurs. These three parameters are involved in the ARNO formulation of subsurface runoff (Franchini and Pacciani, 1991; Todini, 1996), controlling the speed of baseflow release from the third soil layer (Liang et al., 1994) and, specifically, the non-linear part of the baseflow generation function. The sensitivity indices found for these parameters are consistent with the high sensitivity measures reported by Mendoza et al. (2015b), Melsen et al. (2016) and Wi et al. (2017).

The Expt2 parameter is an exponent of the Brooks–Corey relationship (Brooks and Corey, 1964) and controls the hydraulic conductivity between the second and third soil layers. A small value for the Expt2 parameter increases inter-layer drainage for the same soil moisture content and, therefore, increases baseflow generation. The Depth2 parameter is the thickness of the second soil layer. In general, thicker soil layers slow seasonal peak flows and increase water loss due to evapotranspiration (Xie et al., 2007). It should be noted that the parameter Depth2 has been identified as highly sensitive by many authors (Demaria et al., 2007; Mendoza et al., 2015b; Wi et al., 2017; Gou et al., 2020; Lilhare et al., 2020; Yeste et al., 2020; Melsen and Guse, 2021). Finally, Depth3 is the thickness of the third layer of soil, and the large sensitivities obtained here agree with the results reported by Mendoza et al. (2015b) and Wi et al. (2017).

The results in Fig. 6 show that Expt2 and Depth2 also provide large sensitivities for metrics focused on evaporative fluxes (i.e., RR and TRANSP). Other parameters that are relevant for these processes are LAI and the minimum stomatal resistance (Rmin). Indeed, Chaney et al. (2015) reported high sensitivity of annual flow biases to variations in Rmin. LAI is a dimensionless quantity that characterizes intra-annual variations in plant canopies, and it is defined as the one-sided green leaf area per unit ground surface area. On the other hand, Rmin is one of the parameters that control canopy resistance when computing transpiration from each vegetation class, following the formulations proposed by Blondin (1991) and Ducoudré et al. (1993). Both LAI and Rmin provide null sensitivities if the land cover type is bare ground; however, Rmin can also produce null sensitivities in vegetated grid cells.

Figure 6 reveals the large influence of hard-coded parameters on PeakSWE, SnowLength and SUBL – in particular, NEW ALB, ALB THA and ALB AA. The NEW ALB parameter is the new snow surface albedo, which controls the reflection of solar radiation and, therefore, the energy exchange between the atmosphere, forest canopy and the surface layer of the snowpack (Andreadis et al., 2009). Additionally. The ALB AA and ALB THA parameters represent the albedo decay in the accumulation and melting season in the snow albedo curve, respectively (USACE, 1956). These seasons are defined based on the absence or presence of liquid water in the surface snow cover. These results correspond well with the high sensitivities reported by Mendoza et al. (2015b) for these three hard-coded parameters. Finally, the snow surface roughness length (z0SNOW) also affects sublimation rates across Andean subdomains.

4.3 What drives parameter sensitivities across different hydroclimates?

Figure 8 shows the Spearman rank correlation coefficient (rs) between parameter sensitivities and a suite of climatic, topographic, land cover and soil-related grid cell attributes described in Table 2. The magnitude and sign of the correlation quantifies how each sensitivity index varies with a given geophysical attribute. The results show that the magnitude of the correlation varies depending on the metric–parameter combination, with the maximum correlations generally found for soil parameters, such as DsMAX and Ws, with precipitation (rs=0.91) and the aridity index (rs=-0.91) for the FMS function (flashiness of the flow). For this example, the simulation of flow flashiness is highly sensitive to DsMAX and Ws in the wet region but insensitive in the arid area. Conversely, the minimum correlations are found for snow parameters. It should be noted that a weak correlation indicates that there is less spatial pattern in sensitivity; however, the magnitude of the sensitivity index can be high or low. For example, NEW ALB is a highly sensitive parameter across the domain (Fig. 6).

Figure 8The Spearman rank correlation coefficient between integrated first-order DELSA sensitivity indices ISjL and grid cell characteristics. Results are displayed only for the four most sensitive parameters affecting each metric. The crosses indicate correlations with p values larger than 0.05.


The results in Fig. 8 also indicate that high correlations (either positive or negative) are mainly associated with climate indices, which exert a stronger influence compared with the remaining attribute classes. These strong dependencies of the parametric sensitivity on climate variables are somewhat expected, as some combinations of hydrological signatures and parameters inherit strong spatial climate patterns (Addor et al., 2018); compare, for example, the aridity index, shown in Fig. 1d, with DsMAX–FMS and DsMAX–RMSE, shown in Fig. 7. Among climate descriptors, the aridity index, mean annual precipitation and relative humidity yield the highest correlations, and temperature exhibits a relatively lower influence on parametric sensitivity; this result was confirmed with additional correlation analyses including only grid cells with mean annual temperatures below 5 and 2 C (not shown). The lowest correlations are obtained for mean slope (topographic attribute), shrub fraction (land cover attribute) and mean clay content of soil (soil attribute). The key influence of climatic conditions on hydrological behavior is not new, as aridity is commonly regarded as the main driver of water partitioning at the land surface (Budyko, 1974; Hrachowitz et al., 2013).

Figure 8 shows that the extent to which parametric sensitivities are related to grid cell attributes depends on the target evaluation metric (i.e., runoff, evaporative processes and snow processes), with the three distinct groups containing the same influential parameters. In the following subsections, we discuss the results based on these groups, with emphasis on spatial patterns and process interpretation across our study basins. Table 6 summarizes, for each evaluation metric (i.e., physical process to be represented) and climatic zone (using the classification from Table 5), the three most important parameters. Hence, the lists contained therein can be used to guide the selection of parameters for hydrologic model calibration, based on the hydroclimatic regime and target process that modelers would like to represent.

4.3.1 Runoff-oriented metrics

The results presented in Fig. 7 (see RMSE and TRMSE) and Fig. 8 (RMSE) show a direct relationship between the sensitivities provided by INFILT, as well as the degree of aridity, especially in semiarid to hyperarid subdomains. The runoff-oriented metrics are also sensitive to baseflow generation parameters Ds, Ws and DsMAX in most basins – with a similar spatial distribution of ISjL values – except for those located in the north and in some areas of southern Patagonia, where climatic conditions are arid or hyperarid. In basins located in the extreme north, small sensitivities can be attributed to local climate characteristics: most precipitation events in that area occur in summer (i.e., December–March) due to orographic rains caused by air masses coming from the Amazon region, and there is usually little recharge to the aquifers. Additionally, the third soil layer in these basins generally does not reach saturation; therefore, runoff simulations in those areas are insensitive to variations in Ds and Ws because the non-linear part of the baseflow function is only activated when the moisture storage in the third layer exceeds a threshold (Gou et al., 2020). Because of the dependency of DsMAX on precipitation, this parameter could be playing a key role in the baseflow generation processes over Andean regions (Fig. 7). Finally, a similar spatial distribution of integrated first-order sensitivities for Ds, Ws and DsMAX is expected, as they all focus on baseflow generation (see DsMAX–FMS, DsMAX–RMSE, Ws–TRMSE, Ws–FMS and Ds–TRMSE in Fig. 7).

Expt2 is identified as sensitive for runoff-oriented metrics in basins with semiarid to hyperarid climates, characterized by small annual precipitation amounts and permanent water stress. In these hydroclimatic regimes, there is usually not enough water to reach the third soil layer, so water is stored in the second layer and drainage is mainly controlled by Expt2, affecting the vertical redistribution of soil moisture (FMS) and low flows (TRMSE), as shown in Fig. 7. Depth2 provides large runoff sensitivities in dry-subhumid to hyperarid hydroclimatic regimes, for the same reasons as Expt2. Variations in the depth of the second soil layer change the soil moisture of the layer, and higher (lower) values of Depth2 for the same volumetric water content produce lower (higher) soil moisture, affecting drainage between soil layers.

Table 6Summary of the most sensitive VIC parameters found for each metric (rows) and climatic type. The three most important parameters are determined based on the median of the integrated first-order DELSA sensitivity indices and are sorted by ranking (i.e., first, second and third most sensitive).

Download Print Version | Download XLSX

Finally, Depth3 provides large sensitivities for all runoff-oriented metrics, with similar spatial patterns to Ws, Ds and DsMAX, although to a smaller degree (Fig. 7). Depth3 is particularly sensitive in humid-subhumid and humid catchments, suggesting a direct relationship with mean annual precipitation or with the size and intermittency of storms (Abdulla and Lettenmaier, 1997). In these climatic domains, periodic heavy-rainfall events enable a continuous recharge of the second and third soil layers – which may reach saturation – and, thus, a constant baseflow generation that affects runoff response and the retention time of soil moisture, producing higher baseflow during wet seasons (Shi et al., 2008). Interestingly, in humid subdomains, baseflow parameters yield high sensitivities in both rainfall- and snowfall-dominated grid cells, although ALB THA emerges as the most important parameter for RMSE and TRMSE in snowfall-dominated grid cells (not shown).

4.3.2 Evaporative processes

The evaluation metrics associated with these processes are RR (a measure of the overall water balance) and TRANSP (plant transpiration), with LAI, Rmin, Expt2 and Depth2 being the most important parameters.

Figure 7 shows a pronounced spatial variability in LAI sensitivities across a large domain that comprises very different land cover types. One can note that LAI yields high sensitivities for nearly all hydroclimatic regimes, as this parameter controls the evaporation from the canopy layer and canopy transpiration. In hyperarid climates, the LAI is usually less important, given the permanent water stress common for grid cells with bare soil. In summary, LAI is influential wherever vegetation exists, regardless of the prevailing hydroclimatic regime.

The parameter Rmin yields parametric sensitivities across humid-subhumid and humid areas (Fig. 7 and Table 6). In the canopy resistance process, there is a stomatal resistance multiplier, gsm[n], defined as a soil moisture stress factor that depends on the water in the root zone for the nth surface cover class. Thus, when the soil moisture in layer n is less than the fraction of the moisture content at the wilting point, the value of gsm[n] is zero, whereas when the soil moisture is greater than the fractional content of soil moisture at the critical point (∼70 % of field capacity), the value of gsm[n] is one. For the intermediate condition, gsm[n] values vary linearly with soil moisture in that layer, which explains why Rmin provides high sensitivities in very humid (i.e., large precipitation) climates.

Finally, our results show that Expt2 and Depth2 yield large sensitivities to RR and TRANSP in all hydroclimatic regimes, as they affect the soil moisture content in layer 2, which indirectly affects the gsm[n] factor in the canopy resistance formulation. These parameters show a lower relative sensitivity in humid/subhumid and humid climates, as the Rmin parameter becomes more relevant when there is no soil moisture stress (i.e., gsm[n]∼1).

4.3.3 Snow processes

Figure 7 shows that NEW ALB, ALB THA and ALB AA yield high sensitivities throughout the study domain, especially in areas where snow processes dominate hydrological responses. In particular, the NEW ALB parameter is important throughout the domain and reaches the highest values for snow-oriented evaluation metrics. Additionally, the results in Fig. 7 show that ALB AA and ALB THA dominate snow responses in different domains: the ALB THA parameter yields large sensitivities in humid and subhumid mountainous areas located south of 34 S, with large effects on the snow season length and the maximum SWE accumulation, whereas ALB AA shows greater sensitivity for the other climatic regimes, affecting SnowLength and sublimation in semiarid, colder environments in northern Chile (26–29 S).

5 Conclusions

In this study, we have revisited parameter sensitivities in the Variable Infiltration Capacity hydrological model. To this end, we have implemented the DELSA method at every 0.05×0.05 (∼5 km) grid cell contained in 101 basins across continental Chile (i.e., a total of 5574 grid cells), spanning a broad diversity of hydroclimatic (from hyperarid to humid) and physiographic (e.g., topography and land cover) conditions. Our experiments consider a suite of 43 parameters included in soil, vegetation and snow process representations, with 3 of these corresponding to hard-coded parameters (i.e., not exposed to model users). We use eight model evaluation metrics that account for runoff components, evapotranspiration and snow processes, and we conduct correlation analyses to disentangle relationships between parametric sensitivities and pixel-scale attributes. The main findings of this study are as follows:

  • A total of 31 of the 43 (i.e., 72 %) parameters yield little or no sensitivity, most of which correspond to soil and vegetation processes. Therefore, calibrating such parameters will lead to minimal improvements in system representations with considerable computational costs.

  • The three model evaluation metrics focused on snow accumulation and ablation processes were found to be highly sensitive to hard-coded parameters. Exposing these parameters will certainly expand our abilities to perform extensive analysis and increase our opportunities to improve model fidelity and characterize model uncertainty.

  • For some evaluation metrics, the climate attributes examined here are highly correlated with parameter sensitivities, which therefore inherit spatial patterns observed in climate variables across the territory. In particular, mean annual precipitation and the aridity index are highly correlated with Ds, Ws and DsMAX sensitivities when examining RMSE, TRMSE and FMS. Unexpectedly, temperature yields a relatively lower influence among climate descriptors, even for metrics and parameters associated with snow processes. The rest of the attributes (topographic, soil and land cover) provided generally low correlations and, therefore, small predictive power on parameter sensitivities.

  • Parametric sensitivities are strongly related to the climate types in the case study basins. In humid environments, the most important parameters are related to the third soil layer (Ws, Ds, DsMAX and Depth3) and vegetation (Rmin); in arid regimes, the most influential parameters are associated with the first soil layers (INFILT, Expt2 and Depth2).

  • In snow-dominated areas, the hard-coded parameters NEW ALB, ALB THA and ALB AA provide large sensitivities to maximum SWE, snow season length and sublimation.

  • The leaf area index (LAI) is a crucial parameter wherever there is vegetation on the ground. Although such conditions are more frequent in humid environments, the relevance of this parameter depends on vegetation characteristics rather than the underlying climatic conditions.

This study reaffirms overparameterization issues in the VIC model (e.g., Demaria et al., 2007; Gou et al., 2020) as well as the fact that relative parameter importance varies depending on the specific metric or variable analyzed (Chaney et al., 2015; Bennett et al., 2018; Yeste et al., 2020; Melsen and Guse, 2021) and the physiographic or climatic site characteristics (Liang and Guo, 2003; Demaria et al., 2007; Lilhare et al., 2020). However, the results and conclusions reported here are not directly comparable to previous studies due to differences in the experimental designs and domains of interest. In particular, our study considers (1) a large number of soil, vegetation and snow parameters – only comparable to Bennett et al. (2018), who included 46 parameters (excluding snow processes) over the semiarid Colorado River basin; (2) a larger number of process-oriented metrics (compared with the previous efforts listed in Table 1); and (3) a very large sample of grid cells at a relatively high (∼5 km) horizontal resolution, spanning very different physiographic characteristics and hydroclimatic conditions. Hence, our study contributes to the existing literature by providing guidance on relevant VIC parameters for a suite of target processes and climate types.

Future studies aiming at improving spatial calibration density or parameter regionalization techniques using VIC – or any similar hydrology or land surface model – could incorporate this information to define spatially varying target parameters and to examine the extent to which the spatial patterns in parameter sensitivities relate to calibrated parameter fields. Finally, the strong correlations found here between parameter sensitivities and hydroclimatic properties reaffirm the need to incorporate periods with contrasting climate characteristics in sensitivity analysis and calibration strategies in order to achieve more credible hydrologic model simulations under changing climatic conditions.

Code availability

The codes used in this study are available from the corresponding author upon reasonable request.

Data availability

Land cover descriptors, meteorological forcings and reference VIC model outputs for all grid cells were obtained from the national water balance database (DGA, 2018, 2019a, b); this information may be requested from (last access: 13 June 2022). Other grid cell attributes were obtained from the United States Geological Survey dataset (; USGS, 2022), the CR2MET dataset (; Boisier et al., 2018) and SoilGrids250m 2.0 ( (last access: 13 June 2022) (Hengl et al., 2017). Catchment-scale hydrometeorological data were obtained from the CAMELS-CL dataset (Alvarez-Garreton et al., 2018).

Author contributions

All authors were involved in the conceptualization of this study. UMS and PAM designed the methodology and analysis framework and drafted the paper. UMS configured the VIC model, conducted simulations, analyzed the results and created the figures. AJN and NM provided insights into the sensitivity analysis results. All authors discussed the results and contributed to writing, reviewing and editing the manuscript.

Competing interests

The contact author has declared that neither they nor their co-authors have any competing interests.


Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.


We thank Eduardo Muñoz-Castro and Nicolás Vásquez for their advice and assistance in setting up model simulations as well as Ximena Vargas and Miguel Lagos for their suggestions on earlier versions of this paper. Finally, we thank the editor (Nunzio Romano), Neil Grigg and one anonymous reviewer for their constructive comments, which helped to improve this paper.

Financial support

This research has been supported by the Fondo Nacional de Desarrollo Científico y Tecnológico (project Fondecyt 11200142) and by CONICYT/PIA (project AFB180004). This research was partially supported by the supercomputing infrastructure of the NLHPC (ECM-02). The National Center for Atmospheric Research is a major facility sponsored by the National Science Foundation under cooperative agreement no. 1852977.

Review statement

This paper was edited by Nunzio Romano and reviewed by Neil Grigg and one anonymous referee.


Abbaspour, K. C., Rouholahnejad, E., Vaghefi, S., Srinivasan, R., Yang, H., and Kløve, B.: A continental-scale hydrology and water quality model for Europe: Calibration and uncertainty of a high-resolution large-scale SWAT model, J. Hydrol., 524, 733–752,, 2015. 

Abdulla, F. A. and Lettenmaier, D. P.: Development of regional parameter estimation equations for a macroscale hydrologic model, J. Hydrol., 197, 230–257,, 1997. 

Addor, N. and Melsen, L. A.: Legacy, Rather Than Adequacy, Drives the Selection of Hydrological Models, Water Resour. Res., 55, 378–390,, 2019. 

Addor, N., Nearing, G., Prieto, C., Newman, A. J., Le Vine, N., and Clark, M. P.: A Ranking of Hydrological Signatures Based on Their Predictability in Space, Water Resour. Res., 54, 8792–8812,, 2018. 

Al Nakshabandi, G. and Kohnke, H.: Thermal conductivity and diffusivity of soils as related to moisture tension and other physical properties, Agric. Meteorol., 2, 271–279,, 1965. 

Alvarez-Garreton, C., Mendoza, P. A., Boisier, J. P., Addor, N., Galleguillos, M., Zambrano-Bigiarini, M., Lara, A., Puelma, C., Cortes, G., Garreaud, R., McPhee, J., and Ayala, A.: The CAMELS-CL dataset: catchment attributes and meteorology for large sample studies – Chile dataset, Hydrol. Earth Syst. Sci., 22, 5817–5846,, 2018. 

Andreadis, K. M. and Lettenmaier, D. P.: Trends in 20th century drought over the continental United States, Geophys. Res. Lett., 33, 1–4,, 2006. 

Andreadis, K. M., Storck, P., and Lettenmaier, D. P.: Modeling snow accumulation and ablation processes in forested environments, Water Resour. Res., 45, W05429,, 2009. 

Arheimer, B., Pimentel, R., Isberg, K., Crochemore, L., Andersson, J. C. M., Hasan, A., and Pineda, L.: Global catchment modelling using World-Wide HYPE (WWH), open data, and stepwise parameter estimation, Hydrol. Earth Syst. Sci., 24, 535–559,, 2020. 

Bastidas, L. A., Gupta, H. V., Sorooshian, S., Shuttleworth, W. J., and Yang, Z. L.: Sensitivity analysis of a land surface scheme using multicriteria methods, J. Geophys. Res., 104, 19481–19490, 1999. 

Bastidas, L. A., Hogue, T. S., Sorooshian, S., Gupta, H. V., and Shuttleworth, W. J.: Parameter sensitivity analysis for different complexity land surface models using multicriteria methods, J. Geophys. Res., 111, D20101,, 2006.  

Beck, H. E., van Dijk, A. I. J. M., de Roo, A., Miralles, D. G., McVicar, T. R., Schellekens, J., and Bruijnzeel, L. A.: Global-scale regionalization of hydrologic model parameters, Water Resour. Res., 52, 3599–3622,, 2016. 

Bennett, K. E., Urrego Blanco, J. R., Jonko, A., Bohn, T. J., Atchley, A. L., Urban, N. M., and Middleton, R. S.: Global Sensitivity of Simulated Water Balance Indicators Under Future Climate Change in the Colorado Basin, Water Resour. Res., 54, 132–149,, 2018. 

Berg, P., Donnelly, C., and Gustafsson, D.: Near-real-time adjusted reanalysis forcing data for hydrology, Hydrol. Earth Syst. Sci., 22, 989–1000,, 2018. 

Blondin, C.: Parameterization of Land-Surface Processes in Numerical Weather Prediction, L. Surf. Evaporation, Springer, 31–54,, 1991. 

Bohn, T. J. and Vivoni, E. R.: Process-based characterization of evapotranspiration sources over the North American monsoon region, Water Resour. Res., 52, 358–384,, 2016. 

Boisier, J. P., Alvarez-Garretón, C., Cepeda, J., Osses, A., Vásquez, N., and Rondanelli, R.: CR2MET: A high-resolution precipitation and temperature dataset for hydroclimatic research in Chile, ADS [data set], (last access: 13 June 2022), 2018. 

Brooks, R. H. and Corey, A. T.: Hydraulic Properties of Porous Media and Their Relation to Drainage Design, T. ASAE, 7, 0026–0028,, 1964. 

Budyko, M. I.: Climate and Life, Academic Press, ISBN 0121394506, ISBN 13 978-0121394509, 1974. 

Casper, M. C., Grigoryan, G., Gronz, O., Gutjahr, O., Heinemann, G., Ley, R., and Rock, A.: Analysis of projected hydrological behavior of catchments based on signature indices, Hydrol. Earth Syst. Sci., 16, 409–421,, 2012. 

Cayan, D. R., Das, T., Pierce, D. W., Barnett, T. P., Tyree, M., and Gershunov, A.: Future dryness in the southwest US and the hydrology of the early 21st century drought, P. Natl. Acad. Sci. USA, 107, 21271–21276,, 2010. 

Chaney, N. W., Herman, J. D., Reed, P. M., and Wood, E. F.: Flood and drought hydrologic monitoring: The role of model parameter uncertainty, Hydrol. Earth Syst. Sci., 19, 3239–3251,, 2015. 

Chawla, I. and Mujumdar, P. P.: Isolating the impacts of land use and climate change on streamflow, Hydrol. Earth Syst. Sci., 19, 3633–3651,, 2015. 

Chegwidden, O. S. S., Nijssen, B., Rupp, D. E. E., Arnold, J. R. R., Clark, M. P. P., Hamman, J. J. J., Kao, S. C. S. C., Mao, Y., Mizukami, N., Mote, P. W., Pan, M., Pytlak, E., and Xiao, M.: How do modeling decisions affect the spread among hydrologic climate change projections? Exploring a large ensemble of simulations across a diversity of hydroclimates, Earth's Future, 7, 623–637,, 2019. 

Chen, F., Barlage, M., Tewari, M., Rasmussen, R., Jin, J., Lettenmaier, D., Livneh, B., Lin, C., Miguez-Macho, G., Niu, G., Wen, L., and Yang, Z.: Modeling seasonal snowpack evolution in the complex terrain and forested Colorado Headwaters region: A model intercomparison study, J. Geophys. Res.-Atmos., 119, 13795–13819,, 2014. 

Cherkauer, K. A., Bowling, L. C., and Lettenmaier, D. P.: Variable infiltration capacity cold land process model updates, Global Planet. Change, 38, 151–159,, 2003. 

Chipman, H. A., George, E. I., and McCulloch, R. E.: BART: Bayesian additive regression trees, Ann. Appl. Stat., 4, 266–298,, 2010. 

Clark, M. P., Nijssen, B., Lundquist, J. D., Kavetski, D., Rupp, D. E., Woods, R. A., Freer, J. E., Gutmann, E. D., Wood, A. W., Brekke, L. D., Arnold, J. R., Gochis, D. J., and Rasmussen, R. M.: A unified approach for process-based hydrologic modeling: 1. Modeling concept, Water Resour. Res., 51, 2498–2514,, 2015. 

Clark, M. P., Bierkens, M. F. P., Samaniego, L., Woods, R. A., Uijlenhoet, R., Bennett, K. E., Pauwels, V. R. N., Cai, X., Wood, A. W., and Peters-Lidard, C. D.: The evolution of process-based hydrologic models: historical challenges and the collective quest for physical realism, Hydrol. Earth Syst. Sci., 21, 3427–3440,, 2017. 

Cosby, B. J., Hornberger, G. M., Clapp, R. B., and Ginn, T. R.: A Statistical Exploration of the Relationships of Soil Moisture Characteristics to the Physical Properties of Soils, Water Resour. Res., 20, 682–690,, 1984. 

C3S and Copernicus Climate Change Service (C3S): ERA5: Fifth generation of ECMWF atmospheric reanalyses of the global climate, C3S,!/home (last access: 20 January 2018), 2017. 

Cuntz, M., Mai, J., Samaniego, L., Clark, M., Wulfmeyer, V., Branch, O., Attinger, S., and Thober, S.: The impact of standard and hard-coded parameters on the hydrologic fluxes in the Noah-MP land surface model, J. Geophys. Res.-Atmos., 121, 10676–10700,, 2016. 

DeChant, C. M. and Moradkhani, H.: Toward a reliable prediction of seasonal forecast uncertainty: Addressing model and initial condition uncertainty with ensemble data assimilation and Sequential Bayesian Combination, J. Hydrol., 519, 2967–2977,, 2014. 

Dee, D. P., Uppala, S. M., Simmons, A. J., Berrisford, P., Poli, P., Kobayashi, S., Andrae, U., Balmaseda, M. A., Balsamo, G., Bauer, P., Bechtold, P., Beljaars, A. C. M., van de Berg, L., Bidlot, J., Bormann, N., Delsol, C., Dragani, R., Fuentes, M., Geer, A. J., Haimberger, L., Healy, S. B., Hersbach, H., Hólm, E. V., Isaksen, L., Kållberg, P., Köhler, M., Matricardi, M., Mcnally, A. P., Monge-Sanz, B. M., Morcrette, J. J., Park, B. K., Peubey, C., de Rosnay, P., Tavolato, C., Thépaut, J. N., and Vitart, F.: The ERA-Interim reanalysis: Configuration and performance of the data assimilation system, Q. J. Roy. Meteorol. Soc., 137, 553–597,, 2011. 

Demaria, E. M., Nijssen, B., and Wagener, T.: Monte Carlo sensitivity analysis of land surface parameters using the Variable Infiltration Capacity model, J. Geophys. Res., 112, 1–15,, 2007. 

DGA: Actualización del Balance Hídrico Nacional, SIT No. 417, (last access: 13 June 2022), 2017. 

DGA: Aplicación de la metodología de actualización del balance hídrico nacional a las macrozonas Norte y Centro, SIT No. 435, (last access: 13 June 2022), 2018. 

DGA: Aplicación de la metodología de actualización del balance hídrico nacional en las cuencas de la macrozona Sur y parte de la Macrozona Austral, SIT No. 441, (last access: 13 June 2022), 2019a. 

DGA: Aplicación de la metodología de actualización del balance hídrico nacional en las cuencas de la parte sur de la macrozona Austral e Isla de Pascua, SIT No. 444, (last access: 13 June 2022), 2019b. 

Do, H. X., Gudmundsson, L., Leonard, M., and Westra, S.: The Global Streamflow Indices and Metadata Archive (GSIM) – Part 1: The production of a daily streamflow archive and metadata, Earth Syst. Sci. Data, 10, 765–785,, 2018. 

Dorman, J. and Sellers, P.: A global climatology of albedo, roughness length and stomatal resistance for atmospheric general circulation models as represented by the simple biosphere model (SiB), J. Appl. Meteorol., 28, 833–855, 1989. 

Ducoudré, N. I., Laval, K., and Perrier, A.: SECHIBA, a New Set of Parameterizations of the Hydrologic Exchanges at the Land-Atmosphere Interface within the LMD Atmospheric General Circulation Model, J. Climate, 6, 248–273,<0248:sansop>;2, 1993. 

Fang, H., Wei, S., Jiang, C., and Scipal, K.: Theoretical uncertainty analysis of global MODIS, CYCLOPES, and GLOBCARBON LAI products using a triple collocation method, Remote Sens. Environ., 124, 610–621,, 2012. 

Fang, H., Jiang, C., Li, W., Wei, S., Baret, F., Chen, J. M., Garcia-Haro, J., Liang, S., Liu, R., Myneni, R. B., Pinty, B., Xiao, Z., and Zhu, Z.: Characterization and intercomparison of global moderate resolution leaf area index (LAI) products: Analysis of climatologies and theoretical uncertainties, J. Geophys. Res.-Biogeo., 118, 529–548,, 2013. 

Foglia, L., Hill, M. C., Mehl, S. W., and Burlando, P.: Sensitivity analysis, calibration, and testing of a distributed hydrological model using error-based weighting and one objective function, Water Resour. Res., 45, W06427,, 2009. 

Franchini, M. and Pacciani, M.: Comparative analysis of several conceptual rainfall-runoff models, J. Hydrol., 122, 161–219,, 1991. 

Friedman, J. H.: Multivariate Adaptive Regression Splines, Ann. Stat., 19, 590–606,, 1991. 

Gates, D. M. and Evans, L. T.: Environmental Control of Plant Growth, Bull. Torrey Bot. Club, 91, 235,, 1964. 

Gharari, S., Clark, M. P., Mizukami, N., Wong, J. S., Pietroniro, A., and Wheater, H. S.: Improving the representation of subsurface water movement in land models, J. Hydrometeorol., 20, 2401–2418,, 2019. 

Ghiggi, G., Humphrey, V., Seneviratne, S. I., and Gudmundsson, L.: GRUN: an observation-based global gridded runoff dataset from 1902 to 2014, Earth Syst. Sci. Data, 11, 1655–1674,, 2019. 

Göhler, M., Mai, J., and Cuntz, M.: Use of eigendecomposition in a parameter sensitivity analysis of the Community Land Model, J. Geophys. Res.-Biogeo., 118, 904–921,, 2013. 

Gou, J., Miao, C., Duan, Q., Tang, Q., Di, Z., Liao, W., Wu, J., and Zhou, R.: Sensitivity Analysis-Based Automatic Parameter Calibration of the VIC Model for Streamflow Simulations Over China, Water Resour. Res., 56, 1–19,, 2020. 

Gou, J., Miao, C., Samaniego, L., Xiao, M., Wu, J., and Guo, X.: CNRD v1.0: A High-Quality Natural Runoff Dataset for Hydrological and Climate Studies in China, B. Am. Meteorol. Soc., 102, E929–E947,, 2021. 

Gupta, H. V., Kling, H., Yilmaz, K. K., and Martinez, G. F.: Decomposition of the mean squared error and NSE performance criteria: Implications for improving hydrological modelling, J. Hydrol., 377, 80–91, 2009. 

Hamman, J. J., Nijssen, B., Bohn, T. J., Gergel, D. R., and Mao, Y.: The variable infiltration capacity model version 5 (VIC-5): Infrastructure improvements for new applications and reproducibility, Geosci. Model Dev., 11, 3481–3496,, 2018. 

Hengl, T., De Jesus, J. M., Heuvelink, G. B. M., Gonzalez, M. R., Kilibarda, M., Blagotić, A., Shangguan, W., Wright, M. N., Geng, X., Bauer-Marschallinger, B., Guevara, M. A., Vargas, R., MacMillan, R. A., Batjes, N. H., Leenaars, J. G. B., Ribeiro, E., Wheeler, I., Mantel, S., and Kempen, B.: SoilGrids250m: Global gridded soil information based on machine learning, PLoS ONE, 12, e0169748,, 2017. 

Hogue, T. S., Bastidas, L., Gupta, H., Sorooshian, S., Mitchell, K., and Emmerich, W.: Evaluation and transferability of the Noah land surface model in semiarid environments, J. Hydrometeorol., 6, 68–84, 2005. 

Hornberger, G. M. and Spear, R. C.: Approach to the preliminary analysis of environmental systems, J. Environ. Manage., 12, 7–18, 1981. 

Hou, Z., Huang, M., Leung, L. R., Lin, G., and Ricciuto, D. M.: Sensitivity of surface flux simulations to hydrologic parameters based on an uncertainty quantification framework applied to the Community Land Model, J. Geophys. Res.-Atmos., 117, 1–18,, 2012. 

Hrachowitz, M., Savenije, H. H. G., Blöschl, G., McDonnell, J. J., Sivapalan, M., Pomeroy, J. W., Arheimer, B., Blume, T., Clark, M. P., Ehret, U., Fenicia, F., Freer, J. E., Gelfan, A., Gupta, H. V., Hughes, D. A., Hut, R. W., Montanari, A., Pande, S., Tetzlaff, D., Troch, P. A., Uhlenbrook, S., Wagener, T., Winsemius, H. C., Woods, R. A., Zehe, E., and Cudennec, C.: A decade of Predictions in Ungauged Basins (PUB) – a review, Hydrolog. Sci. J., 58, 1198–1255,, 2013. 

Huang, M. and Liang, X.: On the assessment of the impact of reducing parameters and identification of parameter uncertainties for a hydrologic model with applications to ungauged basins, J. Hydrol., 320, 37–61,, 2006. 

Lawrence, D. M., Fisher, R. A., Koven, C. D., Oleson, K. W., Swenson, S. C., Bonan, G., Collier, N., Ghimire, B., van Kampenhout, L., Kennedy, D., Kluzek, E., Lawrence, P. J., Li, F., Li, H., Lombardozzi, D., Riley, W. J., Sacks, W. J., Shi, M., Vertenstein, M., Wieder, W. R., Xu, C., Ali, A. A., Badger, A. M., Bisht, G., van den Broeke, M., Brunke, M. A., Burns, S. P., Buzan, J., Clark, M., Craig, A., Dahlin, K., Drewniak, B., Fisher, J. B., Flanner, M., Fox, A. M., Gentine, P., Hoffman, F., Keppel-Aleks, G., Knox, R., Kumar, S., Lenaerts, J., Leung, L. R., Lipscomb, W. H., Lu, Y., Pandey, A., Pelletier, J. D., Perket, J., Randerson, J. T., Ricciuto, D. M., Sanderson, B. M., Slater, A., Subin, Z. M., Tang, J., Thomas, R. Q., Val Martin, M., and Zeng, X.: The Community Land Model Version 5: Description of New Features, Benchmarking, and Impact of Forcing Uncertainty, J. Adv. Model. Earth Syst., 11, 4245–4287,, 2019. 

Liang, X. and Guo, J.: Intercomparison of land-surface parameterization schemes: sensitivity of surface energy and water fluxes to model parameters, J. Hydrol., 279, 182–209,, 2003. 

Liang, X., Lettenmaier, D. P., Wood, E. F., and Burges, S. J.: A simple hydrologically based model of land surface water and energy fluxes for general circulation models, J. Geophys. Res., 99, 14415–14428,, 1994. 

Liang, X., Wood, E. F., and Lettenmaier, D. P.: Surface soil moisture parameterization of the VIC-2L model: Evaluation and modification, Global Planet. Change, 13, 195–206,, 1996. 

Liang, X., Wood, E. F., and Lettenmaier, D. P.: Modeling ground heat flux in land surface parameterization schemes, J. Geophys. Res.-Atmos., 104, 9581–9600,, 1999. 

Lilhare, R., Pokorny, S., Déry, S. J., Stadnyk, T. A., and Koenig, K. A.: Sensitivity analysis and uncertainty assessment in water budgets simulated by the variable infiltration capacity model for Canadian subarctic watersheds, Hydrol. Process., 34, 2057–2075,, 2020. 

Lohmann, D., Raschke, E., Nijssen, B., and Lettenmaier, D. P.: Hydrologie à l'échelle régionale: II. Application du modèle VIC-2L sur la rivière Weser, Allemagne, Hydrolog. Sci. J., 43, 143–158,, 1998. 

Marks, D. and Dozier, J.: Climate and Energy Exchange at the Snow Surface in the Alpine Region of the Sierra Nevada 2. Snow Cover Energy Balance, Water Resour. Res., 28, 3043–3054, 1992. 

Massoud, E. C., Xu, C., Fisher, R. A., Knox, R. G., Walker, A. P., Serbin, S. P., Christoffersen, B. O., Holm, J. A., Kueppers, L. M., Ricciuto, D. M., Wei, L., Johnson, D. J., Chambers, J. Q., Koven, C. D., McDowell, N. G., and Vrugt, J. A.: Identification of key parameters controlling demographically structured vegetation dynamics in a land surface model: CLM4.5(FATES), Geosci. Model Dev., 12, 4133–4164,, 2019. 

McCabe, M. F., Rodell, M., Alsdorf, D. E., Miralles, D. G., Uijlenhoet, R., Wagner, W., Lucieer, A., Houborg, R., Verhoest, N. E. C., Franz, T. E., Shi, J., Gao, H., and Wood, E. F.: The future of Earth observation in hydrology, Hydrol. Earth Syst. Sci., 21, 3879–3914,, 2017. 

Melsen, L., Teuling, A., Torfs, P., Zappa, M., Mizukami, N., Clark, M., and Uijlenhoet, R.: Representation of spatial and temporal variability in large-domain hydrological models: Case study for a mesoscale pre-Alpine basin, Hydrol. Earth Syst. Sci., 20, 2207–2226,, 2016. 

Melsen, L. A. and Guse, B.: Climate change impacts model parameter sensitivity-implications for calibration strategy and model diagnostic evaluation, Hydrol. Earth Syst. Sci., 25, 1307–1332,, 2021. 

Melsen, L. A., Teuling, A. J., Torfs, P. J. J. F., Zappa, M., Mizukami, N., Mendoza, P. A., Clark, M. P., and Uijlenhoet, R.: Subjective modeling decisions can significantly impact the simulation of flood and drought events, J. Hydrol., 568, 1093–1104,, 2019. 

Mendoza, P. A., Clark, M. P., Barlage, M., Rajagopalan, B., Samaniego, L., Abramowitz, G., and Gupta, H.: Are we unnecessarily constraining the agility of complex process-based models?, Water Resour. Res., 51, 716–728,, 2015a. 

Mendoza, P. A., Clark, M. P., Mizukami, N., Newman, A., Barlage, M., Gutmann, E., Rasmussen, R., Rajagopalan, B., Brekke, L., and Arnold, J.: Effects of hydrologic model choice and calibration on the portrayal of climate change impacts, J. Hydrometeorol., 16, 762–780,, 2015b. 

Misirli, F., Gupta, H. V, Sorooshian, S., and Thiemann, M.: Bayesian recursive estimation of parameter and output uncertainty for watershed models, in: vol. 6, edited by: Duan, Q., Gupta, H. V., Sorooshian, S., Rousseau, A. N., and Turcotte, R., American Geophysical Union, Washington, DC, 113–124,, 2003. 

Mitchell, K. E., Lohmann, D., Houser, P. R., Wood, E. F., Schaake, J. C., Robock, A., Cosgrove, B. A., Sheffield, J., Duan, Q., Luo, L., Higgins, W., Pinker, R. T., Tarpley, J. D., Lettenmaier, D. P., Marshall, C. H., Entin, J. K., Pan, M., Shi, W., Koren, V., Meng, J., Ramsay, B. H., and Bailey, A. A.: The multi-institution North American Land Data Assimilation System (NLDAS): Utilizing multiple GCIP products and partners in a continental distributed hydrological modeling system, J. Geophys. Res., 109, D07S90,, 2004. 

Mizukami, N., Clark, M., Slater, A., Brekke, L., Elsner, M., Arnold, J., Gangopadhyay, S., Clark, M. P., Slater, A. G., Brekke, L. D., Elsner, M. M., Arnold, J. R., Gangopadhyay, S., Clark, M., Slater, A., Brekke, L., Elsner, M., and Arnold, J.: Hydrologic Implications of Different Large-Scale Meteorological Model Forcing Datasets in Mountainous Regions, J. Hydrometeorol., 15, 474–488,, 2014. 

Mizukami, N., Clark, M. P., Gutmann, E. D., Mendoza, P. A., Newman, A. J., Nijssen, B., Livneh, B., Hay, L. E., Arnold, J. R., and Brekke, L. D.: Implications of the Methodological Choices for Hydrologic Portrayals of Climate Change over the Contiguous United States: Statistically Downscaled Forcing Data and Hydrologic Models, J. Hydrometeorol., 17, 73–98,, 2016. 

Montgomery, D. C.: Design and analysis of experiments, Wiley, ISBN 9781118146927, 1991. 

Morris, M. D.: Factorial sampling plans for preliminary computational experiments, Technometrics, 33, 161–174,, 1991. 

Myneni, R. B., Ramakrishna, R., Nemani, R., and Running, S. W.: Estimation of global leaf area index and absorbed par using radiative transfer models, IEEE T. Geosci. Remote, 35, 1380–1393,, 1997. 

Nash, J. and Sutcliffe, J.: River ?ow forecasting through conceptual models part I – A discussion of principles, J. Hydrol., 10, 282–290,, 1970. 

Niu, G.-Y., Yang, Z.-L., Mitchell, K. E., Chen, F., Ek, M. B., Barlage, M., Kumar, A., Manning, K., Niyogi, D., Rosero, E., Tewari, M., and Xia, Y.: The community Noah land surface model with multiparameterization options (Noah-MP): 1. Model description and evaluation with local-scale measurements, J. Geophys. Res., 116, D12109,, 2011. 

Pi, H. and Peterson, C.: Finding the Embedding Dimension and Variable Dependencies in Time Series, Neural Comput., 6, 509–520,, 1994. 

Pokhrel, P. and Gupta, H. V.: On the use of spatial regularization strategies to improve calibration of distributed watershed models, Water Resour. Res., 46, W01505,, 2010. 

Prihodko, L., Denning, A. S. S., Hanan, N. P. P., Baker, I., and Davis, K.: Sensitivity, uncertainty and time dependence of parameters in a complex land surface model, Agr. Forest Meteorol., 148, 268–287,, 2008. 

Rakovec, O., Hill, M. C., Clark, M. P., Weerts, A. H., Teuling, A. J., and Uijlenhoet, R.: Distributed Evaluation of Local Sensitivity Analysis (DELSA), with application to hydrologic models, Water Resour. Res., 50, 1–18,, 2014. 

Rawls, W. J., Ahuja, L. R., Brakensiek, D. L., and Shirmohammadi, A.: Infiltration and soil water movement, in: Handbook of hydrology, edited by: Maidment, D. R., McGraw-Hill Inc., New York, 5.1–5.51, ISBN 0070397325, 1992. 

Razavi, S. and Gupta, H. V.: What do we mean by sensitivity analysis? The need for comprehensive characterization of “global” sensitivity in Earth and Environmental systems models, Water Resour. Res., 51, 3070–3092,, 2015. 

Razavi, S. and Gupta, H. V.: A new framework for comprehensive, robust, and efficient global sensitivity analysis: 2. Application, Water Resour. Res., 52, 423–439,, 2016. 

Reba, M. L., Marks, D., Link, T. E., Pomeroy, J., and Winstral, A.: Sensitivity of model parameterizations for simulated latent heat flux at the snow surface for complex mountain sites, Hydrol. Process., 28, 868–881,, 2014. 

Reynolds, C. A., Jackson, T. J., and Rawls, W. J.: Estimating soil water-holding capacities by linking the Food and Agriculture Organization soil map of the world with global pedon databases and continuous pedotransfer functions, Water Resour. Res., 36, 3653–3662,, 2000. 

Rosero, E., Yang, Z.-L., Wagener, T., Gulden, L. E., Yatheendradas, S., and Niu, G.-Y.: Quantifying parameter sensitivity, interaction, and transferability in hydrologically enhanced versions of the Noah land surface model over transition zones during the warm season, J. Geophys. Res., 115, 1–21,, 2010. 

Rosolem, R., Gupta, H. V., Shuttleworth, W. J., Zeng, X., and de Gonçalves, L. G. G.: A fully multiple-criteria implementation of the Sobol' method for parameter sensitivity analysis, J. Geophys. Res., 117, 1–18,, 2012. 

Roudier, P., Andersson, J. C. M., Donnelly, C., Feyen, L., Greuell, W., and Ludwig, F.: Projections of future floods and hydrological droughts in Europe under a +2C global warming, Climatic Change, 135, 341–355,, 2016. 

Saltelli, A., Ratto, M., Andres, T., Campolongo, F., Cariboni, J., Gatelli, D., Saisana, M., and Tarantola, S.: Global Sensitivity Analysis The Primer,, 2008. 

Schmied, H. M., Eisner, S., Franz, D., Wattenbach, M., Portmann, F. T., Flörke, M., and Döll, P.: Sensitivity of simulated global-scale freshwater fluxes and storages to input data, hydrological model structure, human water use and calibration, Hydrol. Earth Syst. Sci., 18, 3511–3538,, 2014. 

Sheikholeslami, R., Gharari, S., Papalexiou, S. M., and Clark, M. P.: VISCOUS: A Variance-Based Sensitivity Analysis Using Copulas for Efficient Identification of Dominant Hydrological Processes, Water Resour. Res., 57, 1–24,, 2021. 

Shi, X., Wood, A. W., and Lettenmaier, D. P.: How essentialis hydrologic model calibration to seasonal stream flow forecasting?, J. Hydrometeorol., 9, 1350–1363,, 2008. 

Shukla, S., Sheffield, J., Wood, E. F., and Lettenmaier, D. P.: On the sources of global land surface hydrologic predictability, Hydrol. Earth Syst. Sci., 17, 2781–2796,, 2013. 

Sobol', I. M.: Global sensitivity indices for nonlinear mathematical models and their Monte Carlo estimates, Math. Comput. Simul., 55, 271–280, 2001. 

Sobol', I. M. and Kucherenko, S.: A new derivative based importance criterion for groups of variables and its link with the global sensitivity indices, Comput. Phys. Commun., 181, 1212–1217,, 2010. 

Tang, G., Clark, M. P., Papalexiou, S. M., Newman, A. J., Wood, A. W., Brunet, D., and Whitfield, P. H.: Emdna: An ensemble meteorological dataset for north america, Earth Syst. Sci. Data, 13, 3337–3362,, 2021. 

Tian, S., Tregoning, P., Renzullo, L. J., van Dijk, A. I. J. M., Walker, J. P., Pauwels, V. R. N., and Allgeyer, S.: Improved water balance component estimates through joint assimilation of GRACE water storage and SMOS soil moisture retrievals, Water Resour. Res., 53, 1820–1840,, 2017. 

Tian, Y., Woodcock, C. E., Wang, Y., Privette, J. L., Shabanov, N. V., Zhou, L., Zhang, Y., Buermann, W., Dong, J., Veikkanen, B., Häme, T., Andersson, K., Ozdogan, M., Knyazikhin, Y., and Myneni, R. B.: Multiscale analysis and validation of the MODIS LAI product I. Uncertainty assessment, Remote Sens. Environ., 83, 414–430,, 2002. 

Todini, E.: The ARNO rainfall-runoff model, J. Hydrol., 175, 339–382, 1996. 

Tonkin, M. J. and Doherty, J.: A hybrid regularized inversion methodology for highly parameterized environmental models, Water Resour. Res., 41, 1–16,, 2005. 

UNEP: World atlas of desertification, ed. 2, edited by: Middleton, N. and Thomas, D., Arnold, Hodder Headline, PLC, (last access: 13 June 2022), 1997. 

USACE: Snow hydrology: Summary report of the snow investigations, North Pacific Division, Corps of Engineers, US Army, (last access: 13 June 2022), 1956. 

USGS: EarthExplorer, USGS [data set],, last access: 13 June 2022. 

Vano, J. A. and Lettenmaier, D. P.: A sensitivity-based approach to evaluating future changes in Colorado River discharge, Climatic Change, 122, 621–634,, 2014. 

Vásquez, N., Cepeda, J., Gómez, T., Mendoza, P. A., Lagos, M., Boisier, J. P., Álvarez-Garretón, C., and Vargas, X.: Catchment-Scale Natural Water Balance in Chile, in: Water Resources of Chile, Springer, 189–208,, 2021. 

Verbist, K., Santibañez, F., Gabriels, D., and Soto, G.: Documento Técnico No. 25, Atlas de Zonas Áridas de América Latina y el Caribe, (last access: 13 June 2022), 2010. 

Vicuña, S., Vargas, X., Boisier, J. P., Mendoza, P. A., Gómez, T., Vásquez, N., and Cepeda, J.: Impacts of Climate Change on Water Resources in Chile, in: Water Resources of Chile, vol. 13, Springer, 347–363,, 2021. 

Vořechovský, M.: Hierarchical Refinement of Latin Hypercube Samples, Comput. Civ. Infrastruct. Eng., 30, 394–411,, 2015. 

Wi, S., Ray, P., Demaria, E. M. C., Steinschneider, S., and Brown, C.: A user-friendly software package for VIC hydrologic model development, Environ. Model. Softw., 98, 35–53,, 2017. 

Wood, A. W., Kumar, A., and Lettenmaier, D. P.: A retrospective assessment of National Centers for Environmental prediction climate model-based ensemble hydrologic forecasting in the western United States, J. Geophys. Res.-Atmos., 110, 1–16,, 2005. 

Wood, E. F., Lettenmaier, D. P.., and Zartarian, V. G.: A Land-Surface Hydrology Parameterization With Subgrid Variability for General Circulation Models, J. Geophys. Res., 97, 2717–2728,, 1992. 

Woodward, J. L.: Estimating the Flammable Mass of a Vapor Cloud, John Wiley & Sons, Inc., Hoboken, NJ, USA,, 1999. 

Xia, Y., Mitchell, K., Ek, M., Sheffield, J., Cosgrove, B., Wood, E., Luo, L., Alonge, C., Wei, H., Meng, J., Livneh, B., Lettenmaier, D., Koren, V., Duan, Q., Mo, K., Fan, Y., and Mocko, D.: Continental-scale water and energy flux analysis and validation for the North American Land Data Assimilation System project phase 2 (NLDAS-2): 1. Intercomparison and application of model products, J. Geophys. Res., 117, D03109,, 2012.  

Xie, Z., Yuan, F., Duan, Q., Zheng, J., Liang, M., and Chen, F.: Regional parameter estimation of the VIC land surface model: Methodology and application to river basins in China, J. Hydrometeorol., 8, 447–468,, 2007. 

Yadav, M., Wagener, T., and Gupta, H.: Regionalization of constraints on expected watershed response behavior for improved predictions in ungauged basins, Adv. Water Resour., 30, 1756–1774,, 2007. 

Yang, Y., Pan, M., Beck, H. E., Fisher, C. K., Beighley, R. E., Kao, S. C., Hong, Y., and Wood, E. F.: In Quest of Calibration Density and Consistency in Hydrologic Modeling: Distributed Parameter Calibration against Streamflow Characteristics, Water Resour. Res., 55, 7784–7803,, 2019. 

Yang, Y., Pan, M., Lin, P., Beck, H. E., Zeng, Z., Yamazaki, D., David, C. H., Lu, H., Yang, K., Hong, Y., and Wood, E. F.: Global reach-level 3-hourly river flood reanalysis (1980–2019), B. Am. Meteorol. Soc., 102, E2086–E2105,, 2021. 

Yang, Z.-L., Niu, G.-Y., Mitchell, K. E., Chen, F., Ek, M. B., Barlage, M., Longuevergne, L., Manning, K., Niyogi, D., Tewari, M., and Xia, Y.: The community Noah land surface model with multiparameterization options (Noah-MP): 2. Evaluation over global river basins, J. Geophys. Res., 116, 1–16,, 2011. 

Yeste, P., García-Valdecasas Ojeda, M., Gámiz-Fortis, S. R., Castro-Díez, Y., and Esteban-Parra, M. J.: Integrated sensitivity analysis of a macroscale hydrologic model in the north of the Iberian Peninsula, J. Hydrol., 590, 125230,, 2020. 

Zegers, G., Mendoza, P. A., Garces, A., and Montserrat, S.: Sensitivity and identifiability of rheological parameters in debris flow modeling, Nat. Hazards Earth Syst. Sci., 20, 1919–1930,, 2020. 

Zhao, Q., Ye, B., Ding, Y., Zhang, S., Yi, S., Wang, J., Shangguan, D., Zhao, C., and Han, H.: Coupling a glacier melt model to the Variable Infiltration Capacity (VIC) model for hydrological modeling in north-western China, Environ. Earth Sci., 68, 87–101,, 2013. 

Zhao, R.-J., Zuang, Y.-L., Fang, L.-R., Liu, X.-R., and Zhang, Q.-S.: Xinanjiang Model, IAHS-AISH Publ., 129, 351–356, 1980. 

Zink, M., Kumar, R., Cuntz, M., and Samaniego, L.: A high-resolution dataset of water fluxes and states for Germany accounting for parametric uncertainty, Hydrol. Earth Syst. Sci., 21, 1769–1790,, 2017. 

Short summary
This paper characterizes parameter sensitivities across more than 5500 grid cells for a commonly used macroscale hydrological model, including a suite of eight performance metrics and 43 soil, vegetation and snow parameters. The results show that the model is highly overparameterized and, more importantly, help to provide guidance on the most relevant parameters for specific target processes across diverse climatic types.