Development, sensitivity and uncertainty analysis of LASH model

Desenvolvimento, sensibilidade e análise de incertezas do modelo LASH

Abstracts

Many hydrologic models have been developed to help manage natural resources all over the world. Nevertheless, most models have presented a high complexity regarding data base requirements, as well as, many calibration parameters. This has brought serious difficulties for applying them in watersheds where there is scarcity of data. The development of the Lavras Simulation of Hydrology (LASH) in a GIS framework is described in this study, which focuses on its main components, parameters, and capabilities. Coupled with LASH, sensitivity analysis, parameter range reduction, and uncertainty analysis were performed prior to the calibration effort by using specific techniques (Morris method, Monte Carlo simulation and a Generalized Likelihood Uncertainty Estimation -GLUE) with a data base from a Brazilian Tropical Experimental Watershed (32 km²), in order to predict streamflow on a daily basis. LASH is a simple deterministic and spatially distributed model using long-term data sets, and a few maps to predict streamflow at a watershed outlet. We were able to identify the most sensitive parameters which are associated with the base flow and surface runoff components, using a reference watershed. Using a conservative threshold, two parameters had their range of values reduced, thus resulting in outputs closer to measured values and facilitating automatic calibration of the model with less required iterations. GLUE was found to be an efficient method to analyze uncertainties related to the prediction of mean daily streamflow in the watershed.

GLUE methodology; Morris method; hydrologic modeling; parameter ranges; automatic calibration

AGRICULTURAL ENGINEERING

Development, sensitivity and uncertainty analysis of LASH model

Desenvolvimento, sensibilidade e análise de incertezas do modelo LASH

Samuel Beskow^I,II; Carlos Rogério de MelloII,^*; Lloyd Darrell Norton^III

^IUFPel – Centro de Desenvolvimento Tecnológico/Engenharia Hídrica – 96060-290 – Pelotas, RS – Brasil

^IIUFLA – Depto. de Engenharia, C.P. 3037 – 37200-000 – Lavras, MG – Brasil

^IIIUSDA/ARS – National Soil Erosion Research Laboratory, 275 South Russell Street, Purdue University, 47907-2077 – West Lafayette, IN – USA.

ABSTRACT

Many hydrologic models have been developed to help manage natural resources all over the world. Nevertheless, most models have presented a high complexity regarding data base requirements, as well as, many calibration parameters. This has brought serious difficulties for applying them in watersheds where there is scarcity of data. The development of the Lavras Simulation of Hydrology (LASH) in a GIS framework is described in this study, which focuses on its main components, parameters, and capabilities. Coupled with LASH, sensitivity analysis, parameter range reduction, and uncertainty analysis were performed prior to the calibration effort by using specific techniques (Morris method, Monte Carlo simulation and a Generalized Likelihood Uncertainty Estimation -GLUE) with a data base from a Brazilian Tropical Experimental Watershed (32 km²), in order to predict streamflow on a daily basis. LASH is a simple deterministic and spatially distributed model using long-term data sets, and a few maps to predict streamflow at a watershed outlet. We were able to identify the most sensitive parameters which are associated with the base flow and surface runoff components, using a reference watershed. Using a conservative threshold, two parameters had their range of values reduced, thus resulting in outputs closer to measured values and facilitating automatic calibration of the model with less required iterations. GLUE was found to be an efficient method to analyze uncertainties related to the prediction of mean daily streamflow in the watershed.

Key words: GLUE methodology, Morris method, hydrologic modeling, parameter ranges, automatic calibration

RESUMO

Diversos modelos hidrológicos têm sido desenvolvidos no intuito de auxiliar na gestão de recursos naturais em todo o mundo. Porém, a maioria desses modelos apresenta um alto grau de complexidade em relação tanto à necessidade de base de dados, quanto ao número de parâmetros de calibração. Em virtude desses fatores, se torna difícil a aplicação em bacias hidrográficas que têm bases de dados reduzidas. Neste artigo é descrito o desenvolvimento do modelo Lavras Simulation of Hydrology (LASH) em uma estrutura de SIG, buscando enfatizar seus principais componentes e parâmetros, bem como suas potencialidades. Além da descrição do modelo, também foram realizadas a análise de sensibilidade, a redução do intervalo de parâmetros e a análise de incertezas, anteriormente à fase de calibração, utilizando metodologias específicas (método de Morris, simulação de Monte Carlo e o método Generalized Likelihood Uncertainty Equation (GLUE)), com a base de dados de uma bacia hidrográfica experimental tropical brasileira (32 km²), a fim de simular a vazão total média diária. O LASH é um modelo classificado como determinístico e distribuído, que utiliza dados de longo termo e poucos mapas para predizer vazão na seção de controle de bacias hidrográficas. Foi possível identificar os parâmetros mais sensíveis do modelo para a bacia hidrográfica de referência, os quais estão associados com os componentes de escoamento de base e superficial direto. Em função do limiar conservador utilizado neste estudo, foram reduzidos os intervalos de dois parâmetros, dessa forma gerando resultados simulados mais realísticos e também facilitando a calibração automática do modelo com um menor número de iterações necessárias. O método da GLUE mostrou ser eficiente frente à análise de incertezas relacionadas à predição de vazão na bacia de estudo.

