Multiobjective sensitivity analysis and optimization of a distributed hydrologic model

Introduction Conclusions References


Introduction
With the development of information technology (e.g., high performance computing cluster and remote sensing technology), there has been a prolific development of integrated, distributed and physically-based watershed models (e.g., MIKE-SHE, Refsgaard and Storm, 1995) over the past two decades, which are increasingly being used to support decisions about alternative management strategies in the areas of land use change, climate change, water allocation, and pollution control.Though in principle parameters of distributed and physically based models should be assessable from catchment data (in traditional conceptual rainfall-runoff models, parameters are obtained through a calibration process), these models still need a parameter calibration process in practice due to scaling problems, experimental constraints, etc. (Beven and Binley, 1992;Gupta et al., 1998;Madsen, 2003).Problems, arising in calibrating distributed hydrologic models, include how to handle large number of distributed parameters and optimization problems with multiple but often conflicting objectives.
In the literature, to deal with large number of distributed model parameters, this is often done by aggregating distributed parameters (e.g., Yang et al., 2007), or screening out the unimportant parameters through a sensitivity analysis (e.g., Muleta and Nicklow, 2005;Yang, 2011).Sensitivity analysis can be used to not only screen out the most insensitive parameters, but also study the system behaviors identified by parameters and their interactions, qualitatively or quantitatively.However, most of applications in environmental modelling are based on the one-at-a-time (OAT) local sensitivity analysis, which is "predicated on assumptions of model linearity which appear unjustified in the cases reviewed" (Saltelli and Annoni, 2010) fashion in hydrologic modelling.This is not only due to the increasing availability of multi-variable (e.g., flow, groundwater level, etc.) or multi-site measurements, but also due to the intrinsic different system responses (e.g., peaks and baseflow in the flow series).Instead of finding a single optimal solution in the single objective optimization (SOO), the task in the multiobjective optimization (MOO) is to identify a set of optimal trade-off solutions (called a Pareto set) between conflicting objectives.In hydrology, the traditional method to solve multiobjective problems is to form a single objective, e.g., by giving different weights to these multiple objectives or applying some transfer function.Over the past decade, several MOO algorithms approaches have been applied to the conceptual rainfall-runoff models (e.g., Yapo et al., 1998;Gupta et al., 1998;Madsen, 2000;Boyle et al., 2000;Vrugt et al., 2003;Liu and Sun, 2010), and now increasing applied to distributed hydrologic models (e.g., Madsen, 2003;Bekele and Nicklow, 2007;Shafii and Smedt, 2009;MacLean et al., 2010).And there are some papers (Tang et al., 2006;Wöhling et al., 2008) to comparatively study their strengths with the application in hydrology.It is worth noting that the multiobjective calibration is different from statistical uncertainty analysis which is based on the concept (or similar concept) of "equifinality" (see discussion in Gupta et al., 1998;Boyle et al., 2000).This paper applies two sensitive analysis techniques (Morris method and State Dependent Parameter method) and ε-NSGAII in the multiobjective sensitive analysis and calibration framework.This was implemented to calibrate a distributed hydrological and discussed; and finally the main results are summarized and conclusions are drawn in "conclusions" section.
2 Hydrologic model MOBIDIC MOBIDIC (MOdello di Bilancio Idrologico DIstribuito e Continuo; Castelli et al., 2009;Campo et al., 2006) is a distributed and raster-based hydrological balance model.MO-BIDIC simulates the energy and water balances on a cell basis within the watershed.
Figure 1 gives a schematic representation of MOBIDIC.The energy balance is approached by solving the heat diffusion equations in multiple layers in the soil-vegetation system, while the water balance is simulated in a series of reservoirs (i.e., boxes in Fig. 1) and fluxes between them.
For each cell, water in the soil is simulated by where W g [L] and W c [L] are the water contents in the soil gravitational storage and capillary storage, respectively, and and S as [LT −1 ] are infiltration, percolation, interflow, evaporation, and adsorption from gravitational to capillary storage, which are modeled through following equations: (2) Full Once the surface runoff (Q h and R d ) and baseflow are calculated, three different methods can be used for river routing, i.e., the lag method, the linear reservoir method, Muskingum-Cunge method (Cunge, 1969).Muskingum-Cunge method was used in this study.MOBIDIC uses either a linear reservoir or the Dupuit approximation to simulate the groundwater balance which relates the groundwater change to the percolation, water loss in aquifers and baseflow.In this case study, the linear reservoir method was used.
Although there are many distributed parameters in MOBIDIC, normally these distributed parameters are calibrated through the "aggregate" factors (e.g., the multiplier for hydraulic conductivity) based on their initial estimations.And hereafter we use the term "factor" (instead of "model parameter") when we conduct the sensitivity analysis and optimization, to avoid the confusion with the term "model parameter" used in model description.A factor can be a model parameter or a group of model parameters, and in this paper it is a change to be applied to a group of model parameters.In MO-BIDIC, normally nine factors (i.e., nine groups of parameters) need to be calibrated.These factors, their explanations, and their corresponding model parameters are listed in Table 1.

