Evaluating two GIS-based semi-distributed hydrological models in the Bhagirathi-Alkhnanda River catchment in India

Uttarakhand, a Himalayan state of India, has been suffering from several natural calamities, such as flash floods, landslides, drought, etc. The increased frequency and magnitude of flash floods, as well as landslides in different parts of the state, are creating havoc for the people of the state. It has destroyed the infrastructure along with the lives and property of the people living in the potentially hazardous areas of the state. The availability of water in the state is primarily due to the monsoon rainfall, as well as the snowmelt discharge from the Himalayan glacier. The Alkhnanda and Bhagirathi river basin of the state receives around 80% of the water from the snowmelt of the Himalayan glacier (G. Government of Uttarakhand, 2014; Ministry & Government, 2014). These rivers also receive ample water in the rainy season. As a result, frequent flash flood events occur in different parts of the state (Rana et al., 2013). As reported, the flash floods have been causing severe destruction in the upper Himalayan region of the state. The state of Uttarakhand is also known for small hydropower plants, which have altered the natural flow of the rivers (Pranuthi et al., 2014). Further, climate change will also alter the spatial and temporal distribution of water in the basin. The impacts of climate change on the water potential of different basins of India have been studied and explained by different researchers (Salvi et al., 2011; Parth Sarthi et al., 2015; Mahmood & Jia, 2016; Pichuka et al., 2017).

The hydrological model is generally used to simulate the hydrologic processes of a basin. The model can be used to understand the present spatio-temporal distribution water in a river basin as well as its future distribution. The hydrological models are categorized into three different classes, i.e., the lumped, distributed, and semi-distributed model (Gebre, 2015). These classifications are based on the spatial variability considered in modeling the rainfall-runoff process of a catchment (Sok & Oeurng, 2016). In the present study, we employed two semi-distributed hydrologic models, i.e., HEC-HMS and SWAT, to simulate the rainfall-runoff process of the three catchments of the Bhagirathi and Alaknanda river basin. The three catchments are Tehri, Jyosimath, and Uttarkashi. The SWAT (Soil and Water Assessment Tool) model was developed by the USDA-Agricultural Research Service. It is a physically based hydrological model that needs basic knowledge about the climate of the catchment, soil properties of the area, topography, cropping pattern, and the information-related land-use practices going on within the specified watershed (Gitau & Chaubey, 2010). A review of the literature has shown that HEC-HMS and SWAT models have been used for streamflow modeling in different river basins in India as well as in the world (Tan et al., 2017). Both the models, i.e., HEC-HMS and SWAT, are semi-distributed, but the operating functions of both models are different. In the SWAT model, five different hydrological inputs such as rainfall, temperature, relative density, wind speed, and sunshine hour are used for simulation of runoff. In the case of the HMS model, the runoff is generated based on the corresponding rainfall data. The model is based on the curve number method and land-use property of the catchment. The SWAT model is generally used to find the impact of different land management practices on water, sediment, etc. In contrast, the HEC-HMS is designed to simulate the complete hydrologic processes of dendritic watershed systems (Yener et al., 2007; Malagò et al., 2015; Sok & Oeurng, 2016).

In the present paper, we have carried out a parallel study of the results obtained from these two different models and their applicability to the hilly catchment of Uttarakhand. For developing the hydrologic model using HEC-HMS, Watershed Management System (WMS) is used to prepare the input files required for HEC-HMS. The catchment was split into three sub-catchments to obtain the different geographical data sets at the sub-watershed level. The hydrological analysis was then performed in HEC-HMS for these sub-watersheds applying all the parameters acquired from the WMS modeling. Generally, any rainfall-runoff model considers the impacts of physical properties of the watershed, such as climatic, topological, land-use data, and soil data across the catchment to predict the corresponding discharge (Halwatura & Najim, 2013; Mahmood et al., 2016). HEC-HMS uses the curve number (CN) method for generating the runoff of a watershed (Song et al., 2011). As such, we have prepared the CN map of the Uttarakhand basin, which is derived from the soil map and land use data of Uttarakhand.

