Abstract

The river situation in a dry or semi-dry area is extremely affected by climate change and precipitation patterns. In this study, the impact of climate alteration on runoff in Kashafrood River Basin (KRB) in Iran was investigated using the Soil and Water Assessment Tool (SWAT) in historical and three future period times. The runoff was studied by MIROC-ESM and GFDL-ESM2G models as the outputs of general circulation models (GCMs) in the Coupled Model Intercomparison Project Phase 5 (CMIP5) by two representative concentration pathway (RCP) scenarios (RCP2.6 and RCP8.5). The DiffeRential Evolution Adaptive Metropolis (DREAM-ZS) was used to calibrate the hydrological model parameters in different sub-basins. Using DREAM-ZS algorithm, realistic values were obtained for the parameters related to runoff simulation in the SWAT model. In this area, results show that runoff in GFDL-ESM2G in both RCPs (2.6 and 8.5) in comparing future periods with the historical period is increased about 232–383% and in MIROC-ESM tends to increase around 87–292%. Furthermore, GFDL-ESM2G compared to MIROC-ESM in RCP2.6 (RCP8.5) in near, intermediate, and far future periods shows that the value of runoff increases 59.6% (23.0%), 100.2% (35.1%), and 42.5% (65.3%), respectively.

INTRODUCTION

The impact of climate change on the availability of water resources has been attracting the attention of researchers throughout the world (Faramarzi et al. 2009; Rodrigues et al. 2014; Zuo et al. 2015). In the context of precipitation frequency changes and climate change, numerous studies indicate that global warming leads to ecosystem degradation and water crisis (Falkenmark & Rockström 2006; Zuo et al. 2015) and exacerbated water scarcity in dry and semi-dry regions. About 600 million people live in areas with less than 500 m3 water per capita, therefore enduring a severe shortage of water (Pereira et al. 2009). Thus, there is a great deal of concern among researchers about the increasing threat of water scarcity not only in the region but also at global scale (Vörösmarty et al. 2010). As a result, in the global warming context, comprehensive assessment of water resources is a vital part of understanding the availability of water and improving water management towards maintainable, effective, and instrumental usage of scarce freshwater resources. Furthermore, climate change has a significant impact on water availability in river basins and there is an interaction between water resources and climate change (Faramarzi et al. 2009). Spatial and temporal variabilities of water resources in basins are severely impacted by climate change (Vörösmarty et al. 2010). Due to the high variation and low rate of precipitation in dry and semi-dry areas, many sectors can suffer direct or indirect effects on the economy and can lead to loss of people (Anyamba et al. 2014) and, most likely, effects on runoff variation (Chiew et al. 2009). At present, general circulation models (GCMs) are the most frequently used models for the projection of different climatic change scenarios. The Intergovernmental Panel on Climate Change (IPCC) gathers and reviews global climate models as a part of the international climate change assessment reports for policymakers and public users (IPCC 1990). Up to now, IPCC has released five different versions of GCM. The IPCC Fifth Assessment Report of the Intergovernmental Panel on Climate Change (AR5) has introduced a new way of developing scenarios. Numerous studies have been reported about runoff assessment using representative concentration pathway (RCP) scenarios of AR5 on future climate conditions at global and regional scales.

Tan et al. (2017) used five GCMs under RCP2.6, RCP4.5, and RCP8.5 for two future periods (2015–2044 and 2045–2074) and concluded an increase of future annual runoff in Kelantan River Basin in Malaysia. Sun & Peng (2019) using ten GCMs of the CMIP5 under RCP4.5 estimated the future annual runoff in the north–south transect of eastern China for a future period (2011–2100). Yang et al. (2019) simulated precipitation and runoff in Luanhe River Basin in north China using GCMs under four RCPs (i.e., RCP2.6, RCP4.5, RCP6.0, and RCP8.5) for a future period (2020–2030). The result showed that runoff may increase between 39 and 58% compared to the observed runoff depth.

Iran, with a dry and semi-dry climate in the Middle East, is recognized to have low precipitation, high temperature, and high evapotranspiration demand. The water leveling is relatively fragile and endangered by water scarcity that can significantly exacerbate because of climate change (IPCC 2007). Arguably, climate change is one of the big challenges facing Iran; it has an especially harmful impact on water resources which creates problems for the environment and the economy, especially the agricultural sectors. There is a need to predict the potential climate change impacts on the interval and amount of precipitation that can be used for supporting and controlling water resources suitably and decreasing the water scarcity problem (Ternik et al. 2013). A semi-distributed hydrological model, Soil and Water Assessment Tool (SWAT), was used to supply more insights into the mechanisms of hydrological procedures and to assess hydrological processes and effects of climate change on water resources spatiotemporally (developed by United States Department of Agriculture (USDA; Arnold et al. 1998). In northeast Iran, the Kashafrood River Basin (KRB) is the longest river and supplies surface water for Mashhad city. Owing to its extreme dry weather, this basin has been suffering from water scarcity and pollution (Afshar et al. 2017). Up to now, water scarcity related to climate change in the KRB has not been well measured in the analysis of climate change and climate policy construction. Thus, the main purpose of this study was to assess the potential climatic changes on future runoff in KRB according to the CMIP5; Taylor et al. 2012).

The correct values of the parameters in the hydrological model (i.e., SWAT) cannot be directly measured (Li et al. 2010) and should be computed through a calibration process (Laloy et al. 2010; Laloy & Vrugt 2012; Leta et al. 2015; Schoups & Vrugt 2010). The calibration process in a modeling scheme plays an important role in extracting the optimal parameter values which accord reasonably with reality (Leta et al. 2016).

Formal and informal Bayesian frameworks have been developed to calibrate and analyze uncertainty in different river basin models. Formal Bayesian approaches apply likelihood with employing the real assumption to measure the total predictive uncertainty and the posterior probability densities of the parameters (Kuczera & Parent 1998; Schoups & Vrugt 2010; Muñoz et al. 2014; Rivera et al. 2015). The Markov Chain Monte Carlo (MCMC) algorithms can be mentioned as a good example of formal Bayesian approaches. DiffeRential Evolution Adaptive Metropolis (DREAM-ZS), as the most commonly used MCMC algorithm, originated from the DREAM (Vrugt et al. 2009a), which is used as a robust and effective sampler (Vrugt et al. 2009b; Laloy & Vrugt 2012; Laloy et al. 2012). The DREAM-ZS relies on the samplings from the past state and produces candidate points for individual chains. DREAM-ZS provides some advantages: (i) posterior sampling requires only three parallel chains, leading to reduced time for burn-in; (ii) requirement of less function assessments than DREAM, leading to suitable limiting dispensation; (iii) comprising a snooker updater which induces jumps beyond the parallel direction updates and boosts the candidate points' variety (Ter Braak & Vrugt 2008). Although many techniques exist to analyze the uncertainty (Makowski et al. 2002; Mantovan & Todini 2006; Yang et al. 2008; Zhang et al. 2015; Samadi et al. 2017), only a small number of studies have used the DREAM-ZS algorithm to analyze and calibrate the uncertainty in different hydrological models (Laloy & Vrugt 2012; Joseph & Guillaume 2013; Surfleet & Tullos 2013; Laloy et al. 2015; Leta et al. 2015; Zahmatkesh et al. 2015; Nourali et al. 2016; Engeland et al. 2016; Liu et al. 2017; Athira et al. 2018; Han & Zheng 2018).