Palavras-chave: metodologia GLUE, método de Morris, modelagem hidrológica, intervalos de parâmetros, calibração automática

Introduction

Studies concerning models applied to watersheds usually consider the following steps: calibration, validation, Most hydrologic models are too complex to be used in and prediction (Arabi et al., 2007). Calibration is a fundaareas with limited data (Beskow et al., 2009b). Under this mental step towards application of hydrologic models, be-aspect, a model with a simple approach which makes use of cause it enables one to fit a set of parameters which are less data may be preferable for better water resource man-unique for a specific watershed. It is very important to deagement. The LASH model is a simple hydrologic model termine which parameters may cause the most significant efdesigned for prediction of streamflow in watersheds where fect on the output of interest through sensitivity analysis there are limited physical data. (Benaman and Shoemaker, 2004). Sensitivity analysis evalu-

ates which parameters should be taken into account in the calibration phase (Blasone et al., 2008). This kind of analysis has been successfully applied in several studies concerning hydrologic modeling prior to calibration (Muleta and Nicklow, 2005; White and Chaubey, 2005; Van Griensven et al., 2006).

A model may contain parameters not based on field data (Benaman and Shoemaker, 2004). Therefore, to these parameters are given suggested ranges of values found for other sites. Nevertheless, some parameters may have very broad ranges instead of narrow and acceptable ranges, thus resulting in not only many unrealistic estimations, but also inefficiency in the optimization method (Blasone et al., 2008). Methodologies based on sensitivity and uncertainty analyses including Monte Carlo simulations have been widely used in studies concerning reduction of parameter ranges which both have broad suggested intervals and generate too many inaccurate predictions (Arabi et al., 2007; Wei et al., 2008; Huang and Liang, 2006).

The objectives of this paper were to: (i) describe the LASH model as well as its parameters and capabilities; (ii) identify the most sensitive parameters of the model and to reduce their suggested ranges; and (iii) evaluate uncertainties with respect to streamflow predictions prior to model calibration.

Material and Methods

A new hydrological computational model was developed for this study, named as Lavras Simulation of Hydrology (LASH). This model is similar to the one presented by Mello et al. (2008), however, this new model uses a distributed approach instead of either lumped or semi-distributed. LASH was designed to take into account temporal and spatial variability of all the variables included in the hydrologic components by dividing the watershed into homogeneous grid cells.

This is a semi-physically based and continuous simulation model, using the following components on a daily basis: evapotranspiration, interceptation of precipitation by land cover, capillary rise, soil water availability, surface runoff, subsurface flow, and base flow. The model is divided into three basic modules: (i) its first module is designed to compute surface runoff flow, subsurface flow, base flow, and capillary rise, which are drained from the soil layer considered in the water balance; (ii) the second module generates flow within each cell to the stream network; this module takes into account the lag effect using the concept of linear soil reservoir (Collischonn et al., 2007; Tucci, 2005); (iii) in the last module, LASH employs the Muskingum-Cunge Linear Model to propagate the flows through the channel network.

Soil water balance is calculated at each time step for each grid cell in the watershed (Equation 1). The number of grid cells depends on both the cell size and how large the study watershed is.

where j and i are indexes associated to time step and grid cell, respectively; A_t^j is the soil water availability (mm) for the grid cell i at the end of the time step j; represents the soil water availability (mm) for the grid cell i at the start of the time step j; Δt is the time step (daily); P_i corresponds to the precipitation (mm Δt^–1) minus the interceptation of precipitation by land cover; ET_i is the evapotranspiration (mm Δt^–1); Ds_i is the surface runoff (mm); Dss_i represents the subsurface flow (mm Δt^–1); DB_i is the base flow (mm Δt^–1); and DCR_i corresponds to the capillary rise depth (mm Δt^–1). The variable is computed for each time step for each cell.

Once precipitation begins, it is stored on the vegetation cover until maximum interceptation storage (I_max) is reached, which is calculated (Equation 2) for each grid cell as a linear function of Leaf Area Index (LAI) (Almeida et al., 2007; Collischonn et al., 2007; Zhou et al., 2006). LASH uses the Penman-Monteith equation (Allen et al., 1998) to compute how much from intercepted water is evaporated in each time step.

Imax_i = α. LAI_i (2)

where α is the interceptation coefficient, assumed to be 0.2 mm (Collischonn et al., 2007); LAI is the leaf area index (m² m^–2), which can be obtained from either literature or field trials. Since LAI values may have a great variation over time, an option was implemented in LASH so that users can input a separate file (linked to the land-use map) in order to represent such time dependent variation for each land-use.

The Curve-Number Modified Mishra-Singh (CN-MMS) model (Mishra et al., 2003) was employed in LASH for estimating the surface runoff component (D_S), in mm (Equation 3). All of the variables with a subscript i in equations mean that such variable is computed for each grid cell at time j.

where I _a is the initial abstraction (mm); M corresponds to the antecedent soil moisture (mm); and S is the soil potential maximum retention (mm).