HEC-HMS is a physically based semi-distributed hydrologic model developed by the US Army Corps of Engineers to generate the flow of a basin based on various inputs like soil type, land use property, and rainfall data (Halwatura & Najim, 2013; Gebre, 2015; Sok & Oeurng, 2016). HEC-HMS model structure consists of a basin model, meteorological model, control specifications, and input data (time series data) (Martin et al., 2012; Roy et al., 2013; Sampath et al., 2015). For simulating the outflow, generally, the deficit and constant loss model is used to compute the losses from the catchment, SCS unit hydrograph model is used to transform the flows, and the usual monthly base flow method is used to account for the base flow (Song et al., 2011; Ashish et al., 2012; Pichuka et al., 2017). In contrast, the SWAT model is a physically based semi-distributed model intended to predict the influence of land management practices on water, sediment, and agricultural chemical yields in large heterogeneous watersheds with varying soil, land-use, and management conditions over long periods. In this model, a catchment is split into several sub-watersheds or sub-basins (Gitau & Chaubey, 2010; Arnold et al., 2012; Malagò et al., 2015). Sub-basins are further partitioned into hydrological response units (HRUs) based on soil types, land-use types, and slope classes that provide a tremendous level of spatial simulation (Birhanu et al., 2016; Tan et al., 2017). The model predicts the outflow at various HRU using the water balance equation. As input, the SWAT model requires a digital elevation map (DEM), soil map, land-use map, and climate data for modeling a watershed. After data preparation, the model is set up by performing the following four main steps: (i) watershed delineation, (ii) hydrologic response unit (HRU) definition, (iii) model run, and (iv) calibration and validation of the model (Fleming & Doan, 2009; Dhami & Pandey, 2013; Fiseha, 2013). The model calibration and validation were done using the SUFI-2 algorithm, which was executed in the SWAT-CUP 2012.

The present study is conducted on Uttarakhand, an Indian Himalayan state. We have considered two river basins, namely, the Bhagirathi river basin covering the Uttarkashi and Tehri Garhwal catchment, and Alkhnanda river basin in Chamoli and Josh math of Uttarakhand. The river basins are shown in Figure 1. As per the historical flood data, the majority of flash flood events have happened in this area. Figure 2 shows the major cloudburst events that have occurred in the state.

The WMS is used to generate all the hydrological inputs needed for rainfall-runoff modeling. The digital elevation model (DEM) is used for the extraction of physical features of the catchment. The physical characteristics are the total watershed area, the perimeter of the catchment, length of the river, slope, lag time, time of concentration, etc. Figure 3 shows the different processes performed by WMS. Based on the requirement, the 30 m resolution DEM data for the study area has been downloaded from the Shuttle radar topographic mission (SRTM).

Elevation plays an essential role in the hydrological analysis of any catchment. DEM is used to obtain the geographical as well as topographical features of an area (Tan et al., 2017). This DEM downloaded for the study area is shown in Figure 4 (Karim et al., 2016). The DEM was analyzed in WMS to extract basin parameters like slopes, areas of sub-basins, longest flow paths, data related to the elevation of the catchments, water flow direction, streamlines of the river, etc.

After processing of the DEM, we have to define the outlet position. Once we identify the location of outlets in the DEM, the watershed delineation is performed which creates the sub-watersheds that belong to the particular stream portion. The delineation process is done for an outlet, and the WMS forms the watershed as per the contributing area of an outlet. Since three outlets were defined, as per the availability of discharge data, we get three sub-watersheds. Figure 5 shows the location of different outlets, and Figure 6 shows the stream network of Alkhnanda and Bhagirathi river basins, respectively.

The daily precipitation data of different stations of the state were collected from the Indian Meteorological Department (IMD). Also, the daily discharge data at three gauging locations, namely, Tehri, Uttarkashi, and Jyosimath, were obtained from 1985 to 2010, from the Central Water Commission, Lucknow. The primary data of the meteorological stations selected in the present study are presented in Table 1.

Uttarakhand state encounters massive rain in the form of cloudburst, particularly during the monsoon season, which is from June to August, due to its locality at high altitude. Generally, the climate of Uttarakhand is cold, with a high wind velocity throughout the year. The heavy rainfall has usually been seen during July and August. Flash flood generally occurs in the Upper Himalayan region during this period. The average annual precipitation is about 1,900 mm. The rain gauge stations available at Uttarkashi, Jyosimath, and Tehri are examined to study the fluctuations of mean monthly rainfall patterns. In all the selected stations, the maximum rainfall is observed in July and August (Figure 7). The monthly average rainfall during the period from July to September for the watershed of Uttarkashi, Jyosimath, and Tehri is 560, 381, and 580 mm, respectively.