According to the literature review, the assessment of uncertainty analysis on the SWAT model with DREAM-ZS algorithm for predicting runoff with CMIP5 models in multisite calibration has not been reported yet. Therefore, to fill this gap, the current study was conducted to survey the uncertainty prediction capabilities of the DREAM-ZS algorithm to extract posterior ranges of parameters to predict runoff in future time periods.

In this research, the DREAM-ZS algorithm is applied as a formal Bayesian method for uncertainty analysis, calibration, and measurement of SWAT model parameters in the KRB. First, the main goal of this research is to adjust and verify the SWAT model for the KRB as a large mountainous watershed in Iran with a new uncertainty analysis algorithm (DREAM-ZS). Second, the quantity of monthly runoff is calculated at the sub-basin level considering the effect of water resources management activities in the KRB in future times until the end of the 21st century and with new climatic models (GCMs models) of CMIP5 and RCP scenarios.

MATERIALS AND METHODS

The study area

The study area, Kashafrood Basin (KB), is located in Razavi Khorasan Region (RKR), in northeastern Iran. KB is a massive watershed in RKR with an area of 16,870 km2 and is known as a large-scale watershed. The longest river in this area is Kashafrood River (KR) and its total length is about 285 km. Geographically, KB is located approximately between longitudes 58° 15′ E and 61° 13′ E and latitudes 35° 35′ N and 37° 07′ N (Figure 1). The topographical elevation of the basin varies from 390 to 3,302 m above sea level. The average slope is 4.7% and the average elevation is 1,846 m above sea level. Figure 1 shows the topography of KB, and the location of gauging and meteorological stations. Due to mountainous topographical conditions, KB has a cold climate, semi-dry with high evapotranspiration in summer and low annual precipitation. The annual precipitation in this area is about 340 mm with relatively high spatial fluctuations. In addition, the maximum and minimum temperatures in a year are 20.6 and 7.1 °C, respectively (Afshar et al. 2018). The data set used in the study consists of a digital elevation model (DEM), land-use map (see Figure 2) and soil map from IRS satellite imagery, daily temperature and precipitation data and five monthly runoff data of the basin for 1997–2011 (see Figure 1). Torogh and Kardeh dams are located respectively south and north of Mashhad and provide water for agricultural and drinking sectors (Figure 1). The main soil type in the center and north of the basin is a heterogeneous blend of silt, sand, and clay, and in the southern and eastern parts of the basin it is clay and silt (Afshar et al. 2018).

Figure 1

DEM (digital elevation model) map and the location of gauging and meteorological stations in KB (adopted from Afshar et al. 2018).

Figure 1

DEM (digital elevation model) map and the location of gauging and meteorological stations in KB (adopted from Afshar et al. 2018).

Figure 2

Land-use map of KB (adopted from Afshar et al. 2018).

Figure 2

Land-use map of KB (adopted from Afshar et al. 2018).

SWAT model and data sets

SWAT, a semi-distributed conceptual model, was improved in United States Agriculture Department–Agricultural Research Services to set a rate and predict the effect of management design scenarios on groundwater resources, water quality, pollution loading, and soil erosion (Arnold et al. 1998; Narsimlu et al. 2015; Afshar & Hassanzadeh 2017; Afshar et al. 2018). The SWAT model was used at a daily time step for simulating the runoff in this study. To delineate the basin, the DEM of the watershed was used in the SWAT model. Each sub-basin is also divided into unique land-use/slope/soil hydrologic response units (HRUs). A description of each HRU is prepared by superimposing the land use, slope, and soil maps (Neitsch et al. 2011). In this research, KRB was divided into 217 sub-basins and 635 HRUs with 20, 10, and 20% threshold values for land-use, slope, and soil classes, respectively. The hydrometric stations are located at the outlet of the sub-basins and are selected in terms of area and topographic status to obtain a homogeneity between the sub-basins. By examining the overlapping of the maps, we concluded that by defining these thresholds, we can integrate small polygons with an area of less than 25 hectares, depending on the scale of the maps, into their larger adjacent polygons. Afterwards, a multiple slope option with four classes of 0–5, 5–10, 10–15 and >15% was selected for slope discretization (Afshar et al. 2018). Moreover, the modified Soil Conservation Service Curve Number (SCS-CN) was used to evaluate the surface runoff. The potential evapotranspiration (PET) was also calculated by the Hargreaves method. The Muskingum method was employed to route runoff. Flow was simulated from 1997 to 2011, and the first three years were selected as the model warm-up period. Finally, the monthly data set was used for calibration purposes, which related to an 11-year time span (2001–2011).

DREAM-ZS algorithm

The DREAM algorithm has been employed as an effective and robust sampler in numerous studies (Joseph & Guillaume 2013; Surfleet & Tullos 2013; Zahmatkesh et al. 2015; Engeland et al. 2016; Liu et al. 2017; Athira et al. 2018; Han & Zheng 2018) and exhibited high efficacy and merits for calibration and uncertainty analysis of hydrological models. This algorithm exhibits a significant sampling capability in solving the calibration complications, including high-dimensional search problems and multimodal target distributions, and can be linked to parallel computation by conducting several chains. It has been developed based on shuffled complex evolution metropolis algorithm (Vrugt et al. 2003) via differential evolution-Markov chain approach (Ter Braak 2006). The DREAM algorithm searches the posterior distributions concerning the values of parameters and estimates the parameter uncertainty by Bayesian updating. Moreover, the effect of parameters, structure of model uncertainties, and input data are segregated from total uncertainty (Vrugt et al. 2008; Hassanzadeh et al. 2019). DREAM-ZS, established to measure the parameter posterior probability function (pdf), is a new version of the DREAM algorithm. In DREAM-ZS, ‘Z’ and ‘S’ stand respectively for samples of pervious steps and snooker updater (Laloy & Vrugt 2012). Convergence of the chain with posterior distribution is checked by the R index as (Gelman & Rubin 1992): 
formula
(1)
where, g and q are length and number of the Markov chains, respectively; W is average q within sequences' variances, and B/g is variance between q sequence means.
Given that the values of posterior distribution are less than 1.2 for all of the parameters, it turns out to be stationary (Vrugt et al. 2009a, 2009b). Following the convergence, upon reaching sufficiently higher sampling parameters, the final 20% from samples could be selected to evaluate the uncertainty, and statistical values for posterior distribution (variance and average values) are then estimated (Gelman & Rubin 1992; Laloy & Vrugt 2012). Statistical analysis should be done and uncertainties and parameters need to be calculated via an appropriate likelihood function. However, results become unreliable when likelihood function is used with no model error distribution. In this research, a combination of DREAM-ZS and standard least square (SLS) method was employed as a simplified, classic function of likelihood. In SLS, the errors may not be dependent and homoscedastic (with a fixed variance), having normal probability distribution (PD) with the zero as the mean value as (Box & Tiao 1992): 
formula
(2)
where Y and θ are the vectors of model output and unknown model parameter sets with d dimension, respectively. Parameter set (θ) in Bayesian framework can be calibrated by minimization of the errors e(θ), including distinction between corresponding observed output and the model prediction. However, the errors involved in the hydrologic models (data, structure, output, and parameter calibration errors) could be kept independent having PD with a constant variance and zero mean (Schoups & Vrugt 2010). For details about DREAM-ZS see Vrugt et al. (2009a).