where γ is an initial abstraction coefficient (dimensionless); and P₅ represents the 5-day antecedent precipitation (mm). Even though γ can vary between 0 and ∞, Mishra et al. (2003) and Mishra et al. (2006) suggest employing values from 0 to 0.5 to calibrate this parameter.

where θ_S is the saturation soil water content (m³ m^–3); θ₀ represents the current soil water content (m³ m^–3); and h is the where j and i are indexes associated to time step and grid control layer of the water budget (rooting depth, mm).

A_mi = (θ_Si – θ_PWPi) . h_i (7)

where A_m is the maximum soil water availability (mm); and θPWP corresponds to the permanent wilting point soil water content (m³ m^–3).

S_i = Am_i – At_i (8)

where A_t represents the soil water availability at time j.

The Brooks and Corey equation (Rawls et al., 1993) was adapted in the LASH to simulate the subsurface (D_SS) e base flow (D_B) components (mm per day).

where K_SS corresponds to the hydraulic conductivity of the subsurface reservoir (mm per day), which is a calibration parameter; A_CC is the minimum soil water availability to generate subsurface flow (mm); PS represents the pore-size index which may be assumed constant and equal to 0.4 due to its low sensitivity (Collischonn et al., 2007).

The following equation was implemented in LASH in order to simulate the base flow component (D_B), in mm per day (Collischonn et al., 2007):

where K_B corresponds to the hydraulic conductivity of the shallow saturated zone reservoir (mm per day), which has been considered as a parameter of calibration due to difficulty in obtaining it through field trials; and A_C is the minimum soil water availability to generate base flow (mm).

The capillary rise component (D_CR), in mm per day, was implemented in LASH to allow simulation of situations in which some areas of the watershed have low soil water availability, thus occurring rise of water from shallow aquifer to the soil layer and making it available for evapotranspiration.

where K_CR corresponds to the maximum flow returning to soil by capillary rise (mm per day), being considered a parameter of calibration; and A_CR is the soil water availability limit (mm) so that capillary rise occurs.

The evapotranspiration (ET) module, in mm per day, was incorporated in LASH making use of the Penman-Monteith equation, described in Allen et al. (1998). The variables for ET calculation can be input into LASH taking into account their spatial and temporal variability, since some of them may have a considerable variation throughout the year. However, if soil moisture is less than a given limit of soil water availability, actual evapotranspiration is less than or equal to the crop evapotranspiration. The relationship between crop evapotranspiration and actual evapotranspiration can be expressed by the coefficient K_S which depends on the soil moisture in time j (Allen et al., 1998).

Once the model computes D_S, D_SS and D_B, it converts each flow component to discharge. It is necessary to account for the delay of the inflow to the stream network, therefore, a flow routing method has to be used. The method of linear reservoirs was chosen to route flow through each cell, since it is a simple approach and it has been successfully used in many other studies (Zhou et al., 2006; Collischonn et al., 2007; Mello et al., 2008). The following equations are used to compute outflow from the three reservoirs as previously stated:

where, QS_i, QSS_i and QB_i are the outflows from surface, subsurface and groundwater reservoirs (m³ s^–1) of the cell i, respectively; , , and are the water volumes in the surface, subsurface and groundwater reservoirs (m³) of the cell i at time j, respectively; TC_i corresponds to the time of concentration (s); C_B represents the recession time (s) which can be calculated from a hydrograph previously monitored in the watershed; and C_S and C_SS are response time parameters. It is worthwhile to point out that (C_S . T_C) is less than (C_SS . T_C) and the latter is less than C_B due to different delays in each reservoir. The time of concentration can be estimated by different equations which are described in literature. Mean daily streamflow at the outlet is obtained by summing the three outflow components as stated before, which are propagated through the drainage network using the Muskingum-Cunge Linear Model.

There are many variables which may have spatial variation. All of the maps used by the model are derived from the DEM, land-use, soil or channel network maps. In addition to the maps, LASH also needs two other files in table format. The first one contains information on climatic data as well as observed discharge (m³ s^–1) over the time. The following climatic data are necessary in the model to compute evapotranspiration on a daily basis according with Penman-Monteith equation.

A second file in table format is used to inform the model of the variation in parameters connected to the land-use over time, for instance, leaf area index (m² m^–2), height (m), albedo (dimensionless), surface resistance (s m^–1), rooting depth (mm), and crop coefficient (dimensionless). The optimization algorithm is based on the Shuffled Complex Evolution (SCE-UA) method (Duan et al., 1992). The SCEUA is a global optimization method that has been extensively used and accepted in the field of hydrology for several years. Complete details about this method can be found in Duan et al. (1992) and Duan et al. (1994).

For this study, all the procedures are based on the methodology described by Benaman and Shoemaker (2004) fol-lowing four steps, however, some adaptations were made in accordance with Arabi et al. (2007). Basically, the methodology uses two sampling methods: (i) the Morris method (Morris, 1991) which can be classified as One-factor-At-a-Time (OAT) sensitivity analysis (Saltelli et al., 2004); and (ii) Generalized Likelihood Uncertainty Estimation – GLUE (Beven and Binley, 1992). A twelve month period was used for all the analyses conducted in this study.

