Abstract

The variation in land use/land cover (LULC) and climate have a direct impact on the accuracy of any hydrological prediction. However, quantification of the effect of these two factors in an ungauged catchment setting is less discussed. Soil and Water Assessment Tool (SWAT) in combination with two regionalization techniques, i.e., Inverse Distance Weighted (IDW) and Kriging were applied on 32 catchments in India where each catchment was considered as ungauged at least once. The combined and isolated impacts of LULC change (LULCC) and climate variability on streamflow for the period of 1990–2011 were quantified at an annual scale through four different cases. Satisfactory results were obtained from SWAT for the analysis of both the gauged and ungauged set-up. The overall outcomes suggest that, due to the influence of the combined effects of LULCC and climate variability, there was a decrease in the annual streamflow volume by more than 21% from the first period (1990–2000) to the second period (2001–2011) in the selected catchment treated as ungauged. The variable climate factor overshadowed the effect of LULCC. The result may be correlated with the increase in temperature and the decrease in rainfall volume, which is distinctive in a monsoon-dominated country like India.

INTRODUCTION

Accurate streamflow measurement records are very much necessary for various events of water resources planning and management (Blöschl 2016). One way of getting such streamflow records is to measure the streamflow on a particular gauging site for a long period. However, in most cases, this is not feasible due to a variety of reasons such as financial, logistics or purely due to lack of intent (Blöschl 2016). This leads to the fact that a large portion of the total number of river basins present around the globe are either poorly gauged or completely ungauged (Sivapalan et al. 2003; Young 2006; Blöschl et al. 2013). In addition, the hydrological parameter observation networks are continuously weakening (Mishra & Coulibaly 2009). Furthermore, the number of gauging stations in comparison to the large river basins are very low. Hence, small watersheds are devoid of water resources assessment data (Ibrahim et al. 2015).

Prediction of hydro-climatological variables for a drainage basin requires modelling, where hydrological models simplify the complex hydrological processes occurring within the basin and present them in a simplified manner. Estimation of streamflow in gauged and ungauged catchments is being carried out using either distributed physically-based models, conceptual and semi-distributed models or data-driven models (Razavi & Coulibaly 2012). This is only possible where long-term hydrological records are available, i.e., gauged catchments. In the case of ungauged catchments, regionalization is the method for estimation of different variables like streamflow, sediment load and water quality parameters. The term regionalization can be explained as follows: it is the transfer of an optimized model parameter set from one or more gauged (donor) catchments to an ungauged (receiver) catchment (Blöschl & Sivapalan 1995; Oudin et al. 2010). Regionalization is categorized into two types: hydrologic model-dependent and hydrologic model-independent (Razavi & Coulibaly 2012). Spatial proximity, physical similarity and regression analysis are some of the widely used regionalization techniques across the globe in both categories (Nathan & McMahon 1990; Sefton & Howarth 1998; Yu et al. 2002; Merz & Blöschl 2004; Parajka et al. 2005; Götzinger & Bárdossy 2007; Oudin et al. 2008; Samuel et al. 2011; Zhang et al. 2015; Müller & Thompson 2016; Swain & Patra 2017a, 2017b).

Considerable progress has been observed in the PUB (Prediction in Ungauged Basins) initiative during the last decade. However, most of the studies were case studies conducted in well-gauged catchments situated in well-developed countries (Visessri & McIntyre 2015). On the other hand, very few studies address the issue of the non-stationary environment during regionalization. Specifically, in developing countries, where the land surface is undergoing rapid change, the issue of non-stationarity holds great importance. Land use/land cover change (LULCC) and climate variability are important elements out of the many active environmental factors affecting catchment hydrology (Chen et al. 2012; Zhang et al. 2016). The reason for this may be that land use change affects evapotranspiration (ET), interception and infiltration. Thus, ultimately surface and sub-surface flows are affected (Niraula et al. 2015). Whereas climate variability brings distinctive changes in hydrological regimes altering spatial and temporal patterns of water resources in an area (Khoi & Suetsugi 2014). A few cases dealing with the non-stationarity during regionalization are discussed below.

The process of predicting the hydrologic response and streamflow response to land use change provided good results [Nash-Sutcliffe efficiency (NSE) ≤ 0.85] when Croke et al. (2004) applied the ‘Identification of unit Hydrographs And Component flows from Rainfall, Evaporation and Streamflow data’ (IHACRES) and CATCHCROP model to the ungauged catchments in Thailand using a simple regionalization and scaling technique. Singh et al. (2011) introduced the idea of trading-space-for-time to counter the changing climate during continuous streamflow prediction and was successful in implementing it during the study of 394 US watersheds. The transferability of hydrological model parameters under non-stationary conditions was investigated by Li et al. (2012) using differential split-sample tests and Monte Carlo simulation in combination with Dynamic Water Balance Model and SIMHYD models in 30 catchments in Australia. Transferring model parameter values obtained from dry periods to wet periods resulted in smaller errors in streamflow prediction. Development of a multiple regression model using few catchment characteristics led to efficient estimation of the non-stationary regional flood quantiles [Q5(t) and Q100(t)], with relative root mean square error (RMSEr) of 38.2% and 60.8% respectively during the application of a non-stationary regional model in 29 catchments of Canada and Northern USA by Leclerc & Ouarda (2007). Lima & Lall (2010) applied a hierarchical Bayesian model using the scaling relationship between catchment properties which resulted in a competent estimation of monthly flood flow probability distribution parameters for 44 ungauged sites in Brazil in non-stationary settings. Singh et al. (2014) tested the combination of streamflow signature regionalization with the trading of space for time approach for robust hydrologic projections in the Olifants basin in southern Africa under changing climate condition. The statistically downscaled General Circulation Model against observed climate data led to minor deterioration in reliability of runoff projections for the historical period. Visessri & McIntyre (2015) demonstrated the use of the regression method against catchment properties for regionalization of three streamflow signatures (runoff coefficient, Base Flow Index, and seasonal elasticity of flow) under land use change and variable data quality in 44 sub-catchments in Thailand. The regression equations derived from regionalization failed to predict the streamflow indices except for runoff coefficient under non-stationary conditions. Halder et al. (2015) and Tian et al. (2014) are the two recent studies on land use change in the Indian context but are deprived of streamflow analysis in ungauged catchments. Li et al. (2009) and Yin et al. (2016) used the Soil and Water Assessment Tool (SWAT) for the Heihe catchment and Jinghe River basin, respectively, to analyze the impacts of land use/land cover (LULC) and climate variability on streamflow, producing satisfactory results. Even though the non-stationary issue was well addressed in the above two cases, the catchments are gauged. Hence, very few documents are available at the moment related to streamflow estimation in ungauged catchments considering the non-stationary environment.