Methodology
The procedure applied here consists of two-step analyses, i.e., a multiobjective sensitivity analysis generally characterizing the basic hydrologic processes and single out the most insensitive parameters, and a multiobjective calibration aiming at trade-offs between different objective functions.

Sensitivity analysis techniques
Sensitivity analysis is to assess how variations in model out can be apportioned, qualitatively or quantitatively, to different sources of variations, and how the given model depends upon the information fed into it (Saltelli et al., 2008).In the literature, a lot of sensitivity analysis methods are introduced and applied, e.g., Yang (2011) applied and compared five different sensitivity analysis methods.Here we adopted an approach which combines two global sensitivity analysis techniques, i.e., the Morris method (Morris, 1991) and SDP method (Ratto et al., 2007).

Morris method
Morris method is based on replicated and randomized one-factor-at-a-time design (Morris, 1991).For each factor X i , Morris method uses two statistics, µ i and σ i , which measure the degree of factor sensitivity, and the degree of nonlinearity or factor interaction, respectively.The higher µ i is, the more important the factor X i is to the model output; and the higher σ i is, the more nonlinear the factor X i is to the model output or more interactions with other factors (details refer to Morris, 1991;Campolongo et al., 2007).Morris method takes m • (n + 1) model runs to estimate these two sensitivity indices for each of n factors with sample size m.The advantage is it is efficient and effective to screen out insensitive factors.Normally m takes values around 50.And according to Saltelli et al. (2008), the sensitivity measure (µ i ) is a good proxy for the total effect (i.e., S Ti in Eq. (4) below), which is a robust measure in sensitivity analysis.