In this study, a prior set of DREAM-ZS algorithm has been considered via the calibrated parameter set of calibration and uncertainty procedures (SWAT-CUP) and sequential uncertainty fitting program phase 2 (SUFI-2) process (Afshar et al. 2018). The flowchart of SWAT modeling with the DREAM-ZS sampler used for predictive uncertainty and calculation of parameter values is illustrated in Figure 3.

Figure 3

Computational framework of this research.

Figure 3

Computational framework of this research.

Sensitivity and uncertainty assessment

To reach a better understanding about the effect of sensitive parameters on runoff, the sensitivity analysis was done by SWAT-CUP tool. The dominant variables of SWAT (see Table 1) could be determined via Latin hypercube sampling (LH-OAT) sensitivity analyzing method (Van Griensven et al. 2006). LH-OAT was applied to explain the initial population of parameters. Posterior distribution of parameters (which were applied to generate model outputs) with uniform prior distributions was obtained based on the preselected ranges (Table 1). Therefore, 95% estimation uncertainty (95PPU) limits (including parameters, model framework, etc.) were achieved via computing 97.5 and 2.5% of cumulative distribution with regard to the outputs.

Table 1

Sensitive parameters (adopted from Afshar et al. 2018)

ParametersRange of parameter
ParametersRange of parameter
MinimumMaximumMinimumMaximum
CN2 −0.4 0.4 ESCO 0.01 
GW_DELAY 400 SLSUBBSN 10 150 
ALPHA_BF CH_N2 0.3 
SOL_AWC −0.3 0.3 CH_K2 150 
SOL_K −0.8 0.8 SFTMP −5 
SOL_BD −0.3 0.3 SMTMP −5 
GW_REVAP 0.02 0.2 SMFMN 10 
SHALLST 1,000 TIMP 0.01 
RCHRG_DP SURLAG 24 
EPCO 0.01 PCPMM −0.5 0.5 
ParametersRange of parameter
ParametersRange of parameter
MinimumMaximumMinimumMaximum
CN2 −0.4 0.4 ESCO 0.01 
GW_DELAY 400 SLSUBBSN 10 150 
ALPHA_BF CH_N2 0.3 
SOL_AWC −0.3 0.3 CH_K2 150 
SOL_K −0.8 0.8 SFTMP −5 
SOL_BD −0.3 0.3 SMTMP −5 
GW_REVAP 0.02 0.2 SMFMN 10 
SHALLST 1,000 TIMP 0.01 
RCHRG_DP SURLAG 24 
EPCO 0.01 PCPMM −0.5 0.5 

DREAM-ZS algorithm is evaluated by three indices: P-factor shows the percentage of the observed data bracketed by 95PPU (Vrugt et al. 2009b; Abbaspour 2011; Alazzy et al. 2015), d-factor is the average thickness of 95PPU band divided by standard deviation of the related calculated variable (Abbaspour 2011), and Nash–Sutcliffe efficiency (NSE) represents the goodness-of-fit between the observed and simulated data (Nash & Sutcliffe 1970). P-factor, NSE, and d-factor are in the range of 0–100%, −∞–1, and 0–∞, respectively. In the calibration process, P-factor close to 100% (observations bracket via prediction uncertainty width) and d-factor value of zero (to reach a smaller uncertainty limit) show good agreement between the observed and simulated values (Yang et al. 2008). Also, the model prediction evaluation is facilitated when the NSE values are 0.75–1 and it may be satisfactory if NSE is greater than 0.36 (Moriasi et al. 2007). The parameter ranges can be determined at acceptable values of NSE efficiency, d-factor, and P-factor.

CMIP5 models

By CMIP5 outputs of GCMs models, 61 climatic time periods can be simulated all over the world by climate change models. In this research, about 14 models were selected to be used for selecting the appropriate CMIP5 models which are most adapted to the observational data (Chaturvedi et al. 2012; Xu et al. 2012; Jiang & Tian 2013; Ou et al. 2013; Iqbal & Zahid 2014; Salzmann et al. 2014; Woldemeskel et al. 2016; Afshar et al. 2017). Table 2 summarizes the total information and names of associated institutions of the models. CMIP5, unlike CMIP3, uses the new RCPs based on radiative forcing level by 2,100, shown in Table 3. The four RCP scenarios are RCP2.6, RCP8.5, RCP4.5, and RCP6.0 (Van Vuuren et al. 2011a). RCP2.6 is representative of the low end of the scenario literature in terms of emissions and radiative forcing (Van Vuuren et al. 2011b), while RCP8.5 includes high emission scenarios which cover the entire range of radiative forcing (Rao & Riahi 2006; Riahi et al. 2007, 2011). However, RCP2.6 and RCP8.5 were utilized in this research to investigate the future climate impacts in runoff. The bias corrected spatial downscaling method was chosen to downscale the precipitation data (for more information, see Afshar et al. 2017).

Table 2

GCMs (global climate models) used in this research

Model nameOrganization
HadGEM2-ES MOHC, Met Office Hadley Centre, UK 
GFDL-CM3 NOAA, Geophysical Fluid Dynamics Laboratory, USA 
GFDL-ESM2M 
GFDL-ESM2G 
MIROC5 MIROC, Centre for Climate System Research, Japan 
MIROC-ESM 
MIROC-ESM-CHEM 
IPSL-CM5A-LR IPSL, Institute Pierre Simon Laplace, France 
IPSL-CM5A-MR 
NorESM1-M NCC, Norwegian Climate Centre, Norway 
BCC-CSM1.1 BCC, Beijing Climate Centre, China 
CSIRO-Mk3.6.0 CSIRO-QCCCE, Commonwealth Scientific and Industrial Research Organization, Australia 
CCSM4 NCAR, National Center for Atmospheric Research, USA 
CESM1(CAM5) NSF-DOE-NCAR, Community Earth System Model Contributors 
Model nameOrganization
HadGEM2-ES MOHC, Met Office Hadley Centre, UK 
GFDL-CM3 NOAA, Geophysical Fluid Dynamics Laboratory, USA 
GFDL-ESM2M 
GFDL-ESM2G 
MIROC5 MIROC, Centre for Climate System Research, Japan 
MIROC-ESM 
MIROC-ESM-CHEM 
IPSL-CM5A-LR IPSL, Institute Pierre Simon Laplace, France 
IPSL-CM5A-MR 
NorESM1-M NCC, Norwegian Climate Centre, Norway 
BCC-CSM1.1 BCC, Beijing Climate Centre, China 
CSIRO-Mk3.6.0 CSIRO-QCCCE, Commonwealth Scientific and Industrial Research Organization, Australia 
CCSM4 NCAR, National Center for Atmospheric Research, USA 
CESM1(CAM5) NSF-DOE-NCAR, Community Earth System Model Contributors 
Table 3

Information about representative concentration pathways (RCPs) (Van Vuuren et al. 2011a)