The first step of this method regards the initial sensitivity analysis and establishment of upper and lower bounds for each parameter. It assesses all of the calibration parameters used in the model to choose those that are most sensitive and those with large uncertainties. The Morris method (Morris, 1991) was applied to carry out the sensitivity analysis of LASH. Local sensitivity indices (d_k) were computed for each parameter (Table 1) applying Equation 15, while global sensitivity indices (d_g) were obtained by taking the average of these local sensitivities at different points sampled in the respective parameter space.

where y (x) represents the model output of interest; Δ is the elementary effect of a small perturbation of the k^th component. The perturbation is a predetermined multiple of 1/ (p - 1) in which p corresponds to the number of intervals that a parameter range is divided by.

In this study, a rescaled sensitivity index (d_r) was applied, which was recommended by Arabi et al. (2007). The d_r index was determined by dividing the global sensitivity indices by their total sum. The values of d_r can vary between 0 and 1 in such a way that the greater the value, the more sensitive the respective parameter is.

The second step regards the Initial Monte Carlo simulation. It starts once the first step is carried out, it is possible to define all of the parameters that will be used in step 2. For this step, the Monte Carlo method was employed to evaluate the results from the model by varying all the uncertain parameters simultaneously (step 1), taking values within their recommended range (Table 1). Values for the sensitive parameters were selected randomly according to their respective probability distribution. As there is no prior information with respect to the probability distribution of the parameters, uniform distributions were assumed for all the parameters by taking into account minimum and maximum values for each parameter (Table 1). Studies carried out by Beven and Freer (2001), Benaman and Shoemaker (2004), Arabi et al. (2007), and Wei et al. (2008) also used the same assumption.

Once all the runs are performed for step 2, a cumulative density function can be generated (Benaman and Shoemaker, 2004). This function describes the variation in output as a result of all the uncertain parameters, which are varied at the same time within their suggested ranges. In this case, the analyzed output variable was mean daily streamflow.

Finally, it is necessary to verify if the analyzed output variable is close to the median of the results obtained through Monte Carlo simulation as well as to determine if the probability distribution has a reasonable range with respect to the model output variable. If so, it can be assumed that the parameter ranges suggested in Table 1 are realistic for this watershed. Otherwise, some parameters may have their ranges unrealistic for this watershed; therefore, another analysis may have to be performed following steps 3 – 4 described below.

The third step regards the Range adjustment using in-terval-spaced sensitivity. This step is necessary when some parameter bounds need to be narrowed, since the combination of different parameters provides a considerable number of unrealistic results.

Interval-spaced sensitivity was performed to assess the effect of different model parameters on the mean daily streamflow values. For this analysis, each parameter of concern is evaluated while keeping the remaining parameters at their base value. Base values can be selected by different approaches: (i) previous manual or automatic calibration; (ii) values obtained from either literature or field measurements carried out in the same watershed. Each analyzed parameter had its range divided into equal intervals. The choice relative to the number of intervals is subjective; however, Arabi et al. (2007) recommended using 20 - 50 intervals for each parameter range. Following this mathematical process a graph can be produced, which allows evaluating if the parameter range is generating unrealistic results. Since model parameters have different ranges, the following equation was applied to compute normalized values for the x-axis for each parameter (Arabi et al., 2007).

where N_i is the normalized value of the parameter i, which is determined as a function of absolute value α_i and its upper (UB_i) and lower (LB_i) bounds (Table 1).

Nash-Sutcliffe (C_NS) and its logarithmic version (log(C_NS)) (Nash and Sutcliffe, 1970) efficiency coefficients were used in this study as cutoff criteria in order to distinguish between behavior and non-behavior range.

where QO_i corresponds to the observed streamflow at time j, QO_M is the mean observed streamflow, QS_j represents the simulated streamflow at time j, QS_M is the mean simulated streamflow, log(QO_j) is the logarithm of the observed streamflow at time j, log(QS_j) is the logarithm of the simulated streamflow at time j, log_m(Q_O) corresponds to the mean logarithm of the observed streamflows. Once the efficiency coefficients are computed for each parameter, ranges (upper and lower bound) can be adjusted.

The fourth step – Final Monte Carlo Simulation – involves Monte Carlo analysis in order to assess the influence of the narrowed parameter range over the model outputs. This part of the analysis is similar to Step 2; however, a GLUE likelihood measure is computed for each model run using equation 19 according with procedures described in Arabi et al. (2007) and Blasone et al. (2008). In addition, new parameter values were taken simultaneously for all the parameters that had a considerable sensitivity.

where L(p|Y) corresponds to the likelihood measure of the parameter set (p) for the observed data (Y); the values σ_S and σ_O represent the variance of the error between model prediction and observed data, and the variance of the observed data, respectively.

Parameter data sets that produce likelihood values less than a certain threshold are considered "nonbehavioral" and then discarded, whereas, the rest of the parameter sets are termed "behavioral" (Stedinger et al., 2008). The latter parameter sets are then used to compute likelihood weights, thereafter, each weight is divided by their sum and then they are sorted. In this way, a cumulative distribution for the model output parameter of interest is created and used for estimating uncertainty bounds.