The importance of streamflow estimation in ungauged catchments across the globe is a well-documented subject. However, in a developing country like India, the issue of continuous streamflow regionalization is not prevalent. The additional factors like changing climate and LULC make the regionalization studies more challenging. Considering the above situation, the objective of this study is to enumerate the isolated as well as combined effects of LULCC and climate variability on streamflow and associated uncertainty at an annual scale in ungauged catchments using hydrological model-dependent regionalization techniques. The hydrological model in question is SWAT, while the two regionalization techniques are Inverse Distance Weighted (IDW) and Kriging. Application of SWAT for quantifying the non-stationary effect on prediction is acknowledged in previous studies by researchers. However, application of this procedure on ungauged catchments is a relatively new subject which needs further attention.

METHODOLOGY

Study area and data

India is home to a diversity of climatic regions. Thirty-two catchments (Figure 1) from Eastern and Southern India, spreading over nine river basins, were selected as the study area. All these catchments are headwater regions of nine different river basins, where flow regulation is minimal. The selected 32 catchments are spread over the following states: Jharkhand, West Bengal, Chhattisgarh, Odisha, Telangana, Andhra Pradesh, Tamil Nadu, Karnataka, and Maharashtra. Mean annual rainfall in the first four states is 1,400 mm and for the rest, it is around 1,000–1,100 mm (http://www.imdpune.gov.in/). The key reason for rainfall in India is the southwest monsoon, from July to October. Apart from that, the southern part of India receives a good share of annual rainfall during October and November due to the northeast monsoon. The central and eastern parts of India experience a high temperature of 44 °C during May and June, whereas December and January are the coldest months with temperatures as low as 5 °C. In southern India, summer is hot and humid with a high temperature of 45 °C. However, winter is pleasant with a minimum temperature of 15 °C. The topography and geographical area of the catchments in these regions accommodate large variation.

Figure 1

Location of 32 Indian catchments along with specified representative catchment.

Figure 1

Location of 32 Indian catchments along with specified representative catchment.

Daily rainfall data of the study area were collected from India Meteorological Department (IMD) for the period of the year 1990 to 2011. Only those rain gauge stations where a maximum of 25% data were missing were considered. Similarly, daily streamflow measurement records were collected for the same period from the India-WRIS database (http://www.india-wris.nrsc.gov.in/wris.html) which is jointly operated by Central Water Commission, India (CWC) and Indian space research organization. The land use database was collected from the Oak Ridge National Laboratory Distributed Active Archive Center (ORNL DAAC) (https://daac.ornl.gov) and soil data were collected from the Food and Agriculture Organization of the United Nations (http://www.fao.org/soils-portal/).

SWAT (Soil and Water Assessment Tool)

SWAT is a semi-distributed, continuous hydrological model which functions on a daily time scale to predict the impact of land management on water, sediment and water quality of large basins (Arnold et al. 1998). Weather, hydrology, erosion, plant growth, nutrients, pesticides, land management, and stream routing are some of the key modules of SWAT (Arnold et al. 1998). This model divides a large watershed into many sub-watersheds, which are further sub-divided into a series of hydrological response units (HRUs), a unique combination of land use, soil, and slope. When parameters are simulated, they are first estimated for each HRU and then accumulated for the whole watershed as a weighted average (Abbaspour et al. 2007). For weather data, SWAT uses the data from that gauging station which is nearer to the centroid of each sub-watershed. Variable storage or the Muskingum method is used for routing of calculated runoff and sediment yield through the stream to the outlet from each sub-watershed. Surface runoff is estimated based on the modified Soil Conservation Service (SCS) curve number method, while peak runoff is calculated based on the modified rational formula (Chow et al. 1988; Rostamian et al. 2008). Water is stored in four volumes: snow, soil profile (0–2 m), shallow aquifer (2–25 m) and deep aquifer (Abbaspour et al. 2007). Potential ET is estimated using either of these three methods: Penman–Monteith (Monteith 1965), Priestley–Taylor (Priestley & Taylor 1972), or Hargreaves (Hargreaves & Samani 1985) method, subject to the availability of data. The modified soil loss equation (MUSLE) is used to calculate sediment yield and was developed by Williams & Berndt (1977). Nutrient cycles are found to be similar to those of the EPIC model (Williams et al. 1984). The water quality modelling part of the model is based on QUAL2E (Brown & Barnwell 1987). The water balance equation that governs SWAT (Neitsch et al. 2005) is explained below: 
formula
(1)
where is the final soil water content (mm); is the initial soil water content, is the amount of precipitation, is the amount of surface runoff, is the amount of ET, is the amount of water entering the vadose zone from the soil profile, and is the amount of return flow on day i (mm).

Calibration, validation and uncertainty

The need for calibration of a hydrological model is to determine the optimized parameter set which might have an impact on the catchment response. These parameters may or may not have any physical meaning, but they have a high level of uncertainty associated with them (Heuvelmans et al. 2004). The first period (1990–2000) was considered as the calibration period, while the second period (2001–2011) was treated as the validation period. Before proceeding with the calibration, sensitivity analysis was performed to identify the most sensitive parameters out of a group of 25 model parameters. The parameters used for calibration are presented in Table 1 with the physical significance of each parameter and range of parameters used. SWAT-CUP was used to calibrate the model against observed discharge. In the present study, the Sequential Uncertainty Fitting (SUFI-2) tool was used to assess the uncertainty associated with the prediction. In SUFI-2, uncertainty data associated with input parameters are represented as uniform distributions, however model output uncertainty is quantified by the 95% prediction uncertainty (95PPU) calculated at the 2.5% and 97.5% levels of the cumulative distribution of output variables attained through Latin Hypercube sampling (Abbaspour et al. 2007; Yang et al. 2008). SUFI-2 represents uncertainties of all sources through parameter uncertainty in the hydrological model. It is a bracket of good solutions represented by the 95PPU, generated by certain parameter ranges. P-factor and R-factor are the two statistical terms used for quantification of the uncertainty in terms of a single signal. P-factor is the percentage of observed data bracketed by the 95PPU and R-factor is the average thickness of the 95PPU band. Theoretically, the value for the P-factor ranges between 0 and 100%, whereas R-factor ranges between 0 and infinity. For P-factor, a value of greater than 70% is recommended, whereas having an R-factor of about 1 for discharge calibration and validation is good (Abbaspour 2015). 95PPU should envelope most of the observations, at the same time it is desirable to have a small envelope.

Table 1

SWAT parameters used in the study

Parameter Description Range of parameter
 
Fitted value (average) 
Min Max 
*r__CN2 Initial SCS runoff curve number for moisture condition II −0.2 0.2 −0.033 
*v__GW_DELAY Groundwater delay time (days) 30 450 320.979 
v__ALPHA_BF Baseflow alpha factor (1/days) 0.451 
*a__GWQMN Threshold depth of water in the shallow aquifer required for return flow to occur (mm H2O) 25 12.178 
v__GW_REVAP Groundwater revap coefficient 0.2 0.102 
r__REVAPMN Threshold depth of water in the shallow aquifer for revap to occur (mm H2O) 10 4.267 
r__SOL_AWC Available water capacity of the first soil layer (mm H2O/mm soil) −0.2 0.4 0.142 
v__ESCO Soil evaporation compensation factor 0.8 0.879 
r__SOL_K Saturated hydraulic conductivity (mm/h) −0.8 0.8 0.109 
v__ALPHA_BNK Baseflow alpha factor for bank storage (days) 0.375 
v__CH_K2 Effective hydraulic conductivity in main channel alluvium (mm/h) 130 118.208 
v__EPCO Plant uptake compensation factor −0.8 0.103 
r__HRU_SLP Average slope steepness (m/m) 0.2 0.133 
v__CH_N2 Manning's n value for the main channel 0.3 0.195 
r__OV_N Manning's n value for overland flow −0.2 −0.052 
r__SLSUBBSN Average slope length (m) 0.2 0.188 
r__SOL_BD Moist bulk density (Mg/m3−0.5 0.6 0.147 
V_SURLAG Surface runoff lag coefficient 10 5.430 
Parameter Description Range of parameter
 
Fitted value (average) 
Min Max 
*r__CN2 Initial SCS runoff curve number for moisture condition II −0.2 0.2 −0.033 
*v__GW_DELAY Groundwater delay time (days) 30 450 320.979 
v__ALPHA_BF Baseflow alpha factor (1/days) 0.451 
*a__GWQMN Threshold depth of water in the shallow aquifer required for return flow to occur (mm H2O) 25 12.178 
v__GW_REVAP Groundwater revap coefficient 0.2 0.102 
r__REVAPMN Threshold depth of water in the shallow aquifer for revap to occur (mm H2O) 10 4.267 
r__SOL_AWC Available water capacity of the first soil layer (mm H2O/mm soil) −0.2 0.4 0.142 
v__ESCO Soil evaporation compensation factor 0.8 0.879 
r__SOL_K Saturated hydraulic conductivity (mm/h) −0.8 0.8 0.109 
v__ALPHA_BNK Baseflow alpha factor for bank storage (days) 0.375 
v__CH_K2 Effective hydraulic conductivity in main channel alluvium (mm/h) 130 118.208 
v__EPCO Plant uptake compensation factor −0.8 0.103 
r__HRU_SLP Average slope steepness (m/m) 0.2 0.133 
v__CH_N2 Manning's n value for the main channel 0.3 0.195 
r__OV_N Manning's n value for overland flow −0.2 −0.052 
r__SLSUBBSN Average slope length (m) 0.2 0.188 
r__SOL_BD Moist bulk density (Mg/m3−0.5 0.6 0.147 
V_SURLAG Surface runoff lag coefficient 10 5.430 

*v_ means the existing parameter value is to be replaced by a given value, *a_ stands for a given value that is added to the existing parameter value, *r_ means that an existing parameter value is multiplied by 1+ the given value.

Regionalization methods

Inverse distance weighted

IDW is a multivariate interpolation technique. This method works with a known set of points scattered around the unknown point, the value of which is calculated using a weighted average of the values of the known points. In the case of catchment studies, it is based on the inverse spatial distance of the centroid of the catchment. According to Shepard (1968), the IDW equation used in the present study can be demonstrated as 
formula
(2)
where represents the parameter value of the ungauged catchment, is the weight function, is the parameter value of the gauged catchment and n is the number of catchments. Weight function is calculated as: 
formula
(3)

The distance between the two centroids of the gauged and ungauged catchment is represented by .

For estimation of model parameters in ungauged catchments, all of the 32 catchments were treated as ungauged, one by one.

Kriging

Kriging is an advanced and more sophisticated geostatistical technique than IDW, developed to interpolate spatially autocorrelated variables. It is used to interpolate values from sample data onto a grid of points for contouring (Vandewiele & Elias 1995). Unlike IDW, where the weight is exclusively dependent on distance to the unknown point, Kriging is based on the overall spatial organization of the measured points. To quantify the autocorrelation, the weight depends on a fitted model to the measured points, the distance to the prediction site, and the spatial relationships between the measured values around the prediction site. One of the basic principles of geography quantified in spatial autocorrelation is that things that are closer together are more alike than things that are farther apart. Model fitting starts with an empirical semivariogram for all pairs of locations parted by distance h such as: 
formula
(4)
which estimates the difference squared between the values of the paired locations. A semivariogram describes the spatial autocorrelation of the measured sample points. The subsequent step is to fit a model to the points forming the empirical semivariogram. The empirical semivariogram offers facts only on the spatial autocorrelation of datasets but information for all possible directions and distances. This necessitates the model fitting, which follows. To fit a model to the empirical semivariogram, the exponential function was selected. Afterwards, each pair of points were detected, and the model was fit through them. After the autocorrelation was carried out, prediction was made using the fitted model. Out of the many Kriging methods, ordinary Kriging based on an exponential variogram was used for interpolation with a nugget of 10% of observed variance, and a sill equal to the variance of the model parameter values. One of the reasons for selecting Kriging was the encouraging results obtained by previous researchers by applying it on a large number of basins (Merz & Blöschl 2004; Parajka et al. 2005; Parajka et al. 2007; Samuel et al. 2011). Furthermore, this method does not require additional parameter calibration.

Change in LULC and climate

Agricultural intensification, urbanization and deforestation in addition to other human interventions are some of the primary reasons for land use change. The impact assessment of land cover change on streamflow regionalization was analyzed using land use decadal data of India for the years 1995 and 2005 (Roy et al. 2015). Figure 2 shows the two LULC maps of India for the years 1995 and 2005. The supervised image classification approach with the maximum likelihood classification tool was implemented to classify the LULC image into 19 classes. Forest and agricultural lands are the two major land use classes for most of the catchments under study.

Figure 2

LULC map of India for the year (a) 1995 and (b) 2005.

Figure 2

LULC map of India for the year (a) 1995 and (b) 2005.

The trend analysis was carried out for annual rainfall and mean annual temperature to analyze the change in pattern using the Mann–Kendall test, a non-parametric method. The indicators of the Mann–Kendall test are explained as follows (Xu et al. 2003; Partal & Kahya 2006; Li et al. 2009): 
formula
(5)
where 
formula
(6)
 
formula
(7)
 
formula
(8)
are the sequential data values is the sample size, t stands for the extent of any given tie while, represents the summation over all ties. The null hypothesis H0, i.e., there is no trend in the dataset, is accepted if Z. The Kendall slope which quantifies the magnitude of the monotonic trend is expressed as: 
formula
(9)
where 1 < j<i<n. represents the median over all combinations of record pairs for the whole dataset. A positive value of indicates an ‘increasing trend’ while a negative value suggests a ‘decreasing trend’.

Impact assessment of LULCC and climate variability

The impact of LULCC on streamflow was analyzed by generating different cases mentioned in Table 2 for calibration and validation periods. SWAT simulations were carried out using land use data of the year 1995 for period 1 and land use data of the year 2005 for second period. The two selected land use maps were able to represent the scenario of the respective periods quite effectively. The cases mentioned above may be explained as follows: case 1 represents a scenario where simulation was carried out using LULC data of the period 1 along with the climatological data of the same period. However, in the second case, simulation was performed using climate data of period 1, and land use data of period 2. According to Li et al. (2009), during the simulation, one factor is variable while others are kept constant, which explains the effect of the variable factor on runoff simulation. Case 1 represents the ‘base’ scenario where the impact of LULC change along with climatological factors upon runoff can be quantified. The difference between case 1 and 2 will quantify the isolated impact of LULC change on streamflow. Case 3 uses land use data of the first period while the climatological data are from the second period. Nonetheless, case 4 was performed using both climatological and land use data of the second period. As explained above, the contrast between different cases quantifies the effects of the above-mentioned two factors on streamflow. However, the four cases mentioned were repeated for the catchment when it was treated as ungauged. The percentage change in streamflow due to LULC and climate change was estimated on a yearly scale.

Table 2

Streamflow simulation using four generated cases

Events Description 
Case 1 LULC data of year 1995 with climatological data from 1990–2000 
Case 2 LULC data of year 2005 with climatological data from 1990–2000 
Case 3 LULC data of year 1995 with climatological data from 2001–2011 
Case 4 LULC data of year 2005 with climatological data from 2001–2011 
Events Description 
Case 1 LULC data of year 1995 with climatological data from 1990–2000 
Case 2 LULC data of year 2005 with climatological data from 1990–2000 
Case 3 LULC data of year 1995 with climatological data from 2001–2011 
Case 4 LULC data of year 2005 with climatological data from 2001–2011 

Model performance evaluation

Based on a previously published recommendation (Moriasi et al. 2007), NSE, RMSE-observations standard deviation ratio (RSR) and Percent Bias (PBIAS) were chosen as indicators to judge the model performance.

NSE (Nash & Sutcliffe 1970) is a measure of correlation, bias and variability (Gupta et al. 2009), which is one of the most widely used criteria for performance evaluation all over the world. NSE varies from −∞ to 1. A value of 1, points towards a perfect match between observed and simulated value. A value of NSE ≤ 0 suggests that the model is no better than the observed mean values (Gupta et al. 2009). Usually, RSR varies from zero to a large positive value. The optimal value ‘0’ represents null RMSE or residual variation, i.e., a perfect model simulation (Moriasi et al. 2007). RSR combines the benefits of error index statistics which includes a scaling or normalization factor (Moriasi et al. 2007). PBIAS measures the average tendency of the simulated data to be larger or smaller than their observed counterparts, it detects poor model performances (Gupta et al. 2009). A positive PBIAS value represents the bias towards underestimation while a negative value suggests a model bias towards overestimation and zero represents the optimal value (Gupta et al. 1999). A value of ±25% can be considered as satisfactory result in the case of streamflow estimation (Van Liew et al. 2007). The expressions for NSE, RSR and PBIAS are: 
formula
(10)
 
formula
(11)
 
formula
(12)
where and represent the observed and simulated streamflow respectively, is the mean of the observed streamflow and t is the number of observations.

RESULTS

The change in LULC may be attributed to the alteration which occurred in forest cover and cropland as most of the catchments out of the 32, experienced major changes in the above-mentioned two land use classes (Figure 2). Only two catchments located in the southern part of the study area experienced a major change in built-up area. The extent of change in forest cover for some catchments is significant while it is marginal for others. However, there are few catchments where the change is minimal which can be neglected. One catchment that is highlighted in Figure 1, where the land use change is significant was selected for the presentation of the non-stationarity analysis for simplified depiction of results. Forest cover (41.53%) and cropland (40.58%) are the two dominant components of the land cover for the selected catchment. Forest cover and cropland area decreased (2.26% and 12.31%, respectively) from the first period to the second period. In addition, there was an increase in wastelands by 54.02% during the period of study. Social development and population growth might be the factors contributing largely to the change in land use pattern (Yin et al. 2016).

The Mann–Kendall test explained the climatic variability of the catchment for the selected period (1990–2011). The trend analysis indicates a decrease in annual rainfall volume for two of the rain gauges (station 1 and 4), while no significant trends are observed for the other two (Figure 3). However, according to the Kendall slope, mean annual temperature tended to increase for three weather stations during that period (Figure 4). The overall results suggest that the climate in and around the catchment is getting warmer and drier. The specific results of the Mann–Kendall test are presented in Table 3 which suggests that the decreasing trend in annual rainfall is significant for two stations (station 1 and 4).

Table 3

Statistics of the Mann–Kendall test for annual rainfall and temperature

Gauging Station Rainfall
 
Temperature
 
β Ζ p β Ζ p 
Station 1 −1.2 −1.75 0.06 3.89 ** 
Station 2 −0.06 −0.17  0.2 2.93 ** 
Station 3 −0.54 −0.85  0.01 0.17  
Station 4 −1.17 −1.58    
Gauging Station Rainfall
 
Temperature
 
β Ζ p β Ζ p 
Station 1 −1.2 −1.75 0.06 3.89 ** 
Station 2 −0.06 −0.17  0.2 2.93 ** 
Station 3 −0.54 −0.85  0.01 0.17  
Station 4 −1.17 −1.58    

*means significant at P < 0.1.

**means significant at P < 0.05.

Figure 3

Annual rainfall series from the Mann–Kendall test.

Figure 3

Annual rainfall series from the Mann–Kendall test.

Figure 4

Annual temperature series from the Mann–Kendall test.

Figure 4

Annual temperature series from the Mann–Kendall test.

The performance of SWAT for the 32 catchments can be termed as satisfactory during calibration and validation periods as per the model performance criteria set by Moriasi et al. (2007). The minimum, maximum, mean and median values of three model performance judging criteria are listed in Table 4. SWAT largely performed better for catchments present on the northern side of the study area than those catchments located on the southern side. The missing records of hydro-climatological data for a few southern catchments might have resulted in such output. Although SWAT performance deteriorated marginally when catchments were treated as ungauged in turn with the application of IDW and Kriging, it is still within the acceptable limit. As mentioned earlier, the effect of non-stationarity due to changing climate and LULC was demonstrated for a single catchment as the representation of the whole analysis. The results of the modeling for this catchment can be put within the bracket of ‘good’ results. As demonstrated below in Figure 5, the simulated monthly streamflow matched with the observed streamflow except for a few peak flow sections. Overestimation of peak flow is observed. The values of NSE, RSR and PBIAS were 0.75, 0.29 and −15.5% respectively during calibration. While during the validation period these were 0.7, 0.33 and −15% respectively. A P-factor of 0.83 and R-factor of 1.06 indicate the output uncertainty is within the acceptable limits. The values of NSE, RSR and PBIAS were 0.71, 0.33 and −17.2% (calibration), and 0.69, 0.33 and −17.55% (validation), when IDW was applied. The Kriging technique produced similar results as the IDW approach with values of NSE, RSR and PBIAS at 0.72, 0.29, and −16.41% (calibration), and 0.71, 0.3 and −16.13% (validation). P-factor and R-factor for IDW and Kriging techniques are found to be above the threshold value that can be termed as ‘acceptable’.

Table 4

SWAT performance in terms of NSE, RSR and absolute PBIAS for the 32 catchments

Model performance criteria Statistical measures Calibration
 
Validation
 
Gauged IDW Kriging Gauged IDW Kriging 
NSE Minimum 0.59 0.51 0.53 0.5 0.45 0.45 
Mean 0.7 0.61 0.62 0.65 0.55 0.55 
Median 0.7 0.62 0.62 0.65 0.59 0.59 
Maximum 0.82 0.7 0.71 0.77 0.66 0.68 
RSR Minimum 0.07 0.08 0.08 0.07 0.09 0.08 
Mean 0.31 0.36 0.35 0.3 0.32 0.32 
Median 0.3 0.35 0.33 0.3 0.35 0.32 
Maximum 0.85 0.86 0.86 0.89 0.9 0.9 
PBIAS Minimum 3.22 4.55 4.29 4.21 4.69 4.8 
Mean 15.02 21.6 19.5 15.66 22.1 20 
Median 14.55 15.23 15.7 16.28 18.4 17.84 
Maximum 41.22 42.25 43.21 42.26 45 46.2 
Model performance criteria Statistical measures Calibration
 
Validation
 
Gauged IDW Kriging Gauged IDW Kriging 
NSE Minimum 0.59 0.51 0.53 0.5 0.45 0.45 
Mean 0.7 0.61 0.62 0.65 0.55 0.55 
Median 0.7 0.62 0.62 0.65 0.59 0.59 
Maximum 0.82 0.7 0.71 0.77 0.66 0.68 
RSR Minimum 0.07 0.08 0.08 0.07 0.09 0.08 
Mean 0.31 0.36 0.35 0.3 0.32 0.32 
Median 0.3 0.35 0.33 0.3 0.35 0.32 
Maximum 0.85 0.86 0.86 0.89 0.9 0.9 
PBIAS Minimum 3.22 4.55 4.29 4.21 4.69 4.8 
Mean 15.02 21.6 19.5 15.66 22.1 20 
Median 14.55 15.23 15.7 16.28 18.4 17.84 
Maximum 41.22 42.25 43.21 42.26 45 46.2 
Figure 5

95PPU plot of observed and simulated streamflow at monthly time scale for the period 1990–2000 (calibration) to 2001–2011 (validation).

Figure 5

95PPU plot of observed and simulated streamflow at monthly time scale for the period 1990–2000 (calibration) to 2001–2011 (validation).

Impact assessment of LULCC and climate variability

The impact of LULCC and climate variability was analyzed through the four generated cases mentioned earlier. The isolated impact of land use change was estimated by comparing case 1 with case 2. The difference between these two cases indicates an increase in annual runoff of 3.52%. The decrease in forest cover might have played a key role in the increased runoff. The isolated impact of the variable climate was analyzed with comparison between case 1 and 3. A decrease in annual runoff volume of 21.7% is detected. Overall results show that the combined effect of LULCC and climate variability results in a decrease in annual runoff volume of 22.53%. The contribution of climate variability to the change in streamflow overshadows the influence of LULC change. The percentage change in annual streamflow volume due to the combined effect of climate variability and LULC change is greater than that of the summation of the isolated impacts of these two factors. The results indicate that apart from these two factors mentioned above few other hydro-climatological or geomorphological factors might have played a role in the change in streamflow quantity.

The change in simulated streamflow between the two periods was analyzed at an annual scale. IDW and Kriging were applied to transfer the optimized model parameters from nearby gauged catchments to the selected catchment treated as ungauged. An increase in annual runoff volume was observed during IDW and Kriging application due to the isolated effect of LULC change (case 1–case 2). For IDW and Kriging, streamflow increased by 2.56% and 2.6%, respectively, due to the isolated impact of LULCC. However, a decrease in annual runoff volume was observed during application of both regionalization techniques due to changing climate (case 1–case 3). The percentage decreases in annual runoff volume are 21.09% and 21.6% for IDW and Kriging, respectively. The combined effect of these two factors resulted in a decrease in annual runoff volume of 21.8% and 22% for IDW and Kriging, respectively. The results obtained for the ungauged catchment scenario matched closely with that of the gauged scenario.

DISCUSSION

This study is an attempt to address the non-stationarity issue through change in LULC and climate during spatial regionalization of continuous streamflow in an ungauged catchment. The small change in annual streamflow due to LULC change indicates that human intervention in the catchment area is minimal. The conversion of cropland to fallow land in addition to decreased forest cover might have played an important role in the increase of annual runoff. However, the variation in climate had a prominent contribution towards change in streamflow quantity. Although there was an increase in annual rainfall during the period (1993–1995), a decreasing trend was observed afterward. The variability in streamflow during the study period is more related to change in rainfall during the monsoon period: one reason for this might be the variable rainfall in monsoon-dominated areas is more intense and distinct than LULCC. Even though climate variability has greater impact on streamflow in the catchment, the degree of correlation between rainfall and streamflow was non-linear due to the combined effects. When the catchment was treated as ungauged, IDW and Kriging produced more or less similar results to that of the gauged scenario. Moving from the gauged to ungauged case, a decreased model performance is observed which is known as a spatial loss in model accuracy. Likewise, moving from calibration to validation stage, the decrease in model performance is known as a temporal loss in model accuracy (Parajka et al. 2005). The satisfactory performance of IDW and Kriging indicates that the presence of well-gauged catchments in the proximity of an ungauged catchment is useful.

Regional climate, soil types and vegetation are certain factors that make the regionalization case specific. Most of the SWAT parameters are related to the factors mentioned above. The model concerned and the objective of the study decide the suitability of the regionalization technique. A change in model or parameters mentioned above might have produced different results. It is essential to recognize that there is no single signal representing model output, but rather an envelope of good solutions expressed by the 95PPU, produced by certain parameter ranges. Quantifying the effect of different uncertainties on model outputs separately is highly desired. However, it is very hard to achieve. Hence, the combined effect is enumerated on model outputs. According to the overall degree of prediction uncertainty, the Kriging approach produced marginally better results than the IDW method for this study.

As mentioned earlier, most of the regionalization studies were carried out considering the peripheral atmosphere to be stationary. The result of this study suggests that the ever-changing climate and LULC has an impact on streamflow. Hence, it is highly recommended for hydrologists, engineers and policymakers to consider these changing parameters during planning and management of water resources and ecology, with the associated socioeconomic impacts.

CONCLUSIONS

The application of SWAT for streamflow simulation in non-stationary environments is justified by the satisfactory results obtained in terms of NSE, RSR and PBIAS. The performance of the above-mentioned model is equally good when the catchments are treated as ungauged. IDW and Kriging are found to be competitive with each other, as both approaches produced more or less similar results. The presence of well-gauged catchments around the pseudo ungauged catchment is crucial for better performance of these two geographical distance based regionalization techniques. Both LULC and climate varied during the study period constantly. The combined effect of these two factors resulted in a decrease in annual runoff of 22.53%. However, an isolated impact of LULC brought about a moderate increase in annual runoff volume of 3.52%. On the other hand, a decrease of 21.7% in annual runoff volume was observed due to the isolated effect of climate change as the climate is getting warmer and drier. Out of these two factors, variable climate played a predominant role in surface water hydrology. Hence regionalization studies must consider the non-stationary effects during any hydrological predictions for better water resources planning. A comparison between a stationary and non-stationary condition may further highlight the importance of LULCC and climate change studies. Consideration of a few more geomorphological factors in the analysis should be a future objective to identify the impact of other non-stationary parameters on streamflow.

ACKNOWLEDGEMENTS

The authors would like to acknowledge the IMD and CWC, India for providing the hydro-climatological data required for this study.

REFERENCES

REFERENCES
Abbaspour
K. C.
2015
SWAT-CUP: SWAT Calibration and Uncertainty Programs – A User Manual
.
Swiss Federal Institute of Aquatic Science and Technology
,
Duebendorf
,
Switzerland
.
Abbaspour
K. C.
,
Yang
J.
,
Maximov
I.
,
Siber
R.
,
Bogner
K.
,
Mieleitner
J.
,
Zobrist
J.
&
Srinivasan
R.
2007
Modelling hydrology and water quality in the pre-alpine/alpine Thur watershed using SWAT
.
Journal of Hydrology
333
(
2–4
),
413
430
.
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
.
Blöschl
G.
2016
Predictions in ungauged basins – where do we stand?
In:
7th International Water Resources Management Conference of ICWRS
,
7
.
Blöschl
G.
&
Sivapalan
M.
1995
Scale issues in hydrological modelling: a review
.
Hydrological Processes
9
,
251
290
.
Blöschl
G.
,
Sivapalan
M.
,
Wagener
T.
,
Viglione
A.
&
Savenije
H.
2013
Runoff prediction in ungauged basins: synthesis across processes, places and scales
.
Eos, Transactions American Geophysical Union
95
(
2
),
465
.
Brown
L.
&
Barnwell
T.
1987
The Enhanced Stream Water Quality Models QUAL2E and QUAL2E-UNCAS
.
Chow
V. T.
,
Maidment
D. R.
&
Mays
L. W.
, (eds)
1988
Applied Hydrology
.
McGraw-Hill Inc.
,
New York
,
USA
.
Croke
B. F. W.
,
Merritt
W. S.
&
Jakeman
A. J.
2004
A dynamic model for predicting hydrologic response to land cover changes in gauged and ungauged catchments
.
Journal of Hydrology
291
(
1–2
),
115
131
.
Götzinger
J.
&
Bárdossy
A.
2007
Comparison of four regionalisation methods for a distributed hydrological model
.
Journal of Hydrology
333
(
2–4
),
374
384
.
Gupta
H. V.
,
Sorooshian
S.
&
Yapo
P. O.
1999
Status of automatic calibration for hydrologic models: comparison with multilevel expert calibration
.
Journal of Hydrologic Engineering
4
(
2
),
135
143
.
Gupta
H. V.
,
Kling
H.
,
Yilmaz
K. K.
&
Martinez
G. F.
2009
Decomposition of the mean squared error and NSE performance criteria: implications for improving hydrological modelling
.
Journal of Hydrology
377
(
1–2
),
80
91
.
Halder
S.
,
Saha
S. K.
,
Dirmeyer
P. A.
,
Chase
T. N.
&
Goswami
B. N.
2015
Investigating the impact of land-use land-cover change on Indian summer monsoon daily rainfall and temperature during 1951–2005 using a regional climate model
.
Hydrology and Earth System Sciences Discussions
12
(
7
),
6575
6633
.
Hargreaves
G.
&
Samani
Z. A.
1985
Reference crop evapotranspiration from temperature
.
Applied Engineering in Agriculture
1
,
96
99
.
Heuvelmans
G.
,
Muys
B.
&
Feyen
J.
2004
Evaluation of hydrological model parameter transferability for simulating the impact of land use on catchment hydrology
.
Physics and Chemistry of the Earth, Parts A/B/C
29
(
11–12
),
739
747
.
Ibrahim
B.
,
Wisser
D.
,
Barry
B.
,
Fowe
T.
&
Aduna
A.
2015
Hydrological predictions for small ungauged watersheds in the Sudanian zone of the Volta basin in West Africa
.
Journal of Hydrology: Regional Studies
4
,
386
397
.
Leclerc
M.
&
Ouarda
T. B. M. J.
2007
Non-stationary regional flood frequency analysis at ungauged sites
.
Journal of Hydrology
343
(
3–4
),
254
265
.
Li
C. Z.
,
Zhang
L.
,
Wang
H.
,
Zhang
Y. Q.
,
Yu
F. L.
&
Yan
D. H.
2012
The transferability of hydrological models under nonstationary climatic conditions
.
Hydrology and Earth System Sciences
16
(
4
),
1239
1254
.
Merz
R.
&
Blöschl
G.
2004
Regionalisation of catchment model parameters
.
Journal of Hydrology
287
(
1–4
),
95
123
.
Mishra
A. K.
&
Coulibaly
P.
2009
Developments in hydrometric network design: a review
.
Reviews of Geophysics
47
(
2
),
1
24
.
Monteith
J. L.
1965
Evaporation and Environment
.
Symposia of the Society for Experimental Biology
19, 205–234. Cambridge University Press, Cambridge, UK
.
Moriasi
D. N.
,
Arnold
J. G.
,
Van Liew
M. W.
,
Binger
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
.
Nathan
R. J.
&
McMahon
T. A.
1990
Identification of homogeneous regions for the purposes of regionalisation
.
Journal of Hydrology
121
(
1–4
),
217
238
.
Neitsch
S. L.
,
Arnold
J. G.
,
Kiniry
J. R.
&
Williams
J. R.
2005
Soil and Water Assessment Tool User's Manual Version 2005. 494
.
Oudin
L.
,
Kay
A.
,
Andréassian
V.
&
Perrin
C.
2010
Are seemingly physically similar catchments truly hydrologically similar?
Water Resources Research
46
(
11
),
1
15
.
Parajka
J.
,
Merz
R.
&
Blöschl
G.
2005
A comparison of regionalisation methods for catchment model parameters
.
Hydrology and Earth System Sciences Discussions
2
(
2
),
509
542
.
Parajka
J.
,
Blöschl
G.
&
Merz
R.
2007
Regional calibration of catchment models: potential for ungauged catchments
.
Water Resources Research
43
(
6
),
1
16
.
Partal
T.
&
Kahya
E.
2006
Trend analysis in Turkish precipitation data
.
Hydrological Processes
20
(
9
),
2011
2026
.
Priestley
C. H. B.
&
Taylor
R. J.
1972
On the assessment of surface heat flux and evaporation using large-scale parameters
.
Monthly Weather Review
100
(
2
),
81
92
.
Razavi
T.
&
Coulibaly
P.
2012
Streamflow prediction in ungauged basins: review of regionalization methods
.
Journal of Hydrologic Engineering
18
(
8
),
958
975
.
Rostamian
R.
,
Jaleh
A.
,
Afyuni
M.
,
Mousavi
S. F.
,
Heidarpour
M.
,
Jalalian
A.
&
Abbaspour
K. C.
2008
Application of a SWAT model for estimating runoff and sediment in two mountainous basins in central Iran
.
Hydrological Sciences Journal
53
(
5
),
977
988
.
Roy
P. S.
,
Roy
A.
,
Joshi
P. K.
,
Kale
M. P.
,
Srivastava
V. K.
,
Srivastava
S. K.
,
Dwevidi
R. S.
,
Joshi
C.
,
Behera
M. D.
,
Meiyappan
P.
,
Sharma
Y.
,
Jain
A. K.
,
Singh
J. S.
,
Palchowdhuri
Y.
,
Ramachandran
R. M.
,
Pinjarla
B.
,
Chakravarthi
V.
,
Babu
N.
,
Gowsalya
M. S.
,
Thiruvengadam
P.
,
Kotteeswaran
M.
,
Priya
V.
,
Yelishetty
K. M. V. N.
,
Maithani
S.
,
Talukdar
G.
,
Mondal
I.
,
Rajan
K. S.
,
Narendra
P. S.
,
Biswal
S.
,
Chakraborty
A.
,
Padalia
H.
,
Chavan
M.
,
Pardeshi
S. N.
,
Chaudhari
S. A.
,
Anand
A.
,
Vyas
A.
,
Reddy
M. K.
,
Ramalingam
M.
,
Manonmani
R.
,
Behera
P.
,
Das
P.
,
Tripathi
P.
,
Matin
S.
,
Khan
M. L.
,
Tripathi
O. P.
,
Deka
J.
,
Kumar
P.
&
Kushwaha
D.
2015
Development of decadal (1985-1995-2005) land use and land cover database for India
.
Remote Sensing
7
(
3
),
2401
2430
.
Samuel
J.
,
Coulibaly
P.
&
Metcalfe
R. A.
2011
Estimation of continuous streamflow in Ontario ungauged basins: comparison of regionalization methods
.
Journal of Hydrologic Engineering
16
(
5
),
447
459
.
Shepard
D.
1968
A two-dimensional interpolation function for irregularly-spaced data
. In:
23rd ACM National Conference
, pp.
517
524
.
Singh
R.
,
Wagener
T.
,
Van Werkhoven
K.
,
Mann
M. E.
&
Crane
R.
2011
A trading-space-for-time approach to probabilistic continuous streamflow predictions in a changing climate-accounting for changing watershed behavior
.
Hydrology and Earth System Sciences
15
(
11
),
3591
3603
.
Sivapalan
M.
,
Takeuchi
K.
,
Franks
S. W.
,
Gupta
V. K.
,
Karambiri
H.
,
Lakshmi
V.
,
Liang
X.
,
McDonnell
J. J.
,
Mendiondo
E. M.
,
O'Connell
P. E.
,
Oki
T.
,
Pomeroy
J. W.
,
Schertzer
D.
,
Uhlenbrook
S.
&
Zehe
E.
2003
IAHS decade on Predictions in Ungauged Basins (PUB), 2003–2012: shaping an exciting future for the hydrological sciences
.
Hydrological Sciences Journal
48
(
6
),
857
880
.
Swain
J. B.
&
Patra
K. C.
2017a
Streamflow estimation in ungauged catchments using regional flow duration curve: comparative study
.
Journal of Hydrologic Engineering (ASCE)
22
(
7
),
04017010
.
Swain
J. B.
&
Patra
K. C.
2017b
Streamflow estimation in ungauged catchments using regionalization techniques
.
Journal of Hydrology
554
,
420
433
.
Vandewiele
G. L.
&
Elias
A.
1995
Monthly water balance of ungauged catchments obtained by geographical regionalization
.
Journal of Hydrology
170
(
1–4
),
277
291
.
Van Liew
M. W.
,
Veith
T. L.
,
Bosch
D. D.
&
Arnold
J. G.
2007
Suitability of SWAT for the conservation effects assessment project: comparison on USDA agricultural research service watersheds
.
Journal of Hydrologic Engineering
12
(
2
),
173
189
.
Visessri
S.
&
McIntyre
N.
2015
Regionalisation of hydrological responses under land use change and variable data quality
.
Hydrological Sciences Journal
61
(
2
),
150113092349004
.
Williams
J. R.
&
Berndt
H. D.
1977
Sediment yield prediction based on watershed hydrology
.
Transactions of the ASAE
20
(
6
),
1100
1104
.
Williams
J. R.
,
Jones
C. A.
&
Dyke
P. T.
1984
A modeling approach to determining the relationship between erosion and soil productivity
.
Transactions of the ASAE
27
,
129
144
.
Xu
Z. X.
,
Takeuchi
K.
&
Ishidaira
H.
2003
Monotonic trend and step changes in Japanese precipitation
.
Journal of Hydrology
279
(
1–4
),
144
150
.
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–2
),
1
23
.
Yu
P. S.
,
Yang
T. C.
&
Wang
Y. C.
2002
Uncertainty analysis of regional flow duration curves
.
Journal of Water Resources Planning and Management
128
(
6
),
424
430
.