RCPsDescriptions
2.6 Peak in radiative forcing at ∼3:1 Wm2 (∼490 ppm CO2 equivalent) before 2,100 and then decline (the selected pathway declines to 2.6 Wm2 by 2,100) 
4.5 Stabilization without overshoot pathway to 4.5 Wm2 (∼650 ppm CO2 equivalent) at stabilization after 2,100 
6.0 Stabilization without overshoot pathway to 6 Wm2 at (∼850 ppm CO2 equivalent) stabilization after 2,100 
8.5 Rising radiative forcing pathway leading to 8.5 Wm2 (∼1,370 ppm CO2 equivalent) by 2,100 
RCPsDescriptions
2.6 Peak in radiative forcing at ∼3:1 Wm2 (∼490 ppm CO2 equivalent) before 2,100 and then decline (the selected pathway declines to 2.6 Wm2 by 2,100) 
4.5 Stabilization without overshoot pathway to 4.5 Wm2 (∼650 ppm CO2 equivalent) at stabilization after 2,100 
6.0 Stabilization without overshoot pathway to 6 Wm2 at (∼850 ppm CO2 equivalent) stabilization after 2,100 
8.5 Rising radiative forcing pathway leading to 8.5 Wm2 (∼1,370 ppm CO2 equivalent) by 2,100 

The best CMIP5 models with the highest fit to the observed data were selected on the following steps, after downscaling the climatic parameters of the precipitation and temperature: The historical data of four stations which were located on the network nodes (with intervals of 0.5 degrees (longitude and latitude) to each other) around each station were extracted by ArcGIS. Then, the historical data of the station were calculated via inverse distance weight (IDW) method. Thereafter, according to the precipitation data and Thiessen's polygon averaging method, the CMIP5 model was assessed by Nash–Sutcliffe efficiency (NSE), percent of bias (PBIAS), and ratio of the root mean square error to the standard deviation of measured (RSR) criteria and coefficient of determination (R2) (Moriasi et al. 2007). The results showed that GFDL-ESM2G and MIROC-ESM might be chosen as the best climate models for comparing and predicting the effects of future climate change on the runoff component in the KRB. These models showed the highest compliance with observational precipitation data, based on evaluation metrics such as NSE (MIROC-ESM) = 0.95 and NSE (GFDL-ESM2G) = 0.92), PBIAS (MIROC-ESM) = −2.88 and PBIAS (GFDL-ESM2G) = −2.93), R2 (MIROC-ESM) = 0.97 and R2 (GFDL-ESM2G) = 0.94), and RSR (MIROC-ESM) = 0.33 and RSR (GFDL-ESM2G) = 0.37).

RESULTS AND DISCUSSION

Since the observed data are very restricted, only the results related to the calibration data set are illustrated and the results of the validation data set are excluded.

Results of DREAM-ZS based optimization

DREAM-ZS was employed for separating the effects of data, model structure, and parameter set uncertainties from the total predictive uncertainty (Vrugt et al. 2008). To choose candidate points for the individual chains, this classic method uses samples of previous steps. Three parallel chain runs were assigned to the DREAM-ZS algorithm (to compute the posterior distribution) initializing prior parameter set distribution samples for 50,000 runs (Table 4). DREAM-ZS convergence was constantly checked by R statistic for each parameter (see Figure 4) and could be achieved when all the sampling parameters fell under the threshold of 1.2 (see Figure 4, dashed line). In the current research, the algorithm reached the favorite convergence rate after performing a few iterations (about 5,250 runs) for all sensitive parameters.

Table 4

Fitted value of posterior parameters with ranges of parameters

ParameterRange of parameter
Fitted valueParameterRange of parameter
Fitted value
MinimumMaximumMinimumMaximum
CN2 0.197 0.281 0.23 ESCO 0.41 0.501 0.50 
GW_DELAY 132.7 159 133.55 SLSUBBSN 96.5 109.5 99.54 
ALPHA_BF 0.0 0.08 0.042 CH_N2 0.102 0.13 0.12 
SOL_AWC 0.21 0.28 0.28 CH_K2 101.5 114 101.72 
SOL_K 0.15 0.32 −0.16 SFTMP −1.0 −1.81 −1.45 
SOL_BD 0.095 0.167 −0.11 SMTMP −0.95 −1.92 −1.70 
GW_REVAP 0.085 0.115 0.11 SMFMN 7.5 8.95 8.16 
SHALLST 460.8 510.5 485.30 TIMP 0.50 0.70 0.51 
RCHRG_DP 0.229 0.32 0.24 SURLAG 9.60 12.2 10.22 
EPCO 0.29 0.41 0.29 PCPMM −0.177 −0.325 −0.25 
ParameterRange of parameter
Fitted valueParameterRange of parameter
Fitted value
MinimumMaximumMinimumMaximum
CN2 0.197 0.281 0.23 ESCO 0.41 0.501 0.50 
GW_DELAY 132.7 159 133.55 SLSUBBSN 96.5 109.5 99.54 
ALPHA_BF 0.0 0.08 0.042 CH_N2 0.102 0.13 0.12 
SOL_AWC 0.21 0.28 0.28 CH_K2 101.5 114 101.72 
SOL_K 0.15 0.32 −0.16 SFTMP −1.0 −1.81 −1.45 
SOL_BD 0.095 0.167 −0.11 SMTMP −0.95 −1.92 −1.70 
GW_REVAP 0.085 0.115 0.11 SMFMN 7.5 8.95 8.16 
SHALLST 460.8 510.5 485.30 TIMP 0.50 0.70 0.51 
RCHRG_DP 0.229 0.32 0.24 SURLAG 9.60 12.2 10.22 
EPCO 0.29 0.41 0.29 PCPMM −0.177 −0.325 −0.25 
Figure 4

Convergence rate for Gelman–Rubin scale-reduction metric of 20 parameters of SWAT according to the DREAM-ZS method.

Figure 4

Convergence rate for Gelman–Rubin scale-reduction metric of 20 parameters of SWAT according to the DREAM-ZS method.

After these chains converged, the posterior parameter distributions of the SWAT model were produced with the last 20% of the parameters' set (i.e., 10,000 parameters' set). These plots are illustrated by Figure 5. Also, the best values (maximum likelihood) for all parameters are represented by (×) (see Table 5) where the x-axis represents the posterior ranges for all parameters. Results show that some parameters are efficiently optimized by the DREAM-ZS process for a complex model and resemble less uncertainty as a simply identifiable parameter. Furthermore, EPCO, TIMP, and ALPHA_BF values maintain a uniform distribution of shape and show the low uncertainty in the SWAT model. In other parameters, distributions of CH_N2, GW_DELAY, and ESCO are markedly varied to be concentrated around their lower or upper limits. The ranges of other parameters based on the results must be relaxed as well. Some parameters reached similar distributions or very close distributions. For instance, SOL_AWC, SOL_BD, ALPHA_BF, SURLAG, SLSUBBSN, and CN2 had almost similar values of 0.26, −0.14, 0.04, 11.2, 99, and 0.21, respectively.

Table 5

Values of criteria of total uncertainty (U1) and parameter uncertainty (U2)