The Thiessen polygon method is applied to calculate average precipitation. The distribution of the rain gauge network in Uttarakhand for rainfall records is shown in Figure 8.

The discharge data are collected from CWC, Lucknow, for different periods at different locations of Uttarakhand. These locations are marked as an outlet for runoff modeling, as shown in Figure 5.

After getting all the primary data for simulation, we need to find out the hydrological nature of the catchment surface. For the estimation of infiltration using the CN method, the curve number of each of the sub-catchments is required. The SCS-CN loss method was used to determine rainfall losses. The selection of this method is based on the availability of data for the region. In general, for the estimation of the CN of any region, the soil map and the land-use/land cover property of the area is needed. For the Uttarakhand catchments, we generated the soil map and land-use map using ArcGIS 10.4 and then combined them to generate the CN map of Uttarakhand (Merwade, n. d.; Mahmood & Jia, 2016). The detailed methodology of CN generation and working of the HMS model is shown in Figure 9 and discussed in detail below and shown in Table 2.

Categorizations/Preparation of soil map of Uttarakhand: The hydrologic soil groups (HSG) are generally divided into four groups, A, B, C, and D. Group A indicates low runoff capacity, and the area has a high infiltration rate. Group B indicates the soil having an average infiltration rate. The soils of group C show moderately fine to rough textures and have a moderate rate of water transmission, and the soils of group D demonstrate slow infiltration and have a very high runoff. The soil map for the Uttarakhand watershed was acquired from the World Soil Database, with a spatial resolution of 1 km (Figure 10(a)) and is also presented in Table 3. Dominant soil type (DOMSOIL) with FAO soil unit codes, name, and determined hydrologic soil group (HSG) located in Uttarakhand are also presented in Figure 10(a). This soil database was generated by FAO in association with the International Institute of Applied Systems Analysis (IIASA) (Halwatura & Najim, 2013; Gautam et al., 2013).

Land-use/Land class categorizations/preparation of LULC map of Uttarakhand: The land-use pattern/land cover significantly affects the hydrological processes in a watershed. For example, forest density considerably affects runoff generation from a watershed. A report released by the Ministry of Agricultural, Government of India shows that there is not much change in the LULC (land-use/land cover) scenario of Uttarakhand in the last decade (Figure 11). As such, the land-use and land-cover data for Uttarakhand were obtained from the world land cover data having a spatial resolution of 1 km. The world land cover data are generated by the Joint Research Center of the European Commission (Mahmood & Babel, 2014; Mahmood & Jia, 2016).

Table 2 shows the weighted CN of different sub-watersheds of Uttarakhand.

The area-wise distribution of different soil groups is shown in Table 3.

The maximum area in Uttarakhand is covered by green forest (approximately 62%), followed by snow cover. Also, for the time period 2007–2014, there are no significant changes occurring in the land-use pattern of the state, as shown in Figure 11.

Table 4 shows the value of Ia, S, impervious layer of different sub-watersheds of Uttarakhand.

There are many methods incorporated in HEC-HMS mode to calculate the discharge at the outlet of the basin. The Muskingum method is simple in nature and does not need many inputs to calculate the outflow. Hence, in this paper, the Muskingum method is used to predict discharge at outlets of different sub-watersheds. Two major inputs in the Muskingum method are K and X values for the channel reach. The value of K can be calculated by the length of reach divided by the average flow velocity (Song et al., 2011; Scharffenberg, 2013; Sok & Oeurng, 2016), and the value of X is between 0.0 and 0.5. If X is 0.5, that means the channel has minimum attenuation and, if 0.0, the maximum attenuation. A value of 0.2 was used to specify the X value in the Muskingum method in different sub-catchments of Uttarakhand (Ntoanidis et al., 2013; Pichuka et al., 2017).