The Jaguara Experimental Watershed (JEW) was used as a case study area to apply the LASH model and the above-mentioned procedures related to sensitivity analysis, adjustment of parameter ranges, and uncertainty analysis. This experimental watershed is located in southern Minas Gerais State, Brazil (Figure 1), has an area of about 32 km² and its drainage network presents a permanent hydrological regime.

The annual mean temperature in this region is approximately 19ºC, varying from 14 to 22ºC. According with Köppen's classification the climate is characterized as Cwa, with high concentration of precipitation during both spring and summer (from October to March), whereas, autumn and winter are dry and cool. The mean annual precipitation is about 1400 mm.

The used period was a consistent year with no gaps or problems derived from monitoring and, in addition, it contained a great variation on hydrological events including high peak flow as well as many consecutive drought days.

Variables with respect to weather were monitored every 30 min in the JEW through a complete and compact weather station for the whole period. Such variables included temperature, relative humidity, wind speed, and solar radiation. The discharge data set was obtained from an automatic gauge station located at the JEW's outlet and a stage-discharge rating curve.

An image from the satellite ALOS from May, 2008 which provides multispectral and 10-m resolution images, was acquired in order to classify the land-uses in JEW. This allowed us to determine that the watershed was occupied by eucalyptus (7.68%), coffee (4.11%), bare soil (9.63%), maize (22.85%), native vegetation (13.13%), and pasture (42.60%). The digital elevation model (DEM) of JEW was obtained with 30-m resolution, which enabled the generation of a slope map and delineation of the channel network. The slope gradients of this watershed ranged from 0 to 46.9%, with a mean value of 11.6%. The soil map of JEW presented the following percentage distribution: Oxisols (59.8%), Cambisols (23.4%) and Fluvic Neosols (16.8%) (Menezes et al., 2009).

Field trials were carried out in the JEW in order to quantify the spatial variation of saturated soil water content and permanent wilting point water content (Figure 1). Several points of measurement were sampled throughout the JEW and used to apply geostatistical procedures; thereafter, a semivariogram and a kriging map were obtained for each variable.

Results and discussion

The Morris method (Morris, 1991) was applied to determine which model parameters should be considered as most sensitive. The greater the rescaled sensitivity index (d_r), the more sensitive the respective parameter is. In this way, the LASH model parameters are outlined in Table 2 in a descending order with respect to sensitivity in mean daily streamflow.

Even though the methodology suggested by Morris (1991) is unsophisticated, it allowed us to identify the most sensitive parameters for the LASH model.

All parameters in the sensitivity analysis were important except n and Q_R which are part of the Muskingum-Cunge routing method (Table 2). This indicates that the model is not sensitive to the Muskingum-Cunge routing method for this watershed. A reasonable explanation for the non-sensitivity of these two parameters is that the channels in JEW are narrow and steep, thus the channel network does not have a significant influence on water accumulation and propagation (Collischonn et al., 2007).

The parameters K_B and A_C were the most sensitive in the model, proving that the JEW is governed by base flow process. In other studies (Collischon et al., 2007; Mello et al., 2008), A_C was considered to have a low sensitivity and, therefore, it was kept as a constant value. However, in this watershed our study showed that this parameter was greatly sensitive and should be taken into account during the calibration step, especially in watersheds that are governed by base flow. Since A_C is computed as function of A_m and the latter variable is calculated in a spatially distributed approach throughout the watershed (based on different values of saturation soil moisture and wilting point soil moisture), A_C might have been estimated better in this study than in the above-mentioned studies. As JEW is composed mostly by deep soils with slope gradients less than 18% (Oxissols) and the mean annual precipitation is greater than 1,400 mm, the aquifer recharge process is highly significant. By analyzing the contribution of each flow component for the studied period, we found results as follows: base flow (62.51%), surface runoff (25.10%), and subsurface flow (12.39%). Therefore, base flow is predominant in comparison to direct surface runoff and subsurface flow, thus justifying why the parameter (K_B) associated with the former component was the most sensitive.

The parameter related to surface runoff (γ) also had a high sensitivity. This behavior goes along with results presented in other papers (Mello et al., 2008; Beskow et al., 2009a). Variation in values for this parameter is expected to occur when simulating different watersheds, since this parameter is strongly affected by the following factors (Mishra et al., 2003; Beskow et al., 2009a): (a) climatic conditions such as temperature and solar radiation, which are important for evapotranspiration calculation; and (b) pluvial behavior of the study watershed due to its influence on antecedent soil moisture content.

The soil water content (θ₀) represents the hydrological initial condition of watersheds, and was found to have considerable sensitivity. The establishment of realistic ranges for this parameter depends on when the simulation starts. As in this case, the simulation began in January, which is the most humid month, soil moisture was high. Viola et al. (2009) recommend setting constant values of θ₀ for long-term hydrologic simulation models due to the difficulty of soil moisture monitoring at the watershed scale. However, we suggest including θ₀ as a calibration variable for improvement of the model performance in future simulation using LASH if this input is unavailable.