Hydrometry gaugesU1
U2
NSEPdNSEPd
SARASSHA 0.72 0.91 1.96 0.72 0.13 0.19 
ZIRBAGOL 0.70 0.89 1.48 0.70 0.14 0.22 
GOLJAGHR 0.57 0.98 5.00 0.57 0.14 0.15 
HESDEHB 0.54 0.98 3.43 0.54 0.15 0.16 
KARTIAN 0.56 0.93 2.79 0.56 0.11 0.17 
Hydrometry gaugesU1
U2
NSEPdNSEPd
SARASSHA 0.72 0.91 1.96 0.72 0.13 0.19 
ZIRBAGOL 0.70 0.89 1.48 0.70 0.14 0.22 
GOLJAGHR 0.57 0.98 5.00 0.57 0.14 0.15 
HESDEHB 0.54 0.98 3.43 0.54 0.15 0.16 
KARTIAN 0.56 0.93 2.79 0.56 0.11 0.17 
Figure 5

Distributions of marginal posterior parameters.

Figure 5

Distributions of marginal posterior parameters.

To evaluate the uncertainty of overland flow parameters for five considered hydrometry gauges following the convergence, the last 20% of the SWAT samples which fitted the model well were extracted to produce the output values and 95% confidence interval was defined by computing 97.5 and 2.5% levels of the cumulative distribution of the output variables. The DREAM-ZS algorithm plotted the parameter uncertainty (dark shaded area), the time series data related to total uncertainty (light shaded area), the best simulation obtained by SLS likelihood function (black line), and the monthly flow discharge (filled circles) observed in five gauges from 2001 to 2011, as shown in Figure 6. The model for KRB was calibrated according to the data from one downstream (KARTIAN), one upstream (SARASSHA), and three mid-stream (GOLJAGHR, HESDEHB, and ZIRBAGOL) gauging stations. Also, in Figure 5, the x-axis of plots is demonstrated time from 0 to 132 months. The obtained values for three evaluation measures (NSE, d-factor, and P-factor), which show the performance of uncertainty analysis, are presented in Table 5. The uncertainty origins, such as U1 and U2, might be considered as the total (and parameter) uncertainty.

Figure 6

95PPU ranges of five runoff stations.

Figure 6

95PPU ranges of five runoff stations.