State-Dependent Parameter method (SDP)
SDP (Ratto et al., 2007) is based on the ANOVA functional decomposition, which apportions the model output uncertainty (100 %, as 1 in Eq. where S i is the main effect of factor X i representing the average output variance reduction that can be achieved when factor X i is fixed, and S i j is the first-order interaction between X i and X j , and so on.In ANOVA based sensitivity analysis, total effect (S Ti ) is frequently used, which stands for the average output variance that would remain as long as X i stays unknown, SDP method uses the emulation technique to approximate lower order sensitivity indices in Eq. ( 3) (e.g., S i and S i j in this study) by ignoring the higher order sensitivity indices.And we define S Di = S i + j S i j (referred to as "quasi total effect" later) as a surrogate to the total effect.The advantage is that it can precisely estimate lower order sensitivity indices at a lower computational cost (normally 500 model runs, which is independent of number of factors).The disadvantage is that it cannot estimate higher order sensitivity indices.
Since these two methods are computationally efficient, these two methods are applied individually.And then, the sensitivity of each factor and its system behaviour will be discussed, qualitatively by Morris method, and quantitatively by SDP method.And then the most insensitive factors will be screened out and excluded in the calibration.
In the context of multiobjective analysis, sensitivity analysis includes: (1) to examine

Multiobjective calibration and ε-NSGAII
In the literature of hydrologic modelling, most applications are single objective based, which aims at a single optimal solution.However, for example in flow calibration, there is always a case that two solutions, one solution better simulates the peaks and poorly simulates the baseflow while the other solution poorly simulates the peaks while better simulates the baseflow.These two solutions, which are called Pareto solutions, are incommensurable, i.e., better fitting of the peaks will lead to worse fitting of the baseflow, and vice versa.This belongs to the domain of MOO, aiming at finding a set of optimal solutions (Pareto solutions), instead of one single solution.
Generally a MOO problem can be formulated as follows: min Where X is an n-dimensional vector and in this study represents the model factors to be calibrated, f i (X) i th objective function, and g i (X) i th constraint function.
ε-NSGAII is an extension of the NSGAII (Deb et al., 2002), a second generation of multiobjective evolution algorithm.The main characteristics of ε-NSGAII include: (i) selection, crossover, and mutation processes as other genetic algorithm by mimicking the process of natural evolution, (ii) an efficient non-domination sorting scheme, (iii) an elitist selection method that greatly aids in capturing Pareto front, (iv) ε-dominance archiving, (v) adaptive population sizing, and (vi) automatic termination to minimize the Introduction

Conclusions References
Tables Figures

Back Close
Full need for extensive parameter calibration.More details refer to Kollat and Reed (2006).
As MOO is very time-consuming, we parallelized the source code and interfaced it with MOBIDIC.
As a comparison, a single objective function is defined as 2-norm of the multiple objectives F (X), which measures how close to the original point (theoretical optimum And SOO was done with the classic Nelder-Mead algorithm (Nelder and Mead, 1965) which is widely used in MOBIDIC applications.
To analyze the Pareto solution and also compare with the solution from SOO, except for traditional methods, the "Level diagrams" proposed by Blasco et al. (2008) was also used.Compared to traditional methods, it can visualize high dimensional Pareto front and synchronizes the objective and factor diagrams.The procedure and includes two steps.In the first step, the vector of objectives (k-dimension) for each Pareto point is mapped to a real number (one-dimension) according to the proximity to the theoretical Data used in MOBIDIC model include (i) Digital Elevation Model (DEM), (ii) soil data, (iii) land cover data, (iv) climate data (precipitation, minimum and maximum temperature, solar radiation, humidity and wind speed), and (v) flow data.9 m DEM, land cover, SSURGO soil data, one station (Davidson river near Brevard) of flow data are from US Geological Survey, and hourly NLDAS climate data from National Aeronautics and Space Administration (NASA).NLDAS integrates a large quantity of observation-based and model reanalysis data to drive offline (not coupled to the atmosphere) land-surface models (LSMs), and executes at 1/8th-degree grid spacing over central North America, enabled by the Land Information System (LIS) (Kumar et al., 2006;Peters-Lidard et al., 2007).
DEM is used to delineate the watershed and estimate the topographic parameters and river system, Land cover for evaporation parameters, soil data for soil parameters, climate data is used to drive MOBIDIC, and flow data are used to calibrate the model and assess model performance.The climate and flow data used in this study are from 1 January 1996 to 30 September 2006.As NLDAS only has hourly temperature daily instead of hourly minimum and maximum temperature needed by MOBIDIC, we compiled the hourly climate data to daily data and run the model at a daily step.After MOBIDIC setup, the initial parameter values are listed in third column of Table 1.
We split the data into a warm-up period (from 1 January 1996 to 30 September 2000), a calibration period (from 1 October 2000 to 30 September 2003), and a validation period (from 1 October 2003 to 30 September 2006).

Objective function selection
After setting up MOBIDIC in the Davidson watershed, three objective functions were used in the multiobjective sensitivity analysis and optimization: Introduction

Conclusions References
Tables Figures

Back Close
Full (1) Standardized root mean square error between the logarithms of simulated and observed outflows: (2) Water Balance Index, calculated as the mean absolute error between the simulated and observed flow accumulation curves: (3) Mean absolute error between the logarithms of simulated and observed flow duration curves In Eqs. ( 7)-( 9), Q obs SRMSE (Eq.7), WBI (Eq.8), and MARD (Eq.9) are measures of the closeness between simulated and observed flow series, water balance, and closeness between simulated and observed flow frequencies, respectively.The smaller these measures are, the better the simulation is, and the minima are (0, 0, 0) meaning a perfect match Introduction

Conclusions References
Tables Figures

Back Close
Full between the simulation and observation.It is worth noting that we use the logarithms of the flows instead of flows to avoid overfitting flow peaks (Boyle et al., 2000;Shafii and De Smedt, 2009).And for SRMSE, we have NS ≈ 1 − SRMSE 2 when N is large (e.g., > 100), where NS is the Nash-Sutcliffe coefficient (Nash and Sutcliffe, 1970), which is widely used in hydrologic modelling.
And accordingly, the single objective function here is the Euclidean norm (2-norm) of SRMSE, WBI, and MARD: 5 Result and discussion

Multiobjective sensitivity analysis
Morris method and SDP method were applied individually to the initially selected factors (in Table 1).For Morris method, its convergences for three objective functions, monitored using the method proposed in Yang (2011), were achieved around 700 ∼ 800 model simulations.Figure 3 gives the sensitivity results for objective functions SRMSE, WBI, and MARD, respectively.In each plot, the horizontal axis (µ) denotes the degree of factor sensitivity, and the vertical axis (σ) denotes the degree of factor nonlinearity or interaction with other factors.
For SRMSE, the most sensitive parameters are group (pα, pγ, and pκ), followed by pβ and rCH, while other parameters (especially rK s and rK f ) are not so sensitive.This applies to the degree of the factor nonlinearity or interaction.The sensitivities of pα, pγ, and pκ indicate the importance of their corresponding processes (i.e., surface runoff, percolation, and adsorption which is related to evapotranspiration) to SRMSE, while interflow (pβ) is less important and other processes/characteristics (e.g., groundwater flow, rK f ) are not important.Introduction

Conclusions References
Tables Figures

Back Close
Full For WBI, the dominating parameter is pκ, followed by group (pα, pγ, pβ and rCH), while other parameters (especially rK f and rW cmax ) are not so sensitive.WBI measures the water balance between observed and simulated flow series, and it is reasonable that pκ which controls the water supply for evaporation is most sensitive while other factors (pα, pγ, pβ and rCH) are sensitive mainly through interaction with this factor, as indicated by the high "σ"s of these factors.
For MARD, the results are nearly the same to SRMSE. Figure 4 gives the sensitivity results based on SDP method for SRMSE, WBI, and MARD, from top to bottom.In each plot, the grey and black bars are S i and S Di for each factor.
For SRMSE, as indicated by R 2 in the legend, main effects (S i ) contribute to 58.7 % of SRMSE uncertainty, and quasi total effects (S Di ) account for 83 % of SRMSE uncertainty which is quite high, while other 17 % due to higher interactions are not explained.Based on S Di (black bar), the most sensitive factors are pγ and pκ, followed by pα and rCH, and then pβ and rW cmax while other factors are not sensitive.This result quantitatively corroborates the result obtain from Morris method.The main effects (S i ) of (pγ, pκ, and pα) are high (i.e., 0.17, 0.18 and 0.14), which suggests these factors should be determined first in model calibration as they lead to the largest reduction in SRMSE uncertainty.For each factor, the difference between the black bar and grey bar shows the first order interaction with other factors.This interaction is very strong in pγ, pκ, pα, and rCH, and very weak in other factors.
For WBI, as indicated by R 2 in the legend, the total main effects (S i ) contribute to 38.4 % of the WBI uncertainty, quasi total effects (S Di ) only account for 57.6 % of WBI uncertainty, and 40 % due to higher interactions are not explained and can not be ignored.However, by analyzing the result with that from Morris method (top right in Fig. 3), we still can get some valuable results: the dominating sensitive factor is pκ with S Di 0.43 (which is same as Morris method), followed by pγ, pα, and rCH, while other factors are not sensitive; the main effect of pκ is as high as 0.27, and it should be fixed Figures in order to get the maximum reduction in WBI uncertainty; the first interaction is high in pκ, pγ, and pα, not obvious in other factors.Similar to the Morris results for SRMSE and MARD, the result of MARD is nearly the same as SRMSE.The similar result for SRMSE and MARD shows similar characteristic relationship between factors and the objective function.This is explainable: a good simulation measured by SRMSE will more likely result in a good measure of MARD, and vice versa.
As aforementioned, in the context of multiobjective sensitivity analysis, sensitivity analysis is to exclude the factors which are insensitive to all the objective functions considered.Based on the analysis above, four most insensitive factors are rKs, rKf, rW cmax ,and rW gmax .However, as shown in Fig. 4, rW cmax is more sensitive than other three factors, and for objective function WBI, as higher order interactions are strong based on SDP (i.e., explains around 40 % of model uncertainty), and evaporation is the most sensitive process to water balance (as indicated by pκ and rCH) and rW cmax is the only factor related to evaporation storage (W c ), therefore, we only exclude rK f , rK s , and rW gmax for calibration.

Multiobjective optimization
After sensitivity analysis, only six factors were involved in the calibration.MOO converged with 482 Pareto front points after totally 22 000 model runs with the parallelized ε-NSGAII, while SOO converged after 686 model runs with the classic Nelder-Mead algorithm.Apparently, ε-NSGAII took more model simulations than the Nelder-Mead algorithm, but simulation time was compensated by the parallelized code running on high performance clusters.
Figure 5 shows optimized non-dominant sets normalized within [0, 1] and the black line is for the factor set with SOO.It is encouraging that except rW cmax , factor ranges decreased a lot.This corroborates the conclusion in the sensitivity analysis: pγ, pκ, pβ, pα, and rC H are the most sensitive and identifiable factors to these three objective functions, while rW cmax is less sensitive and less identifiable.values of pγ and dispersed rW gmax show that optimized factor sets are scattered in the response surface rather than concentrated in a continuous region.And the factor set with SOO is within the range of non-dominant sets.
Figure 6 shows Pareto solutions scattered in the three-dimensional space (top left), and projections in two-dimensional subspaces with corresponding correlation coefficients (r) in the calibration period, with the black dot in each plot denoting the solution for SOO.Correlation coefficients are high and negative for SRMSE and WBI (−0.54), and WBI and MARD (−0.74), and this indicates strong trade-off interactions along the Pareto surface, i.e., better (lower) WBI will eventually result in worse (higher) SRMSE, and vice versa.The correlation coefficient is low (0.13) between SRMSE and MARD, and it is close to 0 when these two variables approach to their minima regions (i.e., SRMSE < 0.53 and MARD < 0.09).Table 2 lists the statistics of these three objectives associated with Pareto sets and the result of SOO.For Pareto sets, in the calibration period, the average SRMSE is 0.49 ranging from 0.47 to 0.57, which corresponds to the average NS 0.78 ranging from 0.67 to 0.78; the average WBI is 0.05 ranging from 0.02 to 0.11; and the average MEAD is 0.08 ranging from 0.03 to 0.11.In the validation period, the average SRMSE is 0.54 ranging from 0.51 to 0.62, which corresponds to the average NS 0.70 ranging from 0.61 to 0.74; the average WBI is 0.05 ranging from 0.04 to 0.09; and the average MEAD is 0.10 ranging from 0.08 to 0.13.And for SOO, SRMSE, WBI and MEAD are 0.48, 0.06 and 0.07 for the calibration period, and 0.57, 0.06 and 0.10 for the validation, and accordingly the "NS"s are 0.77 and 0.67, respectively.According to Moriasi et al. (2007) which suggests NS > 0.75 and WBI < 10 % as excellent modelling of river discharge, all Pareto solutions with MOO and the solution with SOO are close to "excellent" for both calibration and validation periods.
To better visualize Pareto sets and compare with the result of SOO, the level diagrams are plotted in Fig. 7 by applying Euclidean norm (2-norm) to evaluate the distance of each Pareto point to the ideal origin (0,0,0) (ideal values for all three normalized objectives are 0).In Fig. 7, top three plots are for three objectives and the rest for optimized factors, and the black dot in each plot is the solution for SOO.In the level Introduction

Conclusions References
Tables Figures

Back Close
Full diagrams, each objective and each factor of a point (corresponding to a Pareto solution) is represented with the same 2-norm value for all the plots.Compared with MOO, obviously, SOO was trapped in the local optima as seen in top-left plot, which means optimization with Nelder-Mead algorithm was dependent of starting point.The 2-norm has a close linear relationship with SRMSE due to values of SRMSE are 5 to 10 times of other two objective functions, and it does not have such relationship with other two objectives.The scattering of objectives and factors makes it difficult in decision making to select a single solution because there is not a clear trade-off solution (Blasco et al., 2008).However, compared to SOO, the Pareto solutions from MOO can make decision making easy as it can be converted with expert opinion or some utility function.
Figures 8 and 9 show simulated and observed flow duration curves and time series flows, respectively, with grey lines denoting the simulations with MOO and black lines with SOO.Generally, all simulations match the observation well for both the duration curve and time series flow for both calibration period and validation period.For the duration curve, simulations from MOO show a wide range in the low flows with frequencies from 0.85 to 1.0, which reflects the insensitivity of groundwater process (discussed in the sensitivity analysis, i.e., rK f is insensitive to these three objectives).Except for this, there is a slight overestimation of flows, large flows during the calibration period with frequencies from 0.2 to 0.1, and median to large flows during the validation period with frequencies from 0.5 to 0.1.This might be due to the uncertainty in the reanalyzed climate data.And the extreme flow with frequency around 0 is underestimated, and this is because we chose the logarithm scale of the observed and simulated flows instead of normal scale when computing objectives SRMSE and MARD.With SOO, the deviation from the observed is larger.Similar conclusions can be drawn from the time series simulations in Fig. 9, i.e., the wide ranges of low flow period, and underestimation of flow peaks.Other than this, generally all simulations can mimic the observations.
Figure 10 shows the time series of watershed average storages (soil storage expressed as soil saturation, and groundwater depth), and fluxes (evaporation, surface runoff and baseflow) associated with MOO (shaded)  soil saturation varies from 0.2 to 1.0 and groundwater from 0 to 120 mm.The temporal fluctuation of soil moisture is higher than groundwater, but lower than fluxes in evaporation and surface runoff.And this is true for the solution with SOO except the its ranges of soil saturation and groundwater (groundwater is very close to 0 mm).For fluxes with MOO, evaporation and surface runoff have more temporal variation than baseflow, and their magnitudes are larger than baseflow.This applies to fluxes with SOO, and its baseflow is close to 0. This can be confirmed by the De Finetti diagram in Fig. 11: with MOO, the average contributions of evaporation, surface runoff, and basedflow are 49.3, 46.1, and 4.8 %, respectively while the contribution of baseflow is very insignificant.And the contribution of baseflow is almost 0 with SOO.

Conclusion
This study presents a multiobjective sensitivity and optimization approach to calibrate a distributed hydrologic model MOBIDIC with its application in the Davidson watershed for three objective functions (i.e., SRMSE, WBI, and MAED).Results show: 1.The two sensitivity analysis techniques are effective and efficient to determine the sensitive processes and insensitive parameters: surface runoff and evaporation are very sensitive processes to all three objective functions, while groundwater recession and soil hydraulic conductivity are not sensitive and were excluded in the optimization.average NS is 0.75 ranging from 0.67 to 0.78 in the calibration and 0.70 ranging from 0.61 to 0.74 in the validation period.
4. In the case study, evaporation and surface runoff shows similar importance to watershed water balance while the contribution of baseflow can be ignored.
5. Compared to MOO with ε-NSGAII, the application of SOO with the Neld-Mead algorithm was dependent of initial starting point.Furthermore, the Pareto solution provides a better understanding of these conflicting objectives and relations between objectives and parameters, and a better way in decision making.Introduction

Conclusions References
Tables Figures

Back Close
Full , or simple linear regressions where a lot of uncertainty are not fairly accounted for.The use of global sensitivity analysis techniques is very crucial in distributed modelling.Only recently, global sensitivity analysis techniques started to appear in hydrologic modelling.Although most hydrologic applications are based on the single objective calibration, model calibration with multiple and often conflicting objectives arises in a natural Discussion Paper | Discussion Paper | Discussion Paper | model MOBIDIC with its application to the Davidson watershed, North Carolina.The purpose is to study parameter sensitivity of the hydrologic model MOBIDIC and explore the capability of MOO in calibrating the MOBIDIC compared to the traditional SOO used in MOBIDIC applications.This paper is structured as follows: Sect. 2 gives a description of the MOBIDIC model; Sect. 3 introduces the approach in the multiobjective sensitivity analysis and optimization; Sect. 4 gives a brief introduction of the study site, model setup, objective selection, and sensitivity and calibration procedure; in Sect.5, the results are presented Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | where γ, β and κ are percolation coefficient [T −1 ], interflow coefficient [T −1 ], and soil adsorption coefficient [LT −1 ], respectively, P the precipitation [LT −1 ], Q h and R d Horton runoff and Dunne runoff, K s the soil hydraulic conductivity [L T −1 ], W gmax [L] and W cmax [L] the gravitational and capillary storage capacities.
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | levels of their interactions: the sensitivity of each factor to different objective functions, qualitatively or quantitatively; (2) to single out the most sensitive factors and study the physical behaviours of the system; (3) to exclude the most insensitive factors and therefore simplify the process of calibration.Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | optimum measured with a specific norm of objectives; and in the second step, these norm values are plotted against the corresponding values of each objective or factor.1-norm, 2-norm and ∞-norm are suggested.To compare with SOO, 2-norm was used.4 Davidson watershed and objective selection 4.1 Davidson watershed The Davidson watershed, located in southwest mountain area of North Carolina, drains an area of 105 km 2 above the station "Davidson river near Brevard" (see Fig. 2).The Figures Back Close Full Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | elevation ranges from 645 m to 1820 m a.s.l.Based on the NLDAS climate data, the average annual precipitation is 1900 mm and varies from 1400 mm to 2500 mm, and daily temperature changes from −19 • C to 26 • C. The average daily flow is about 3.68 m 3 s −1 .
Discussion Paper | Discussion Paper | Discussion Paper |

i
and Q sim i are observed and simulated flow series at time step i , N the data length, log Q the average of logarithmic transformed observed flows, Q obs Ci and Q sim Ci i th observed and simulated accumulated flows, and Q obs Pi and Q sim Pi i th percentiles of observed and simulated flow duration curves.
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | and SOO (black line).With MOO, Discussion Paper | Discussion Paper | Discussion Paper | 2. For SRMSE and MAED, though all the factors have almost same sensitivities, the non-dominance of Pareto set shows the trade-off between these two objectives in the response surface.3.Both MOO and SOO achieved acceptable results for both calibration period and validation period, in terms of objective functions and visual match between simulated and observed flows and flow duration curves.For example, with MOO, the Discussion Paper | Discussion Paper | Discussion Paper |

Factora
Figure 1 A schematic representation of MOBIDIC.Boxes denote different water storages (gravitational storage W g , capillary storage W c , groundwater storage H, surface storage W s , and river system), solid arrows fluxes (evaporation E t , precipitation P, infiltration I nf , adsorption A d , percolation P c , surface runoff R, interflow Q d , groundwater discharge Q g , and surface runoff and interflow from upper cells (R+Q d ) up ), dashed arrows different routings, and blue characters major model parameters.

Fig. 1 .
Fig. 1.A schematic representation of MOBIDIC.Boxes denote different water storages (gravitational storage W g , capillary storage W c , groundwater storage H, surface storage W s , and river system), solid arrows fluxes (evaporation E t , precipitation P , infiltration I nf , adsorption A d , percolation P c , surface runoff R, interflow Q d , groundwater discharge Q g , and surface runoff and interflow from upper cells (R + Q d ) up ), dashed arrows different routings, and blue characters major model parameters.

Figure 2 Fig. 2 .
Figure 2 The location of the Davidson watershed, North Carolina, with DEM map, river system (lines), and watershed outlet (the triangle point) Fig. 2. The location of Davidson watershed, North Carolina, with DEM map, river system (lines), and watershed outlet (the triangle point).

Figure 3
Figure 3Multiobjective sensitivity analysis result based on the Morris method

Fig. 3 .
Fig. 3. Multiobjective sensitivity analysis result based on the Morris method.

Figure 6 Fig. 6 .Fig. 8 .
Figure 6The Pareto solutions in the three dimensional space (top left), and the projections in the two dimensional subspace (other plots), with MOO, and the black dot is the solution with SOO

Figure 9
Figure 9 Observed flows (dotted) and simulated flows with MOO (grey) and SOO (black line) for the calibration period (top) and validation period (bottom) Fig. 9. Observed flows (dotted) and simulated flows with MOO (grey) and SOO (black line) for the calibration period (top) and validation period (bottom).

Figure 10
Figure 10Time series of watershed average storages (soil water storage expressed as soil saturation, and groundwater depth), and fluxes (evaporation, surface runoff, and baseflow) with MOO (grey) and SOO (black line).For SOO, the groundwater storage and baseflow are close to 0 and hardly seen.

Fig. 10 .
Fig. 10.Time series of watershed average storages (soil water storage expressed as soil saturation, and groundwater depth), and fluxes (evaporation, surface runoff, and baseflow) with MOO (grey) and SOO (black line).For SOO, the groundwater storage and baseflow are close to 0 and hardly seen.