The parameter A_CR presented a considerable sensitivity for this watershed (Table 2). It is advisable to take into account capillary rise, since this may be an important component if the watershed has savannahs, riparian forests or springs, thus causing rise of water from shallow aquifer to the soil layer. Two statistics of precision are presented in Table 2 (C_NS and log (C_NS)), which are goodness-of-fit measures and refer to the minimum values found in the interval-spaced sensitivity analysis. These statistics can be used as a cutoff criterion in this study to eliminate portions of the parameter ranges in which either C_NS or log (C_NS) is negative.

There were only two parameters (K_B and A_C) that had minimum values of log (C_NS) below zero, which indicates that their parameter ranges are producing unrealistic outputs and then they need to be narrowed.

After performing the prior step (sensitivity analysis), we chose all of the parameters presented in Table 2, except Q_R and n, to carry out this analysis for mean daily streamflow as the output variable studied in the Monte Carlo based simulations. The choice for the number of Monte Carlo simulations is quite subjective; nevertheless, this step simply attempts to find out if the parameter ranges should be narrowed. This way, it is not necessary to run as many iterations as for the uncertainty analysis. For the same kind of analysis, Benaman and Shoemaker (2004) used 500 Monte Carlo based runs by varying simultaneously 36 input parameters for the SWAT model in order to simulate streamflow and sediment transport at the watershed scale. In the present study, only 13 parameters of the LASH model were evaluated, and 500 iterations were found to be a reasonable number.

A 12-month simulation period was used to produce Figure 2. It represents the cumulative probability distribution of daily streamflow using 500 Monte Carlo runs which involved variation of the model input parameters simultaneously (using ranges established in Table 1). The mean observed streamflow is presented in Figure 2 through a solid vertical line. A dashed line was drawn in order to emphasize the streamflow corresponding to the median of the cumulative probability distribution. Comparing these two lines one can observe that the line representing observed streamflow is distant from the median of the cumulative probability distribution. Therefore, we concluded that the streamflow results are biased to low. Most of the runs (about 95%) provided mean daily streamflow less than the observed streamflow. This occurred because of uncertainties concerning some parameters, whose initial ranges should be narrowed to facilitate a better automatic calibration. Otherwise, the automatic calibration process would take many parameter sets and would produce unrealistic results.

Thirteen parameters from Table 2 were used to assess interval-spaced sensitivity based on 50 runs for each parameter within its respective range (Table 1). A manual calibration was applied to determine the base values for the model parameters used in this study. The most difficult decision in this step was to define a cutoff criterion (threshold, T), since it is subjective. We used an assumption recommended by Arabi et al. (2007), in which statistics C_NS and log (C_NS) < 0 were considered unacceptable and used as a threshold to reduce parameter ranges. Based on the results presented in Table 2, through the columns Min (C_NS) and Min (log (C_NS)), it was possible to determine which parameter bounds were reduced. Figure 3 shows the four most sensitive parameters and how their ranges were narrowed using the threshold T. In addition to the threshold used (T < 0), Table 3 also illustrates the parameter ranges that would have been narrowed if other T values had been applied.

Two parameters presented a large range that had a major influence on the model output. These would cause too many parameter sets to be unrealistic for this watershed in case of using suggested ranges for calibration, and in addition, the optimization method would be inefficient. Because only one parameter is changed at a time, we should be careful with the choice of the threshold for reduction of ranges. According to Too high values should not be set as threshold because realistic portions of parameter ranges may be removed if all parameters are considered at once (Benaman and Shoemaker, 2004). Arabi et al. (2007) and Benaman and Shoemaker (2004) used similar methodologies for adjustment of parameter ranges of the SWAT model. In both studies, results allowed them to conclude that the streamflow-related parameters did not need a range reduction. As the sediment-related parameters were found to be more sensitive, they had their ranges reduced, thus eliminating portions of the parameter ranges which were causing many unrealistic predictions.

Figure 4 shows the cumulative GLUE distribution for mean streamflow using both the suggested parameter ranges and the new parameter ranges. Each line represents 5000 Monte Carlo runs. Comparing the difference between the two lines it was found that the reduction of ranges caused a considerable impact on the results even though the mean observed data did not correspond exactly to the 50^th percentile of the cumulative probability distribution. This proves that range adjustment should be applied to make the output more accurate.

It is usual to establish uncertainty bounds for the GLUE analysis (Blasone et al., 2008), which were here defined as 5^th and 95^th percentiles of the distribution. These percentiles were 0.194 and 0.423 m³ s^–1, respectively, when the initial parameter ranges were used. In contrast, 0.242 and 0.472 m³ s^–1 were the values found for the same percentiles when the narrowed ranges were applied instead. These bounds included most of the iterations (58.2% of the iterations for the initial parameter ranges and 70.6% for the narrowed ranges), which indicates that the variation of the model parameters was capable of accounting for the total output uncertainty, therefore, measurement and model structure errors were balanced (Blasone et al., 2008). The final GLUE simulation presented results much better than those from the initial GLUE simulation. The median was equal to 0.352 m³ s^–1, differing only by 9.5% from the mean observed streamflow value, which demonstrates acceptable results from the GLUE methodology.

