2.1. Study Area
The study area is the SefidRood basin (SRB), which is situated in the north-western region of Iran (Fig. 1). Drainage of this basin, one of the most important basins in the country, is done by the SefidRood river. The SefidRood is the second largest river in Iran, with a total length of 670 km. The basin's topography ranges from − 31 m in basin outlet to 4300 m in the highest location above mean sea level. The basin is positioned between coordinate: 46◦ 30⸍ E and 51◦ 13⸍ E longitudes and 34◦ 55⸍ N and 37◦ 52⸍ N latitudes. The watershed covers about 59500 km2, and its perimeter is approximately 1650 km. Due to this basin's extent, the variety of climate, vegetation, and soil type are observed, which is one of the reasons for the importance of the SefidRood catchment. The observed river discharge data from the Gilvan hydrometric station located at the basin outlet in daily and monthly time scales are used to evaluate the VIC 3L hydrological model's performance.
2.2. Macro Scale Hydrological Model: VIC-3L
The Macro-scale hydrological model variable infiltration capacity (VIC-3L) (Liang, Wood, and Lettenmaier 1996) was used at a spatial resolution of 0.25◦×0.25◦ (approximately 27.7 km ×27.7 km). This model, capable of solving water and energy balances. In the VIC-3L model, each cell is modeled separately, and grid cells do not interact or current transfer with each other. Meteorological input data for model: precipitation (mm), humidity (vapor pressure), wind speed (m/s), daily maximum and minimum temperature (◦C) are drawn from local synoptic and rain gauge stations within the catchment, then application of Inverse Distance Weighting (IDW) interpolation method has been extracted for each cell of the basin computing network. Leaf area index, albedo coefficient, minimum stomatal resistance, architectural resistance, roughness length, displacement length, and the relative fraction of roots in each soil layer include the model parameters related to the catchment's vegetation. In this research, after classifying and determining the type of coverage based on Advanced Very High-Resolution Radiometer (AVHRR), these parameters have been estimated according to the values provided by (Su et al. 2005). Determination of some soil parameters (such as porosity fraction, saturated hydraulic conductivity, field capacity, bulk density, wilting point, the slope of retention curve) in order to implement the model is possible using library functions, for which it is necessary to determine the soil texture of the study catchment. In this study, soil types are obtained from the Iranian National Soil Conservation Service (INSCS) soil map. According to the values provided by (Rawls, Gimenez, and Grossman 1998), the characteristics and parameters of each soil texture are extracted. However, most of the model parameters can be derived from remote sensing and in-situ measurement observation. Some of them are conceptual and have non-physical values, so their estimation is feasible during the calibration process. The studies indicate the high importance of some soil parameters in estimating soil layers moisture and surface and subsurface runoff modeling (Guo, Liang, and Leung 2004; Nijssen et al. 1997): binf) The infiltration parameter which controls the partitioning of rainfall (or snowmelt) into infiltration and direct runoff [0.01–0.4]. Dsmax) The maximum base-flow velocity [0.01–30]. Ds) The fraction of maximum base-flow velocity where non-linear base flow begins [0.01–1]. Ws) The fraction of maximum soil moisture where non-linear base flow occurs [0.01–1]. C) Power used in base discharge curve [0.01–2]. d2 and d3) The thickness of the second and third soil layer [0.1–2.5 and 0.1–1.5].
2.3. Remote-sensing and Re-analysis data used
2.3.1. ASCAT Surface Soil Moisture Database
Advanced SCATterometer (ASCAT) is a radar system operating in the C-band (5.255 GHz) onboard the Metop satellite since 2006. The scientific basis for retrieving the soil moisture from ASCAT scatterometer data is based on a change detection method (CDM) (Bartalis et al. 2008). This method is mathematically less complicated than the radiation transmission model can be analyzed in inverse. Therefore, the soil moisture content can be estimated directly by measuring the sensor without needing the correction process. Thus, in order to estimate the recovery error for each cell at ground level is quite straightforward. It is noted that one of the disadvantages of the CDM is the integrated representation of the measurement process. In other words, there is no distinct portion to the retrieving and observation of soil, vegetation, and soil-vegetation interactions. The algorithm for processing the scatterometer data is called soil WAter Retrieval Package (WARP), Which is provided in the WARP6 software package (Naeimi et al. 2009). Soil moisture values of this sensor are associated with the depth of 2 cm of soil (topmost) with 3-days/ ~25 km temporal/ spatial resolution (Wanders et al. 2014)
2.3.2. A Global Water Resources Reanalysis Database: The eartH2Observe
In the re-analysis database using data assimilation techniques and meteorological forcing data, climate data are re-analyzed. In general, a constant model and available observational data are used during the re-analysis process. Climate re-analysis databases are continually being updated by monitoring data measured by satellites, radiosondes, and so on. In these databases, using observation data, data simulated by the models are corrected and debug. One of the features of this data is its spatial scale diversity. The Earth2Observe dataset (http://wci.earth2observe.eu/), recently developed by the EU Water Resources/Remote Sensing Research (2017), is one of these re-analysis databases, which consisted of an ensemble of ten global hydrological and land surface models with a different spatiotemporal resolution for the 1979–2014 period using a meteorological forcing dataset based on re-analysis. In this database, error reduction uses meteorological observation and data assimilation techniques (Dutra, Balsamo, and Calvet 2017). In the next section, the modeling protocol of global hydrological and land surface models in this database, which simulated data has been used in the present study, has been described.
2.3.2.1. The LISFLOOD Model
The LISFLOOD model is a distributed large-scale hydrological rainfall-runoff model that has been developed by the floods group of the Natural Hazards Project of the Joint Research Centre (JRC) of the European Commission to simulate water balance of catchment. The basis of this model is the division of the catchment into several cells. To implement the LISFLOOD model can be used a range of grid cells of as little as 100 m for medium-sized catchments to around 10 km for large scale catchments (Burek, Knijff, and De Roo 2013). LISFLOOD is one of the models that can simulate water balance (on a daily time step) and individual flood events (on hourly time intervals, or even smaller). In this model, two layers are considered for the simulation of groundwater and subsurface flow. Also, using a sub-model, surface runoff is routed to the nearest river channel. With the Arno routing model (Todini 1996; Zhao and Liu 1995), the flow is directed to the basin's outlet. The processes that are simulated by the model include snowmelt, infiltration, interception of rainfall, evapotranspiration, surface runoff, and exchange of soil moisture between the two soil layers, etc. One of the limitations of this hydrological model is its non-application in dry areas due to the lack of simulation of capillary ascent and groundwater at great depths (van der Knijff, Younis, and de Roo 2010). This model's surface soil moisture data have been used to calibrate the VIC 3L model in the present study.
2.3.2.2. The HTESSEL Model
Hydrology Tiled ECMWF Scheme for Surface Exchanges over Land (HTESSEL) is a distributed land surface model that simulates the land surface's response to weather conditions and components related to water and energy balance. The model can also simulate soil temperature changes, soil moisture, surface / sub-surface runoff, and improve the model's performance by considering different vegetation and snow cover conditions (Dutra et al. 2017). In general, in this model, for the surface layer, six tiles including bare ground, low and high vegetation, intercepted water and shaded and exposed snow, and two over water including open and frozen water (Balsamo et al. 2009). The vertical discretization has soil layers below ground at 7 cm for the top layer and 21, 72, and 189 cm as suggested in (Deardorff 1978) that can be covered by a single layer of snow (The presence of a snow layer depends on the study area, which can be zero). In the HTESSEL model, Subsurface water fluxes are determined by Darcy's law. In this research, surface soil moisture data of this model have been used.
2.3.2.3. The HBV Model
The HBV model is a semi-distributed conceptual model generated by the Swedish Meteorological and Hydrological Institute (SMHI) by (Bergström and Forsman 1973) to simulate runoff, forecasting, and water balance studies. The model simulates daily discharge using rainfall and temperature (daily) and potential evaporation (monthly) as input. This model with a simple structure and low inputs simulate the most critical variables relating to water balance components, including the river discharge, precipitation, evapotranspiration, soil moisture, and water storage that contribute to runoff formation. In this model, the basin is divided into hydrological units based on height variations and vegetation type. Snowmelt simulation in the HBV model is computed by a degree-day method. The HBV-light version provides the possibility to include observed groundwater levels and use a different response routine with a delay parameter. In this study, surface soil moisture and evapotranspiration data of this model have been used.
2.3.2.4. The W3RA Model
The World-Wide Water (W3RA) model is one of the distributed large-scale hydrological models based on the landscape hydrology component model of the AWRA system, developed by the Australian Meteorological Agency. This model's snowmelt simulation is based on the HBV model's thought, and the model is implemented at a spatial resolution of 0.5 degrees (Dutra et al. 2017). The vertical discretization is considered a 3-layer soil with different depths. In each grid box, two vegetation types are present: deep-rooted tall vegetation type and shallow-rooted short vegetation type. The energy balance determines actual evaporation. In this model, humans' changes and systems, such as irrigation, are not considered (Van Dijk, 2010). In this study, evapotranspiration data of this model have been used.
2.3.2.5. The GLEAM Model
The Global Land Evaporation Amsterdam Model (GLEAM) was developed by VU Amsterdam in collaboration with Ghent University in Belgium and the European Space Agency to estimate evapotranspiration (Miralles et al. 2011). In each grid cell includes 4 types of coverage: bare soil, short vegetation, tall vegetation, open water. The amount of evapotranspiration in this model is calculated based on the equation presented by (Priestley and Taylor 1972) and according to air temperature and net surface radiation. The actual evaporation for each pixel results from aggregating the total of each land cover type weighted by the percentage of each cover type within the pixel. GLEAM data is included daily with a spatial resolution of 0.25 degrees from 1980 to 2016.
2.3.2.6. The SURFEX-TRIP Model
The SURFEX-TRIP model is a distributed model coupling with the ISBA land surface model to compute water and energy balance. In the SURFEX-TRIP model, 14 layers are considered for the soil. The total runoff results from total surface runoff, subsurface flow in the upper soil layer, and free drainage in the soil's lower hydrological column. In this land surface model, the TRIP routing model collects runoff from the cells and transfers water to the basin outlet. The simulation of hydrological components is carried out at the 0.25◦ resolution and daily time scale (Decharme et al. 2010). In this research, surface runoff data of this model have been used.
2.3.2.7. The Ensemble Model
In addition to the above models, the results of all hydrological and land surface models available in the Earth2Observe database are averaged and presented as a new product called Ensemble or Earth2Observe. In this research, evapotranspiration and surface runoff data (in daily time scale) of this database have been used to calibrate the VIC-3L model.
Table 1
The database used in this research.
Model name
|
Data provider (s)
|
Model type
|
Routing
|
Data used
|
Reference(s)
|
ASCAT
|
IPF/TU Wien
|
-
|
-
|
SSM
|
(Wagner et al. 1999)
|
LISFLOOD
|
Joint Research Centre (JRC)
|
GHM
|
Double kinematic
wave
|
SSM
|
(van der Knijff, Younis, and de Roo 2008)
|
HTESSEL
|
European Centre for Medium-Range
Weather Forecasts (ECMWF)
|
LSM
|
CaMa-Flood
|
SSM
|
(Balsamo et al. 2009)
|
HBV
|
Swedish Meteorological and Hydrological Institute
|
GHM
|
Muskingum
|
SSM & ET
|
(Bergström 1972)
|
W3RA
|
Australian National University (ANU) and Commonwealth Scientific and Industrial Research Organization (CSIRO)
|
GHM
|
Cascading linear
reservoirs
|
ET
|
(Van Dijk 2010)
|
GLEAM
|
VU Amsterdam in collaboration with Ghent University in Belgium and the European Space Agency
|
-
|
-
|
ET
|
(Miralles et al. 2011)
|
SURFEX-TRIP
|
Météo France
|
LSM
|
TRIP with stream
|
Surface Runoff
|
(Decharme et al. 2010)
|
2.3. Rescale Surface Soil Moisture
Usually, due to data characteristics or differences in layer depth, systematic biases, and differences between the values of the surface moisture of the observed (SSM of ASCAT, HTESSEL, HBV, and LISFLOOD) and modeled (SSM of VIC-3L) may exist. Several methods can be used to overcome these errors and match the observed and modeled values, such as the cumulative distribution function (CDF) or a mean-standard deviation (Draper et al. 2009; Reichle and Koster 2004). It was decided in order to overcome these errors and match the SSM of the VIC-3L model and the remote sensing sources, the mean-standard deviation be used. The amount of corrected SSM in this way is calculated based on Eq. 1(Draper et al. 2009).
(1)
\(\theta _{{sim}}^{'}\): the adjusted simulated surface soil moisture value, \({\theta _{sim}}\): the simulated (VIC-3L) soil moisture values, \({\bar {\theta }_{sim}}\)and\({\bar {\theta }_{obs}}\): the means of the simulated and observed (ASCAT, HTESSEL, HBV, and LISFLOOD) soil moisture values, \({\sigma _{{\theta _{sim}}}}\)and\({\sigma _{{\theta _{obs}}}}\): the standard deviations of the simulated and observed soil moisture values.
2.4. Calibration Strategy
Despite the different approaches, the Generalized Likelihood Uncertainty Estimation (GLUE) probabilistic method that had a much simpler operational process, especially in the Linux operating system and was first proposed by (Beven and Binley 1992), was used. This method is a combination of the Monte Carlo thought and the Bayesian paradigm, which has high efficiency in calibration and analysis of uncertainty. The fundamental idea behind the Monte Carlo technique is choosing sample points randomly over the desired range to approximate their results (Goodarzi, Ziaei, and Teang Shui 2013). Independ on the calibration parameters is one of the GLUE's hypotheses. The parameters that are not correlated with each other are selected for calibration using this method. The steps for implementing the GLUE method are as follows:
1) Defining the initial domain of calibration parameters, The calibration parameters and their range of variations are mentioned in Sect. 1–3.
2) Generate 1000 random variables for each parameter with uniform distribution (due to the lack of a prior distribution of a parameter).
3) Use the set of random parameters, the VIC-3L hydrological model was implemented on a daily scale.
4) Calculate the likelihood index such as Nash-Sutcliffe efficiency (Beven and Freer 2001) according to the simulated and observational values in each cell of the basin computing network.
5) Opt the most proper set of parameters with the highest Nash-Sutcliffe value (higher than 0.6) as the optimal value of each cell's calibrated parameters.
6) Optimal surface and sub-surface runoff from each grid cell's of basin routed to the outlet.
2.5. Performance Metrics
The performance of the VIC-3L hydrological model in simulating basin river discharge has been evaluated based on indexes, including Pearson Correlation Coefficient (CC), Nash-Sutcliffe (NS), Kling-Gupta (KGE), and Root Mean Square Error (RMSE), in this research. CC is a dimensionless index that is used for investigating the relationship between two variables. It has a value between + 1 (perfectly positively linearly related) and − 1 (entirely negatively linearly related) and is equal to zero if two variables do not have any linear relation whatsoever. NS is a normalized statistic to determine the residual variance's relative magnitude compared to the observed data variance. NS = 1 means a perfect match of the model to the observed data. RMSE measures how much error there is between two data sets. It is close to zero indicates a more accurate model in the simulation of river discharge. Kling-Gupta is an index based on the NS index and the Mean Squared Error (MSE), which facilitates the analysis of the relative importance of its different components (correlation, bias, and variability) in the context of hydrological modeling. These statistics are calculated in Table 2.
Table 2
The statistical indexes which are used in this study
Statistical Index
|
Formula
|
Reference(s)
|
Correlation Coefficient
|
\({\text{CC}}=\frac{{\sum\limits_{{i=1}}^{n} {\left( {{O_i} - \bar {O}} \right)} \left( {{S_i} - \bar {S}} \right)}}{{\sqrt {\sum\limits_{{i=1}}^{n} {{{\left( {{O_i} - \bar {O}} \right)}^2}} \sum\limits_{{i=1}}^{n} {{{\left( {{S_i} - \bar {S}} \right)}^2}} } }}\)
|
(Pearson 1896)
|
Nash-Sutcliffe
|
\(NS=1 - \frac{{\sum\limits_{{i=1}}^{n} {\mathop {(\mathop O\nolimits_{i} - \mathop S\nolimits_{i} )}\nolimits^{2} } }}{{\sum\limits_{{i=1}}^{n} {\mathop {(\mathop O\nolimits_{i} - \overline {O} )}\nolimits^{2} } }}\)
|
(Nash and Sutcliffe 1970)
|
Root Mean Square Error
|
\(RMSE=\sqrt {\frac{1}{n}\sum\limits_{{i=1}}^{n} {\mathop {(\mathop O\nolimits_{i} - \mathop S\nolimits_{i} )}\nolimits^{2} } }\)
|
(Hyndman and Koehler 2006)
|
Kling-Gupta
|
\(KGE=1 - \sqrt {{{(\beta - 1)}^2}+{{(\alpha - 1)}^2}+{{(CC - 1)}^2}}\)
|
(Gupta et al. 2009)
|
\(\beta =\frac{{\bar {S}}}{{\overline {O} }}\)
|
\(\alpha =\frac{{{\sigma _{{S_i}}}}}{{{\sigma _{{O_i}}}}}\)
|
Oi and Si are observed and simulated variable, \(\bar {O}\)and \(\bar {S}\) are the mean of observed and simulated data, \({\sigma _{{O_i}}}\)and \({\sigma _{{S_i}}}\)are the standard deviations of the observed and simulated variable; n is the total number of observations.
2.6. Scenarios
The use of river discharge recorded in the basin outlet is one of the hydrological model's standard calibration methods. Moreover, remote sensing data and re-analysis datasets are the calibration methods that calibrate the model's cell by cell. Therefore in the present research, in addition to calibrating the VIC 3L model based on the observed stream flows of SRB outlet, soil moisture, evapotranspiration, and surface runoff during single-objective and multi-objective approaches to cell by cell calibration, this hydrological model were evaluated. These experiments were carried out for five different calibration scenarios: S1) in this scenario, which is also the basic scenario, calibration using river discharge time series observed of Gilvan hydrometric station. S2) calibration using the surface soil moisture time series of each source: ASCAT, HTESSEL, HBV, and LISFLOOD. S3) calibration using W3RA, GLEAM, HBV, and Ensemble evapotranspiration. S4) use SURFEX-TRIP and Ensemble surface runoff time series to calibration. S5) simultaneous calibration using ASCAT surface soil moisture GLEAM evapotranspiration time series.
S1 aims to calibrate the hydrological model using river discharge observations, following the traditional calibration approach. In scenarios 2 to 4, the time series of surface soil moisture, evapotranspiration, and surface runoff obtained from remote sensing and re-analysis Earth2observe datasets were considered observational values for each SRB grid cell. In other words, improvement of model performance in estimating the time series of SSM, ET, or surface runoff for each grid cell is considered the objective function. While in S5, the matching of the ET and SSM simulated time series with GLEAM ET and ASCAT SSM at the same time is the objective function. Calibration S5 represents the multi-objective calibration approach that combines the strengths of calibration S2 and S3. In general, the scenarios evaluated in the present study can be categorized according to Table 3. The flowchart of the research process is presented in Fig. 2.
Table 3
The scenarios which are used in this study
Scenario
|
Objective functions
|
Products
|
1
|
River Discharge
|
in-situ streamflow
|
2
|
Surface Soil Moisture (SSM)
|
ASCAT
|
HBV
|
HTESSEL
|
LISFLOOD
|
3
|
Evapotranspiration (ET)
|
W3RA
|
HBV
|
GLEAM
|
Ensemble
|
4
|
Surface Runoff
|
SURFEX-TRIP
|
Ensemble
|
5
|
SSM & ET
|
ASCAT & GLEAM
|