As can be seen in Figure 6, the most observation points (over 90%) of DREAM-ZS take place within 95% predictive uncertainty limits and P-factor demonstrates high performance of the DREAM-ZS algorithm. Also, the results illustrate that the total uncertainty (and parameters' set uncertainty) band covers 89–98% (and 11–34%) of the measurements' data by identical d-factor of 1.36–5 (and from 0.15 to 0.29), respectively. The minimum and maximum P and d-factors are for Gand ZIRBAGOL (1.48 and 89%) and GOLJAGHR (5 and 98%) stations, respectively. Considering only the SWAT parameter uncertainty, 11–15% of the observations are bracketed. On the contrary, limits of the lower and upper parameters' uncertainty involve most of the observation of the band. Also, the gap between observed data and the parameter uncertainty may result from the model structural, input and output insufficiency (Laloy et al. 2010). Although the total prediction uncertainty bounds are quite wide, the parameter uncertainty bounds are quite narrow. On average for all stations, the value of NSE obtained was around 0.62 (Table 5). Also, the NSE values show relatively good agreement between simulated and observed runoff values. Using the DREAM-ZS algorithm, the peak runoff values are mostly underestimated, while in only a few months are they overestimated (especially in the end part of the calibration period) in all runoff stations. The accuracy of the model is good in peak discharge estimations, although not excellent. In wet periods, the rainfall gauge stations are not a very good representative of precipitation over the basin, because the map shows there are no rain gauges in some high and most low altitude parts of the study basin. Also, the KRB is a mountainous watershed and is exposed to high spatial and temporal variability in rainfall distribution that cannot be captured exactly using rain gauges in the basin. Therefore, such limitation could have contributed to some degree of simulation uncertainty (Cho et al. 2009; Memarian et al. 2014). Overall, the results show a good agreement between the simulated and measured data providing a realistic modeling.

Results of runoff predictions for the future via CMIP5 models

The spatial variations of runoff in the KRB during the historical (1997–2011) and future time periods (near: 2014–2042, intermediate: 2043–2071, and far: 2072–2100) (Afshar et al. 2018) in GFDL-ESM2G and MIROC-ESM models under both scenarios (RCP8.5 and RCP2.6) are respectively illustrated in Figures (7)(9).

Figure 7

Spatial variation of runoff in the KRB in historical period.

Figure 7

Spatial variation of runoff in the KRB in historical period.

Figure 8

Spatial variation of runoff in the KRB for both models and under RCP2.6 in future periods (left panels belong to MIROC-ESM, right panels belong to GFDL-ESM2G).

Figure 8

Spatial variation of runoff in the KRB for both models and under RCP2.6 in future periods (left panels belong to MIROC-ESM, right panels belong to GFDL-ESM2G).

Figure 9

Spatial variation of runoff in the KRB for both models and under RCP8.5 in future periods (left panels belong to MIROC-ESM, right panels belong to GFDL-ESM2G).

Figure 9

Spatial variation of runoff in the KRB for both models and under RCP8.5 in future periods (left panels belong to MIROC-ESM, right panels belong to GFDL-ESM2G).

The simulation results show that in the northern part of the watershed, the runoff values would be higher than the southern part with almost the same pattern for different time periods. However, there are some differences in the sub-basin scale in the amount of runoff (Figures (7)–(9)). The mean annual runoff in the near future under RCP2.6 (RCP8.5) in GFDL-ESM2G and MIROC-ESM models, respectively, is about 173 (182) and 108 (147) mm/year. Also, in the intermediate future period, the mean value of annual runoff under RCP2.6 and RCP8.5 in GFDL-ESM2G (MIROC-ESM) is 141 (70) and 160 (119) mm/year, respectively. Simulated annual runoff in the far future period via the GFDL-ESM2G (MIROC-ESM) model under RCP2.6 and RCP8.5 is about 176 (124) and 125 (75) mm/year. The results show that in the near and intermediate periods, by both models under RCP8.5, simulated runoff will be higher than RCP2.6, while in the far future under RCP2.6, the mean annual of the simulated runoff will be higher than that in the other scenario. Furthermore, the average annual runoff obtained by the GFDL-ESM2G model in all future periods and under both scenarios is higher than the other selected model (MIROC-ESM). By comparing the runoff values in historical and future periods it is clear that the hydrology condition of the watershed will be improved due to the increase of runoff in all future periods. More details are presented in Table 6 with regards to comparing runoff in future time periods and via both CMIP5 models with historical time.

Table 6

Comparison of runoff in historical and future time periods in both CMIP5 models (mm/year)

Period timesGFDL-ESM2G
MIROC-ESM
RCP2.6RCP8.5RCP2.6RCP8.5
Near future 172.72 181.53 108.19 147.50 
Mid future 140.65 160.42 70.24 118.71 
Late future 176.51 124.74 123.86 75.44 
Historical 37.54 
Period timesGFDL-ESM2G
MIROC-ESM
RCP2.6RCP8.5RCP2.6RCP8.5
Near future 172.72 181.53 108.19 147.50 
Mid future 140.65 160.42 70.24 118.71 
Late future 176.51 124.74 123.86 75.44 
Historical 37.54 

The minimum, maximum, and mean annual runoff values are shown in Tables 7 and 8, respectively, for upstream (in the higher part of a stream) and downstream (the direction that a river flows) sub-basins in the KRB. In general, sub-basins located in the upstream show higher runoff value, while downstream sub-basins show smaller annual runoff values under both scenarios and models. The average annual runoff values for upstream in the MIROC-ESM model under RCP2.6 for near, intermediate, and far future are about 138.8, 94.8, and 163.6 mm/year, while in the downstream it decreases to around 43.5, 51.5, and 47.9%, respectively. In the mentioned model under RCP8.5 scenario almost the same pattern is seen in the future periods and in downstream runoff the three future periods show 41.4, 38.2, and 45.3% less values than those for the upstream. In the MIROC-ESM model under both scenarios and future periods, the range of mean annual runoff variation in upstream and downstream is between 38.2 and 51.1%, while in the GFDL-ESM2G model the range of this variation is between 14.0 and 44.0%. The lowest difference between upstream and downstream will happen in the far future under RCP8.5 in the GFDL-ESM2G model and under RCP2.6 (RCP8.5), the variation in the near, intermediate, and far future is about 40.9 (44.0), 41.0 (30.2), and 38.7 (14.0)%, respectively. On the other hand, runoff declined about 60.7% from the upstream to downstream in the historical time period.

Table 7

Variations of runoff for upstream and downstream of the KRB in the MIROC-ESM model

RangeHistoricalRCP2.6
RCP8.5
Near futureMid futureLate futureNear futureMid futureLate future
 Upstream sub-basins of KRB 
Min 1.57 46.39 23.87 48.84 73.13 38.08 23.21 
Max 125.95 249.20 179.78 278.57 353.00 467.46 206.22 
Mean 54.22 138.84 94.84 163.58 186.68 147.26 97.91 
 Downstream sub-basins of KRB 
Min 0.70 6.63 4.78 9.54 15.48 10.11 7.46 
Max 96.23 215.18 155.08 254.13 277.05 465.54 203.75 
Mean 21.32 78.38 46.32 85.23 109.40 90.94 53.59 
RangeHistoricalRCP2.6
RCP8.5
Near futureMid futureLate futureNear futureMid futureLate future
 Upstream sub-basins of KRB 
Min 1.57 46.39 23.87 48.84 73.13 38.08 23.21 
Max 125.95 249.20 179.78 278.57 353.00 467.46 206.22 
Mean 54.22 138.84 94.84 163.58 186.68 147.26 97.91 
 Downstream sub-basins of KRB 
Min 0.70 6.63 4.78 9.54 15.48 10.11 7.46 
Max 96.23 215.18 155.08 254.13 277.05 465.54 203.75 
Mean 21.32 78.38 46.32 85.23 109.40 90.94 53.59 
Table 8

Variations of runoff for upstream and downstream of the KRB in the GFDL-ESM2G model

RangeHistoricalRCP2.6
RCP8.5
Near futureMid futureLate futureNear futureMid futureLate future
 Upstream sub-basins of KRB 
Min 1.57 74.61 60.56 83.19 74.67 71.52 49.78 
Max 125.95 397.13 330.62 403.39 434.39 335.27 213.59 
Mean 54.22 217.99 177.61 219.55 233.74 189.47 134.29 
 Downstream sub-basins of KRB 
Min 0.70 47.94 19.94 53.87 36.12 41.31 27.96 
Max 96.23 281.34 236.77 290.47 276.64 268.97 239.37 
Mean 21.32 128.68 104.70 134.64 130.75 132.16 115.46 
RangeHistoricalRCP2.6
RCP8.5
Near futureMid futureLate futureNear futureMid futureLate future
 Upstream sub-basins of KRB 
Min 1.57 74.61 60.56 83.19 74.67 71.52 49.78 
Max 125.95 397.13 330.62 403.39 434.39 335.27 213.59 
Mean 54.22 217.99 177.61 219.55 233.74 189.47 134.29 
 Downstream sub-basins of KRB 
Min 0.70 47.94 19.94 53.87 36.12 41.31 27.96 
Max 96.23 281.34 236.77 290.47 276.64 268.97 239.37 
Mean 21.32 128.68 104.70 134.64 130.75 132.16 115.46 

On a sub-basin scale (Figure 10), the highest rate of runoff in both models, RCPs, and in all future periods occurred in sub-basin number 18 and also under RCP8.5 respectively happened in sub-basin number 100 (MIROC-ESM in the intermediate and far future) and sub-basin 170 (GFDL-ESM2G model in far future). On the other hand, in the GFDL-ESM2G model under both scenarios and all future periods, the lowest value of runoff was seen in sub-basin number 190, while in the MIROC-ESM model in all future periods and RCPs sub-basin number 203 (except in far future and RCP2.6 in sub-basin 216), experienced the minimum value of this component.

Figure 10

Sub-basins' number map of KRB and selected sub-basin for runoff component.

Figure 10

Sub-basins' number map of KRB and selected sub-basin for runoff component.

CONCLUSIONS

In this study, formal Bayesian frameworks (DREAM-ZS) were investigated to simulate runoff discharge over a semi-dry basin system (with an area of 16,870 km2) in Iran by the SWAT model. Twenty parameters out of initial parameters were determined as sensitive parameters to assess uncertainty and calibration processes at five gauging stations. The DREAM-ZS algorithm was used to optimize the parameters by taking SLS as a likelihood function. After 50,000 runs, the best posterior parameters were obtained. Investigation of marginal density implied that both parameter ranges and shape had almost a uniform posterior distribution for some parameters. The distribution of ALPHA_BF, CN2, and CH_N2 could be well-defined, whereas the remaining parameters were either positively or negatively skewed. Results of calibration and uncertainty analysis indicated that at all stations the peak runoff in most months predicts underestimation via the calibrated model, while the falling limb had overestimated values in some months. The base flow had a high correspondence between simulated and observation data for all calibration stations. For the total uncertainty, NSE was obtained in the range of 0.54–0.72, most data were bracketed within 95PPU bounds (89–98%), and various d-factors (1.48–5.00) were noticed. In DREAM-ZS, the peak runoff points were slightly overestimated and underestimated as compared to the observed data in some months. Therefore, uncertainty analysis approved the high efficiency of DREAM-ZS in predicting parameters uncertainty via the KRB as a future research plan. The mean annual runoff in the MIROC-ESM model under RCP2.6 (RCP8.5) from near to far future periods, respectively, was 108 (147), 70 (119), and 124 (75) mm/year while in the other selected model the results of the simulation were 173 (181), 141 (160), and 176 (125) mm/year. Comparison of GFDL-ESM2G and MIROC-ESM models for near, intermediate, and far future periods under RCP2.6 (RCP8.5) showed the mean annual runoff in GFDL-ESM2G was about 59.6% (23.0%), 100.2% (35.1%), and 42.5% (65.3%) more than the other model. The results of this research can provide reliable information about runoff at sub-basin scale for policymakers.

A follow-up work to this study would be an investigation on water resource security (demand, availability, scarcity, reliability, and vulnerability) based on the blue and green water concepts by applying the mentioned models in this region.

REFERENCES

REFERENCES
Abbaspour
K. C.
2011
SWAT-CUP4: SWAT Calibration and Uncertainty Programs–A User Manual
.
Swiss Federal Institute of Aquatic Science and Technology (Eawag)
,
Dübendorf
,
Switzerland
.
Afshar
A. A.
&
Hassanzadeh
Y.
2017
Determination of monthly hydrological erosion severity and runoff in torogh dam watershed basin using SWAT and WEPP models
.
Iranian Journal of Science and Technology, Transactions of Civil Engineering
41
(
2
),
221
228
.
Afshar
A. A.
,
Hasanzadeh
Y.
,
Besalatpour
A. A.
&
Pourreza-Bilondi
M.
2017
Climate change forecasting in a mountainous data scarce watershed using CMIP5 models under representative concentration pathways
.
Theoretical and Applied Climatology
129
(
1–2
),
683
699
.
Anyamba
A.
,
Small
J. L.
,
Britch
S. C.
,
Tucker
C. J.
,
Pak
E. W.
,
Reynolds
C. A.
,
Crutchfield
J.
&
Linthicum
K. J.
2014
Recent weather extremes and impacts on agricultural production and vector-borne disease outbreak patterns
.
PloS One
9
(
3
),
e92538
.
Arnold
J. G.
,
Srinivasan
R.
,
Muttiah
R. S.
&
Williams
J. R.
1998
Large area hydrologic modeling and assessment part I: model development
.
Journal of the American Water Resources Association
34
(
1
),
73
89
.
Athira
P.
,
Nanda
C.
&
Sudheer
K. P.
2018
A computationally efficient method for uncertainty analysis of SWAT model simulations
.
Stochastic Environmental Research and Risk Assessment
32
(
6
),
1479
1492
.
Box
G. E. P.
&
Tiao
G. C.
1992
Bayesian Inference in Statistical Analysis
.
Wiley Interscience
,
New York
,
USA
.
Chaturvedi
R. K.
,
Joshi
J.
,
Jayaraman
M.
,
Bala
G.
&
Ravindranath
N. H.
2012
Multi-model climate change projections for India under representative concentration pathways
.
Current Science
103
(
7
),
791
802
.
Chiew
F. H. S.
,
Teng
J.
,
Vaze
J.
,
Post
D. A.
,
Perraud
J. M.
,
Kirono
D. G. C.
&
Viney
N. R.
2009
Estimating climate change impact on runoff across southeast Australia: method, results, and implications of the modeling method
.
Water Resources Research
45
(
10
),
1
17
.
Cho
J.
,
Bosch
D.
,
Lowrance
R.
,
Strickland
T.
&
Vellidis
G.
2009
Effect of spatial distribution of rainfall on temporal and spatial uncertainty of SWAT output
.
Transactions of the ASABE
52
(
5
),
1545
1556
.
Engeland
K.
,
Steinsland
I.
,
Johansen
S. S.
,
Petersen-Øverleir
A.
&
Kolberg
S.
2016
Effects of uncertainties in hydrological modelling. A case study of a mountainous catchment in Southern Norway
.
Journal of Hydrology
536
,
147
160
.
Falkenmark
M.
&
Rockström
J.
2006
The new blue and green water paradigm: breaking new ground for water resources planning and management
.
Journal of Water Resources Planning and Management
132
,
129
132
.
Faramarzi
M.
,
Abbaspour
K. C.
,
Schulin
R.
&
Yang
H.
2009
Modeling blue and green water resources availability in Iran
.
Hydrological Processes
23
,
486
501
.
Gelman
A.
&
Rubin
D. B.
1992
Inference from iterative simulation using multiple sequences
.
Statistical Science
7
,
457
472
.
Hassanzadeh
Y.
,
Afshar
A. A.
,
Pourreza-Bilondi
M.
,
Memarian
H.
&
Besalatpour
A. A.
2019
Toward a combined Bayesian frameworks to quantify parameter uncertainty in a large mountainous catchment with high spatial variability
.
Environmental Monitoring and Assessment
191
(
1
),
23
.
IPCC
1990
Climate Change: The IPCC Scientific Assessment (1990)
.
Cambridge University Press
,
London
,
UK
.
IPCC
2007
Climate Change 2007: Impacts, Adaptation, and Vulnerability
.
Exit EPA Disclaimer Contribution of Working Group II to the Third Assessment Report of the Intergovernmental Panel on Climate Change
.
Cambridge University Press
,
Cambridge
,
UK
.
Iqbal
W.
&
Zahid
M.
2014
Historical and future trends of summer mean air temperature over South Asia
.
Pakistan Journal of Meteorology
10
(
20
),
67
74
.
Joseph
J. F.
&
Guillaume
J. H.
2013
Using a parallelized MCMC algorithm in R to identify appropriate likelihood functions for SWAT
.
Environmental Modelling & Software
46
,
292
298
.
Leta
O. T.
,
Nossent
J.
,
Velez
C.
,
Shrestha
N. K.
,
Van Griensven
A.
&
Bauwens
W.
2015
Assessment of the different sources of uncertainty in a SWAT model of the River Senne (Belgium)
.
Environmental Modelling & Software
68
,
129
146
.
Li
X.
,
Weller
D. E.
&
Jordan
T. E.
2010
Watershed model calibration using multi-objective optimization and multi-site averaging
.
Journal of Hydrology
380
(
3–4
),
277
288
.
Liu
Y. R.
,
Li
Y. P.
,
Huang
G. H.
,
Zhang
J. L.
&
Fan
Y. R.
2017
A Bayesian-based multilevel factorial analysis method for analyzing parameter uncertainty of hydrological model
.
Journal of Hydrology
553
,
750
762
.
Memarian
H.
,
Balasundram
S. K.
,
Abbaspour
K. C.
,
Talib
J. B.
,
Boon Sung
C. T.
&
Sood
A. M.
2014
SWAT-based hydrological modelling of tropical land-use scenarios
.
Hydrological Sciences Journal
59
(
10
),
1808
1829
.
Moriasi
D. N.
,
Arnold
J. G.
,
Van Liew
M. W.
,
Bingner
R. L.
,
Harmel
R. D.
&
Veith
T. L.
2007
Model evaluation guidelines for systematic quantification of accuracy in watershed simulations
.
Transactions of the ASABE
50
(
3
),
885
900
.
Muñoz
E.
,
Rivera
D.
,
Vergara
F.
,
Tume
P.
&
Arumí
J. L.
2014
Identifiability analysis: towards constrained equifinality and reduced uncertainty in a conceptual model
.
Hydrological Sciences Journal
59
(
9
),
1690
1703
.
Narsimlu
B.
,
Gosain
A. K.
,
Chahar
B. R.
,
Singh
S. K.
&
Srivastava
P. K.
2015
SWAT model calibration and uncertainty analysis for streamflow prediction in the Kunwari River Basin, India, using sequential uncertainty fitting
.
Environmental Processes
2
(
1
),
79
95
.
Neitsch
S. L.
,
Arnold
J. G.
,
Kiniry
J. R.
&
Williams
J. R.
2011
Soil and Water Assessment Tool Theoretical Documentation Version 2009
.
Texas Water Resources Institute
,
College Station, TX
,
USA
.
Nourali
M.
,
Ghahraman
B.
,
Pourreza-Bilondi
M.
&
Davary
K.
2016
Effect of formal and informal likelihood functions on uncertainty assessment in a single event rainfall-runoff model
.
Journal of Hydrology
540
,
549
564
.
Ou
T.
,
Chen
D.
,
Linderholm
H. W.
&
Jeong
J. H.
2013
Evaluation of global climate models in simulating extreme precipitation in China
.
Tellus A: Dynamic Meteorology and Oceanography
65
(
1
),
19799
.
Pereira
L. S.
,
Cordery
I.
&
Iacovides
I.
2009
Coping With Water Scarcity: Addressing the Challenges
.
Springer Science & Business Media
,
Dordrecht, Netherlands
, p.
272
.
Rao
S.
&
Riahi
K.
2006
The role of non-CO2 greenhouse gases in climate change mitigation: long-term scenarios for the 21st century
.
The Energy Journal
27
(
3
),
177
200
.
Riahi
K.
,
Gruebler
A.
&
Nakicenovic
N.
2007
Scenarios of long-term socio-economic and environmental development under climate stabilization
.
Technological Forecasting and Social Change
74
(
7
),
887
935
.
Riahi
K.
,
Rao
S.
,
Krey
V.
,
Cho
C.
,
Chirkov
V.
,
Fischer
G.
,
Kindermann
G.
,
Nakicenovic
N.
&
Rafaj
P.
2011
RCP 8.5—a scenario of comparatively high greenhouse gas emissions
.
Climatic Change
109
(
1–2
),
33
57
.
Rodrigues
D. B.
,
Gupta
H. V.
&
Mendiondo
E. M.
2014
A blue/green water based accounting framework for assessment of water security
.
Water Resources Research
50
,
7187
7205
.
Salzmann
M.
,
Weser
H.
&
Cherian
R.
2014
Robust response of Asian summer monsoon to anthropogenic aerosols in CMIP5 models
.
Journal of Geophysical Research: Atmospheres
119
(
19
),
11321
11337
.
Samadi
S.
,
Tufford
D. L.
&
Carbone
G. J.
2017
Assessing parameter uncertainty of a semi-distributed hydrology model for a shallow aquifer dominated environmental system
.
Journal of the American Water Resources Association
53
(
6
),
1368
1389
.
Sun
G.
&
Peng
F.
2019
Evaluation of future runoff variations in north–south transect of eastern China: effects of CMIP5 models outputs uncertainty
.
Journal of Water and Climate Change
(in press)
,
1
15
.
Taylor
K. E.
,
Stouffer
R. J.
&
Meehl
G. A.
2012
An overview of CMIP5 and the experiment design
.
Bulletin of the American Meteorological Society
93
(
4
),
485
498
.
Ter Braak
C. J.
&
Vrugt
J. A.
2008
Differential evolution Markov chain with snooker updater and fewer chains
.
Statistics and Computing
18
(
4
),
435
446
.
Terink
W.
,
Immerzeel
W. W.
&
Droogers
P.
2013
Climate change projections of precipitation and reference evapotranspiration for the Middle East and Northern Africa until 2050
.
International Journal of Climatology
33
(
14
),
3055
3072
.
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
(
1
),
10
23
.
Van Vuuren
D. P.
,
Edmonds
J.
,
Kainuma
M.
,
Riahi
K.
,
Thomson
A.
,
Hibbard
K.
,
Hurtt
G. C.
,
Kram
T.
,
Krey
V.
,
Lamarque
J. F.
,
Masui
T.
,
Meinshausen
M.
,
Nakicenovic
N.
,
Smith
S. J.
&
Rose
S. K.
2011a
The representative concentration pathways: an overview
.
Climatic Change
109
,
5
31
.
Van Vuuren
D. P.
,
Stehfest
E.
,
den Elzen
M. G.
,
Kram
T.
,
van Vliet
J.
,
Deetman
S.
,
Isaac
M.
,
Goldewijck
K. K.
,
Hof
A.
,
Mendoza Beltran
A.
,
Oostenrijck
A.
&
van Ruijven
B.
2011b
RCP2. 6: exploring the possibility to keep global mean temperature increase below 2 °C
.
Climatic Change
109
(
1–2
),
95
.
Vörösmarty
C. J.
,
Mclntyre
P. B.
,
Gessner
M. O.
&
Duegeon
D.
2010
Global threats to human water security and river biodiversity
.
Nature
467
,
555
561
.
Vrugt
J. A.
,
Gupta
H. V.
,
Bouten
W.
&
Sorooshian
S.
2003
A shuffled complex evolution metropolis algorithm for optimization and uncertainty assessment of hydrologic model parameters
.
Water Resources Research
39
(
8
),
1
18
.
Vrugt
J. A.
,
Ter Braak
C. J.
,
Clark
M. P.
,
Hyman
J. M.
&
Robinson
B. A.
2008
Treatment of input uncertainty in hydrologic modeling: doing hydrology backward with Markov chain Monte Carlo simulation
.
Water Resources Research
44
(
12
),
1
15
.
Vrugt
J. A.
,
Ter Braak
C. J. F.
,
Diks
C. G. H.
,
Robinson
B. A.
,
Hyman
J. M.
&
Higdon
D.
2009a
Accelerating Markov chain Monte Carlo simulation by differential evolution with self-adaptive randomized subspace sampling
.
International Journal of Nonlinear Sciences and Numerical Simulation
10
(
3
),
273
290
.
Vrugt
J. A.
,
Ter Braak
C. J. F.
,
Gupta
H. V.
&
Robinson
B. A.
2009b
Equifinality of formal (DREAM) and informal (GLUE) Bayesian approaches in hydrologic modeling
.
Stochastic Environmental Research and Risk Assessment
23
(
7
),
1011
1026
.
Woldemeskel
F. M.
,
Sharma
A.
,
Sivakumar
B.
&
Mehrotra
R.
2016
Quantification of precipitation and temperature uncertainties simulated by CMIP3 and CMIP5 models
.
Journal of Geophysical Research: Atmospheres
121
(
1
),
3
17
.
Xu
L.
,
Xie
S. P.
&
Liu
Q.
2012
Mode water ventilation and subtropical countercurrent over the North Pacific in CMIP5 simulations and future projections
.
Journal of Geophysical Research: Oceans
117
(
C12
),
1
12
.
Yang
J.
,
Reichert
P.
,
Abbaspour
K. C.
,
Xia
J.
&
Yang
H.
2008
Comparing uncertainty analysis techniques for a SWAT application to the Chaohe Basin in China
.
Journal of Hydrology
358
(
1
),
1
23
.
Zhang
J.
,
Li
Q.
,
Guo
B.
&
Gong
H.
2015
The comparative study of multi-site uncertainty evaluation method based on SWAT model
.
Hydrological Processes
29
(
13
),
2994
3009
.
Zuo
D.
,
Xu
Z.
,
Peng
D.
,
Song
J.
,
Cheng
L.
,
Wei
S.
,
Abbaspour
K. C.
&
Yang
H.
2015
Simulating spatiotemporal variability of blue and green water resources availability with uncertainty analysis
.
Hydrological Processes
29
(
8
),
1942
1955
.