Uncertainty analysis has been successfully applied in many other studies using hydrologic models. Muleta and Nicklow (2005) and Arabi et al. (2007) used the GLUE methodology in order to analyze uncertainties with respect to streamflow and sediment yield estimates simulated through the SWAT model. They found that sediment yield predictions had more uncertainties compared to streamflow estimates. Blasone et al. (2008) performed an uncertainty analysis for the MIKESHE model and concluded that uncertainties in parameters might not address the total uncertainty of spatially-distributed variables. These researchers point out that this occurs when a bias is observed in predictions due to uncertainties which can arise from both model structure and measurement error, as occurred in their study for groundwater elevation. Wei et al. (2008) computed model predictive uncertainties for the RHEM model to assess erosion risk for different scenarios. The authors were able to provide different conservation plans for decision makers instead of a single value concerning predicted soil loss. The uncertainties related to the VIC-3L model parameters and their effect on simulated streamflow values were evaluated successfully by Huang and Liang (2006).

The four most sensitive parameters found in this study cannot be easily measured in the field. Huang and Liang (2006) recommended for parameters like these that they should be estimated through model calibration. Using a conservative threshold (T < 0), only the parameters K_B and A_c had their ranges narrowed. However, if results from GLUE were not satisfactory, other thresholds would have to be used, for instance C_NS and log (C_NS) < 0.15 or < 0.3 (Table 3), in order to reduce parameter ranges and to proceed with a new GLUE analysis.

Before narrowing the suggested parameter ranges, the initial Monte Carlo simulation resulted in a mean streamflow equal to 0.23 ± 0.01 m³ s^–1 for a 95% confidence interval, which means that 95% of the 5,000-iteration sets will result in mean streamflow between 0.22 and 0.24 m³ s^–1. The final Monte Carlo simulation produced a mean streamflow of 0.31 ± 0.01 m³ s^–1 when the reduced parameter ranges were taken. These results clearly indicate that the reduction of bounds influenced considerably the mean streamflow values, thus improving mean streamflow predictions by 35% and reducing uncertainties linked to the input parameters. Although this study did not focus on the calibration of the LASH model, Figure 5 was added to allow readers to analyze the good performance of the model for JEW. This type of analysis can be more efficient if a given watershed has more data available to be applied to (Muleta and Nicklow, 2005). This analysis should be made for other watersheds before assuming the same parameters are the most sensitive as well as their ranges, since parameter ranges depend on specific characteristics of the watershed. Moreover, since this model is being applied for the first time, the data from other watersheds should be used in future for comparison.

Conclusions

The sensitivity analysis was performed for the LASH model with data from a medium-sized Brazilian watershed, indicating that the most sensitive parameters were K_B, A_C, θ₀, and γ. Using a conservative threshold (C_NS or log (C_NS) < 0) the ranges of two parameters were reduced, thus improving considerably the uncertainty analysis through GLUE methodology and making results more accurate for the watershed under consideration. Moreover, it will be easier to choose which parameters and their respective ranges should be considered to perform the calibration of LASH for other watersheds. These factors will help speed up the optimization efforts.

Acknowledgements

To CNPq and FAPEMIG (PPM IV-060/10) for financial support of this research and scholarships to the first and second authors and to the USDA-ARS National Soil Erosion Research Laboratory at Purdue University for supporting this research.

Received March 03, 2010

Accepted September 30, 2010