In the hydrological analysis, the calibration process is a process or methodology by which the model's input parameters are arranged in a way that the model's output in the form of simulated flow is able to match the nature of the observed flow. The historical flow at different time periods of two outlets, namely, Jyosimath (Chamoli) and Tehri, was used for calibration and validation. The use of different time periods for calibration and validation have the advantage that we can visualize both the calibration and validation period separately without any overlapping. As per the availability of data, the time period from 1985 to 2005 was selected as the calibration period for both the sub-watershed, and the period from 2006 to 2009 for Tehri and 2006 to 2010 for Jyosimath (Chamoli) as the validation period. These selections of the time period are based on the availability of data. The characteristics like land-use pattern and soil properties of the sub-watershed were counted as unchanged during the simulation period. In HEC-HMS, Nelder–Mead and the univariate gradient are the two principal algorithms available to optimize the objective function (HEC-HMS Tutorial, n.d.; Gebre, 2015; Wang et al., 2015). In this study, the peak weighted RMS error is taken and is minimized by utilizing the Nelder–Mead algorithm, and the optimal parameters of the model are estimated as shown in Tables 5 and 6. After optimization, the results obtained are used to validate the model.

As discussed earlier, all the inputs have the same source as used in the HEC-HMS model, except in the SWAT model, we use five-time series data. The outlet locations are the same as the HEC-HMS model. For the LULC and soil map, the map used in the HEC-HMS model is also used in the SWAT model. The detailed methodology and the working of the SWAT model are explained in Figure 12, where we can see that first we used the land-use property and soil map along with historical discharge data to calibrate the SWAT model, then the model is validated using some part of the historical discharge data.

In the present study, the Arc SWAT 2012 was used to develop the SWAT hydrological model for Uttarakhand. The advantage of using Arc SWAT 2012 is that it can manage and handle various spatial data sets easily (Sok & Oeurng, 2016). There are five main steps in SWAT modeling, delineation of the catchment, HRU definition, calibration, and validation of the model, and finally, sensitivity analysis of different parameters. In the present study, the Uttarakhand catchment is delineated into three sub-basins based on the outlet locations. These three sub-basins were further divided into four HRUs. The HRUs' threshold values of the land use, soil, and slope were fixed at 3, 3, and 5%, respectively, which are suggested by many authors. The first two years (1989–1990) are considered as a starting period to initiate the hydrological phenomena. The Hargreaves method, which needs only the precipitation and temperature data for the estimation of evapotranspiration (ET) is used in this simulation. The CN method is used for estimating surface runoff (SURQ), and the variable storage method is used for predicting streamflow (SF) routing. After model preparation, further analysis of the model, like parameter sensitivity analysis, calibration, and validation, was conducted with the help of the SWAT-CUP tool, which is a part of the SWAT model. In SWAT the global sensitivity analysis method is incorporated for evaluation of different parameters related to discharge, so this method was used to assess the most significant parameters for monthly streamflow simulations in the different sub-watersheds. In most of the previous work, a sequential uncertainty fitting algorithm (SUFI-2) was used for calibration of the model due to its simplicity and perfection. In the present analysis, we also used the SUFI-2, with 20 different parameters' unions (with 500 iterations) for the period 1991–2003 for Tehri and 1991–2005 for Jyosimath. In a particular repetition, the SUFI-2 includes the goodness of fit and the 95% prediction uncertainty (95PPU) between simulated and observed discharge. The new parameters obtained after the completion of each iteration, which has different ranges, are used to re-calibrate the model until the best parameters' series were received. These best parameters, thus obtained, were then applied to validate the monthly streamflow from 2003 to 2010.

The process of measurement of the rate of variation in model output with respect to model input for pre-defined parameters is called the sensitivity analysis. Generally, we do a sensitivity analysis to identify the required parameters and the accuracy needed for calibration (Von Storch & Zwiers, 1999; Goyal, 2014). The model calibration is defined as the method of determining model parameters by correlating model output with observed data for a given set of conditions. These outputs obtained during the calibration process are further used as input to run the model which is known as validation of the model. The following methods are used to evaluate the sensitivity of HMS and SWAT model results.

Nash–Sutcliffe efficiencies have the range from −∞ to 1. An efficiency of 1 (E = 1) shows a perfect match between model and observations. An efficiency of 0 indicates that the model predictions are as accurate as of the mean of the observed data, whereas an efficiency less than zero (−∞ < E < 0) occurs when the observed mean is a better predictor than the model. If the model efficiency is close to 1, that means the model is accurate.

Figure 13 shows the variation of the objective function, i.e., peak-weighted RMS error in the HMS model, which terminates after 32 iterations of the model run, with a function value of 490. Tables 5 and 6 represent the initial and optimized value of different parameters used in the HMS model during the calibration of the model for Uttarkashi and Tehri, respectively. It can be observed that there is not much difference between the initial and final value of the optimized parameter, and the objective function sensitivity is approaching zero.

The discharge obtained by using HMS and SWAT hydrological models are presented in the form of scatter-plots (Figures 14–17) for both calibration and validation periods. Tables 6 and 7 show the optimized value of different parameters during the HEC-HMS calibration process. Table 6 shows the optimized parameter value during the calibration process in the SWAT model. At both Uttarkashi and Tehri, where the discharge data are available, the orientation of observed flow was well paired by patterns of simulated flow during the calibration period. However, the peak flows during rainy seasons were not reproduced well by the SWAT model. At Tehri, in most of the years, the peak flows were underestimated by the model by 22% during calibration and by 33.1% during validation. However, the results are different in the case of Uttarkashi. The primary cause of the underestimation of peak flow is that there is a significant contribution of snowfall discharge to the rivers and the low density of the rain gauge stations in the state. The over- and underestimation of flow also affects the study of the effect of climate change on future estimated discharge. For example, at Uttarkashi, the model overestimated the flow by 17%. This means that if the projected increase at Uttarkashi is 30%, then the actual projected increase will be 13%. In the case of validation, there is more underestimation of flow at Tehri (33%), whereas at Uttarkashi, the overestimation rate is the same as calibration results in the HMS model. Overall sensitivity analysis shows that the calibration and validation results were satisfactory. However, the SWAT model shows better satisfactory results in comparison to the HMS model. The SWAT model could estimate the peak discharge better than the HEC-HMS model. The sensitivity analysis also verifies the better performance of the SWAT hydrological model.

The values, as shown in Figures 14–17 vary from 0.81 to 0.76 for the SWAT model and 0.78 to 0.73 for the HMS model, i.e., the SWAT model has better prediction capability than the HMS model. The SWAT model has five different hydrological components in comparison with the HMS model. As a result, the model represents the physical processes in a better way than the HMS model.

A comparative analysis was performed for both the models based on the parameters to check the statistical soundness. Table 8 shows the reference value, which is based on the previous study. The analysis shows that both the models are statistically sound for the modeling of discharge in the study area. In Tehri, the SWAT model is more satisfactory then the HMS model, as indicated by the PBIAS values of the model in Table 7.

For visualizing the impact of climate change on the discharge pattern at the Chamoli basin, the downscaled rainfall data under RCP4.5 and RCP8.5 have been used as input in the calibrated HEC-HMS model. The one limitation in projection is that we only considered the change in rainfall pattern and not on the other input parameters of the model. The predicted discharge pattern is shown in Figure 18. The prediction of peak discharge shows an increase of 27% in the case of RCP4.5 and an increase of 48% in the case of RCP8.5. The results, as shown in Figure 18, verify the rainfall prediction for the catchment under RCP8.5 and RCP4.5.

The primary objective of this study was to assess the performance of the HEC-HMS and SWAT models in simulating streamflow to find a suitable model for hydrological modeling in the Alkhnanda-Bhagirathi river basin. Calibration and validation for each catchment were performed. Based on the previous study, we found that the coefficient related to soil moisture storage and the CN are the most sensitive parameters for simulation of runoff. As such, precision has been maintained while evaluating these parameters. Both the models could simulate the monthly streamflow at the outlet of the catchment, i.e., Tehri and Chamoli. However, there is an underestimation in predicting the peak flow during the monsoon season, especially in July and August. However, the results obtained during the lean period or the period other than these two months are satisfactory in nature and acceptable. The statistical study shows that the SWAT model performs better than the HMS model, i.e., the comparative analysis showed the SWAT model to be slightly better at predicting the overall change in streamflow. The future prediction shows that the peak discharge of Alkhnanda may be increased by 27 and 47% under RCP4.5 and RCP8.5, respectively. One of the significant drawbacks of this study is that we did not consider the snowmelt data due to unavailability of the data. Moreover, we did not consider the change in LULC over time. As such, we also recommend further study considering snowmelt data and LULC change patterns.

The flow data have been provided by CWC, Lucknow.

Data cannot be made publicly available; readers should contact the corresponding author for details.