• Allen, R.G.; Pereira, L.S.; Raes, D.; Smith, M. 1998. Crop Evapotranspiration: Guidelines for Computing Crop Water Requirements. FAO, Roma, IT. (FAO Irrigation and Drainage Paper, 56).
• Almeida, A.C.; Soares, J.V.; Landsberg, J.J.; Rezende, G.D. 2007. Growth and water balance of Eucalypitus grandis hybrid plantations in Brazil during a rotation for pulp production. Forest Ecology and Management 251: 10-21.
• Arabi, M.; Govindaraju, R.S.; Hantush, M.M. 2007. A probabilistic approach for analysis of uncertainty in the evaluation of watershed management practices. Journal of Hydrology 333: 459-471.
• Benaman, J.; Shoemaker, C.A. 2004. Methodology for analyzing ranges of uncertain model parameters and their impact on total maximum daily load process. Journal of Environmental Engineering 130: 648-656.
• Beskow, S.; Mello, C.R.; Coelho, G.; Silva, A.M.; Viola, M.R. 2009a. Surface runoff in a watershed estimated by dynamic and distributed modeling. Revista Brasileira de Ciência do Solo 33: 169-178. (in Portuguese, with abstract in English).
• Beskow, S.; Mello, C.R.; Norton, L.D.; Curi, N.; Viola, M.R.; Avanzi, J.C. 2009b. Soil erosion prediction in the Grande River, Brazil using distributed modeling. Catena 79: 49-59.
• Beven, K.J.; Binley, A.M. 1992. The future of distributed models: model calibration and uncertainty prediction. Hydrological Processes 6: 279-298.
• Beven, K.J.; Freer, J. 2001. Equifinality, data assimilation, and uncertainty estimation in mechanistic modeling of complex environmental systems using the GLUE methodology. Journal of Hydrology 249: 11-29.
• Blasone, R.; Madsen, H.; Rosbjerg, D. 2008. Uncertainty assessment of integrated distributed hydrological model using GLUE with Markov chain Monte Carlo sampling. Journal of Hydrology 353: 18-32.
• Collischonn, W.; Allasia, D.G.; Silva, B.C.; Tucci, E.M. 2007. The MGB/ IPH model for large-scale rainfall-runoff modeling. Hydrological Science Journal 52: 878-895.
• Duan, Q.; Sorooshian, S.; Gupta, V. 1992. Effective and efficient global optimization for conceptual rainfall-runoff models. Water Resources Research 28: 1015-1031.
• Duan, Q.; Sorooshian, S.; Gupta, V. 1994. Optimal use of the SCE-UA global optimization method for calibrating watershed models. Journal of Hydrology 158: 265-284.
• Huang, M.; Liang, X. 2006. On the assessment of the impact of reducing parameters and identification of parameter uncertainties for a hydrologic model with applications to ungauged basin. Journal of Hydrology 320: 37-61.
• Mello, C.R. ; Viola, M.R.; Norton, L.D.; Silva, A.M.; Weimar, F.A. 2008. Development and application of a simple hydrologic model simulation for a Brazilian headwater basin. Catena 75: 235-247.
• Menezes, M.D.; Curi, N.; Marques, J.J.G.S.M.; Mello, C.R.; Araújo, A.R. 2009. Pedologic survey and geographic information system for evaluation of land use within a small watershed, Minas Gerais State, Brazil. Ciência e Agrotecnologia (UFLA) 33: 1544-1553. (in Portuguese, with abstract in English).
• Mishra, S.K.; Singh, V.P.; Sansalone, J.J.; Aravamuthan, V. 2003. A modified SCS-CN method: characterization and testing. Water Resources Management 17: 37-68.
• Mishra, S.K.; Sahu, R.K.; Eldho, T.I.; Jain, M.K. 2006. An improved I_a-S relation incorporating antecedent moisture in SCS-CN methodology. Water Resources Management 20: 643-660.
• Morris, M.D. 1991. Factorial sampling plans for preliminary computational experiments. Technometrics 33: 161-174.
• Muleta, M.K.; Nicklow, J.W. 2005. Sensitivity and uncertainty analysis coupled with automatic calibration for a distributed watershed model. Journal of Hydrology 306: 127-145.
• Nash, J.E.; Sutcliffe, J.E. 1970. River flow forecasting through conceptual models. Part I. a discussion of principles. Journal of Hydrology 10: 282-290.
• Rawls, W.J.; Ahuja, L.R.; Brakensiek, D.L.; Shirmohammadi, A. 1993. Infiltration and soil water movement. p. 1-51. In: Maidment, D.R., ed. Handbook of hydrology. McGraw-Hill, New York, NY, USA.
• Saltelli, A.; Tarantola, S.; Campolongo, F.; Ratto, M. 2004. Sensitivity analysis in practice: a guide to assessing scientific models. John Wiley, Chichester, UK.
• Stedinger, J.R.; Vogel, R.M.; Lee, S.U.; Batchelder, R. 2008. Appraisal of the generalized likelihood uncertainty estimation (GLUE) method. Water Resources Research 44: 1-17.
• Tucci, C.E.M. 2005. Hydrological models. Editora Universidade/ UFRGS, Porto Alegre, RS, Brazil. (in Portuguese).
• Van Griensven, A.; Meixner, T.; Grunwald, S.; Bishop, T.; Diluzio, M.; Srinivasan, R. 2006. A global sensitivity analysis tool for the parameters of multi-variable catchment models. Journal of Hydrology 324: 10-23.
• Viola, M.R.; Mello, C.R.; Acerbi Júnior, F.W.; Silva, A.M. 2009. Hydrologic modeling in the Aiuruoca river basin, Minas Gerais State. Revista Brasileira de Engenharia Agrícola e Ambiental 13: 581-590. (in Portuguese, with abstract in English).
• Wei, H.; Nearing, M.A.; Stone, J.J.; Breshears, D.D. 2008. A dual Monte Carlo approach to estimate model uncertainty and its application to the Rangeland Hydrology and Erosion Model. Transactions of the American Society of Agricultural and Biological Engineers 51: 515-520.
• White, K.L.; Chaubey, I. 2005. Sensitivity analysis, calibration, and validations for a multisite and multivariable SWAT model. Journal of the American Water Resources Association 41: 1077-1089.
• Zhou, M.C.; Ishidaira, H.; Hapuarachchi, J.M.; Kiem, A.S.; Takeuchi, K. 2006. Estimating potential evapotranspiration using Shuttleworth-Wallace model and NOAA-AVHRR NDVI data to feed a distributed hydrological model over the Mekong River basin. Journal of Hydrology 327: 151-173.

Publication Dates

• Publication in this collection
  11 July 2011
• Date of issue
  June 2011

History

• Received
  03 Mar 2010
• Accepted
  30 Sept 2010