1 Introduction

Hydrologic models are essential to understand hydrological processes, for quantifying interactions among natural physical factors (Boorman and Sefton 1997) and assessing management strategies (Loukas et al. 2007). However, even with same input data, different hydrological models that are applied in same study basin might generate dissimilar simulation results (Jiang et al. 2007; Maurer et al. 2010; Shi et al. 2011). Even the correct global trend can be attained together, but otherness was still existed in the processes and spatial interactions (Ferrant et al. 2011). So the model selection becomes a priority for a successful hydrology studies in special region (Maurer et al. 2010; Vansteenkiste et al. 2013).

In practice, model selection is limited by common practice of modellers; it is rare that an objective model selection approach is authentic (Najafi et al. 2011; Nasr et al. 2007), which causes trouble for water resources decision makers in continuing water resource management (Wood 2004). In the phase2 of the Distributed Model Intercomparison Project (DMIP2) (Smith et al. 2012, 2013), several distributed models was mixed compared to lumped benchmark, results indicated there was none single model can perform best in all cases. DANUBIA component (Barthel et al. 2012) comprised of 17 model components and discussed the integrated simulation of global change influence on agriculture and groundwater. Najafi et al. (2011) farther suggested that joint application of different models would be significant for water resource management after using four hydrologic models in Tualatin River basin. These previous studies emphasized the importance of model selection and combined application of multiple models, however, there is little focus on the model performance of spatial hydrological processes, but only the discharges, and did not give detail interpretations of output deviations among different models.

The semi-distributed SWAT (Arnold et al. 1998) and fully distributed European Hydrological System (MIKE SHE) (Abbott et al. 1986), are two physically based hydrological models, their simulation including all the water balance components of a watershed. These two models have many good applications for hydrology studies (Kaini et al. 2012; Najafi et al. 2011; Thampi et al. 2010; Thompson et al. 2004). Fontaine et al. (2002) incorporated a modified snowfall-snowmelt routine through elevation bands in SWAT. MIKE SHE is a deterministic hydrological model whose input data and parameters are independent in each grid, and spatial heterogeneities of study basin can be described detailedly. Therefore, both SWAT and MIKE SHE are capable to be applied in alpine basin (Ahl et al. 2008; Debele et al. 2009; Liu et al. 2011, 2012; Rahman et al. 2012; Smerdon et al. 2009).

Since SWAT and MIKE SHE can fulfil their modelling tasks independently, there are only small numbers of studies intended to quantify the differences between two of them. Furthermore, most attentions focused on the goodness of fit indices for the modelled discharges at outlet (El-Nasr et al. 2005; Golmohammadi et al. 2014). While, the calibration solely at basin outlets alone and ignoring other hydrological components was not able to greatly improve model’s reliability and accuracy (Smith et al. 2013).

Mountainous watersheds, headstream of most river basins in arid region, play an important role in water resource management for the downstream region (Rahman et al. 2012). Tarim River basin, which is the longest inland river in world, is located in Xinjiang Province in northwest China. The limited water resources have severely affected sustainable development of this region and caused a vulnerable ecological environment (Chen et al. 2006; Liu et al. 2011). Yarkant River is the largest tributaries and primary water sources of Tarim River. Since most studies of Yarkant River basin have focused on the single effect of snow and glacier melt variation (Chen et al. 2006, 2010; Gao et al. 2010; Zhang et al. 2008), it would be meaningful to perform an integrated modelling study and understand the simulation effect in such region.

The studies of hydrological processes are essential to having the reasonable estimations for water resources utilization and management. In this study, because of extreme topographical condition, there is a strong spatial heterogeneity in Yarkant River basin. Two popular distributed models SWAT and MIKE SHE are applied to assess model response in arid alpine region. The spatiotemporal distributions of runoff, snowmelt and evapotranspiration form two models are used to explain effects of modules’ structure and algorithm on hydrological processes, and reinforce the understanding water cycle processes.

2 Study Area and Models

2.1 Study Area

Yarkant River (Fig. 1) originated from the north slope of Karakoram Mountains and forms a 2.5 × 104 km2 oasis along the lower reach below Kaqun station. This oasis is an important agroeconomic zone and cotton production region, and the largest irrigation region in Xinjiang. A region above Kaqun Station (Fig. 1) was selected as study catchment with an area of 50,248 km2, and a main stream length of 585 km. The distribution of precipitation and temperature is strongly uneven in this study region (Kang and He 1991; Yang 1989). The only one internal metrological station named Tashkurgan (Fig. 1) shows the mean annual precipitation of 95 mm and pan evaporation of 1500 mm during 2000–2009, and the average temperature was below zero from November to next March.

Fig. 1
figure 1

Location of Yarkant river basin and hydrologic station, meteorological stations and channel network

Full size image

There are abundant glaciers in this study region of Yarkant River with a total cover area of 5574 km2, the estimated amount of ice storage is 685 km3 (Yang 1989). The melt water provides rich water resources to Yarkant River, with a mean annual discharge of 6.87 × 1010 m3 detected at Kaqun station. Meanwhile, because of seasonal snow and ice melt, temporal distribution of discharge varies seasonally, and proportion of discharge from June to September was approximately 80 % of total discharge.

The study catchment has very complex terrain with extreme variations in elevation gradient: the altitude from 8611 m down to 1450 m with an average elevation of 4450 m and slope of 30.78°. The alpine meadow and snow-ice were dominant land cover types having proportions of 29.63 and 25.89 %, respectively. Dystric cambisols, lithic leptosols and haplic chernozems are three prominent soil types and accounts for 36.55, 17.46 and 12.67 % of the basin area, respectively.

2.2 Hydrological Models

Considering the modelling approach, both SWAT and MIKE SHE has their unique features from model structure to governing equations. In SWAT, the catchment is subdivided into sub-basins that are connected with the river network; the sub-basin is further subdivided into specific soil/land use/slope characteristic units called Hydrologic Response Units (HRUs). HRUs are the smallest computed unit without any spatial relationship or interaction among them. Semi-empirical Soil Conservation Service (SCS) approach is used to calculate surface flow. After runoff contribution to channels is obtained through water balance control, the variable storage approach or Muskingum methods can be used to compute channel flow. Elevation bands are the important configuration for snowmelt simulation in the SWAT model. Catchment in SWAT can be split into ten elevation bands maximally, within precipitation and temperature adjusted regarding elevation gradient change coefficients.

In MIKE SHE, catchment is split into a number of square grids with explicit geo-information in horizontal direction. In vertical direction, the vertical variations of soil and hydrogeological characteristics are assigned in number of layers with variable depths. Infiltration is calculated depends on the groundwater conditions and the soil physical attributes which parameters have a clear physical meaning. When water enters unsaturated soil zone, one-dimensional Richards equation or Gravity flow equation can be used to calculate the vertical soil interflow. In saturation zone, three-dimensional Darcy equation is employed with an iterative implicit finite difference technique. A dynamic coupling is set between MIKE SHE and MIKE 11 (Thompson et al. 2004) to simulate the dynamic interaction between aquifers and stream channels. For the snowmelt module, a threshold temperature is assigned to determine the beginning time of snowfall and snowmelt.

3 Methodology

3.1 Input Data

Digital Elevation Model (DEM) data with the resolution of 90 m was used to obtain physical attributes of catchment terrain and river network. The catchment was divided into 25 sub-basins (Fig. 1) in SWAT and simulation resolution of MIKE SHE is 2 km × 2 km. The land use types were determined by the land use data of Global Land Cover Network (GLCN) with resolution of 500 m × 500 m. The remote sensing data images from Moderate Resolution Imaging Spectroradiometer (MODIS) were used to achieve the related parameters such as the leaf area index, root depth, crop coefficient, and growing period of each vegetation type. The properties of the different soil types were based on the data of Harmonized World Soil Database (HWSD).

From China Meteorological Data Sharing Service System, daily data including precipitation, maximum/minimum/mean temperature, wind speed, relative humidity and solar radiation at Tashkurgan station were collected from 2000 to 2009. In SWAT model, precipitation and maximum/minimum temperatures in each elevation band were interpolated from the precipitation lapse rate (PLAPS) and temperature lapse rate (TLAPS). In MIKE SHE model, the catchment was divided into ten regions according the elevation of each grid at 700-meter intervals, every grid in one region shares the same precipitation and temperature time series, which got through the same PLAPS and TLAPS values set as 76 mm/km and -7 °C/km, respectively, those are based on some investigates about the amount of precipitation and temperature in the different altitude districts (Gao et al. 2010; Jing 2010; Wang et al. 2009).

3.2 Calibration

The simulated period was parted warm-up period 2000–2002, calibration period 2003–2007 and verification period 2008–2009. The auto-calibration modules, SWAT-CUP of SWAT and Auto Calibration Tool (ACT) in MIKE SHE package, have been used to calibrate parameters to improve the calibration efficiency. Sequential Uncertainty Fitting (SUFI-2) (Abbaspour et al. 2004) and Shuffled Complex Evolution approach (SCE) (Vrugt et al. 2003) were used in SWAT-CUP and ACT. The simulated daily discharge was calibrated by the objective function according the observation record at Kaqun station. Four statistical coefficients were used to determine the model performances: Nash-Sutcliffe efficiency coefficient (Nush and Sutcliffe 1970) NSE, Pearson correlation coefficient R, root-mean-square error RMSE and percentage bias PBIAS. Their formulations are written as:

$$ NSE=1-\frac{{\displaystyle {\sum}_{i=1}^n}{\left({Q}_{obs,i}-{Q}_{sim,i}\right)}^2}{{\displaystyle {\sum}_{i=1}^n}{\left({Q}_{obs,i}-{\overline{Q}}_{sim,i}\right)}^2} $$
(1)
$$ R=\frac{{\displaystyle {\sum}_{i=1}^n}\left({Q}_{obs,i}-{\overline{Q}}_{obs,i}\right)\left({Q}_{sim,i}-{\overline{Q}}_{sim,i}\right)}{\sqrt{{\displaystyle {\sum}_{i=1}^n}{\left({Q}_{obs,i}-{\overline{Q}}_{obs,i}\right)}^2}\sqrt{{\displaystyle {\sum}_{i=1}^n}{\left({Q}_{sim,i}-{\overline{Q}}_{sim,i}\right)}^2}} $$
(2)
$$ RMSE=\sqrt{\frac{{\displaystyle {\sum}_{i=1}^n}{\left({Q}_{sim,i}-{Q}_{obs,i}\right)}^2}{n}} $$
(3)
$$ PBIAS=\frac{{\displaystyle {\sum}_{i=1}^n}\left({Q}_{sim,i}-{Q}_{obs,i}\right)}{{\displaystyle {\sum}_{i=1}^n}{Q}_{obs,i}} $$
(4)

where Q obs,i and Q sim,i are the measured and simulated discharges at ith day (m3/s), respectively; \( {\overline{Q}}_{obs} \) and \( {\overline{Q}}_{sim} \) are the average measured and simulated discharges in the simulation period (m3/s), respectively; and n is the value of the time steps.

4 Results and Discussion

4.1 Simulation Results

Figure 2 presents the comparison of simulated and observed discharges at Kaqun station in 2003–2009 at a daily temporal scale. The discharge hydrography curves from SWAT and MIKE SHE caught the overall trend of rising and recession and reflected the temporal characters of stream flow at outlet. According to the evaluation criteria of calibration and verification, both the applications of two models in Yarkant River basin obtained acceptable performances. SWAT obtained slightly better indices (NSE = 0.76, R 2 = 0.78, RMSE = 144.51, PBAIS = 0.22 %) than MIKE SHE (NSE = 0.71, R 2 = 0.70, RMSE = 169.88, PBAIS = 8.12 %).

Fig. 2
figure 2

The observed and simulated daily discharges at Kaqun station

Full size image

The criteria of discharge calibration are not the overall referential measurements to evaluate the accuracy of hydrologic model application (Vázquez and Feyen 2007), but response of a model to natural hydrologic process can be analysed through the variation of different hydrological elements. Table 1 provides the quantified results of the average annual values of different hydrological components in the whole basin. ETa in SWAT and MIKE SHE was furtherly divided into four parts according to their sources: snow sublimation, canopy interception, river/pound water evaporation and soil evaporation/transpiration. From Table 1, it can be seen there are some notable differences in runoff, snow storage and ETa between SWAT and MIKE SHE. It would be worth to exam the differences in more detailed view regarding the different structures and arithmetic.

Table 1 The average annual values of the water components from SWAT and MIKE SHE inYarkant River basin in 2003–2009
Full size table

4.2 Runoff

Both in SWAT and MIKE SHE, channel flow was derived from surface runoff and base flow, but the constitutions of base flow are quite distinct in two models. In SWAT, the average annual base flow was 82.6 mm which is the sum of the groundwater runoff (16.3 mm) and lateral subsurface runoff (66.5 mm). While, the base flow in MIKE SHE is groundwater runoff that equals to 31.1 mm. The constitution and distribution of runoff in daily scale is presented in Fig. 3 below. Previous studies reported that the rapid lateral subsurface flow provides a dominant contribution to the storm flow in headwater catchments (Kienzler and Naef 2008; Verseveld et al. 2009; Swarowsky et al. 2012). Especially in alpine area, most of ice-snowmelt water infiltrates into soil and contributes to stream as subsurface lateral flow (Fan et al. 2014) . This feature is described reasonably in SWAT model with subsurface runoff contribution rate of 41.4 %, what is anymore, in SWAT, base flow present an obvious seasonality, and little contribution to stream from the November to the following March (Fig. 3).

Fig. 3
figure 3

The daily constitution of the simulation runoff from SWAT and MIKE SHE in Yarkant River basin

Full size image

In MIKE SHE, Richards and Gravity flow equation only take vertical flow into account in unsaturated zone, but when soil moisture is saturated, and head pressure of aquifer water is higher than surface, the infiltrated soil water could flow back to surface, and this part accounted for 39.4 % of total stream runoff in this simulation (close to the subsurface lateral flow of 41.4 % in SWAT). As well as the appropriate soil water recharge to aquifer, MIKE SHE got a persistent and steady contribution to stream in the low water period as base flow (Fig. 3). The proportion of 21.3 % agreed very well with the result of 23 % from Fan et al. (2013), which was obtained through a multiple base flow separate approach.

4.3 Snowmelt

4.3.1 Temporal Distribution of Snowmelt

The simulated daily snowmelt in Yarkant River basin from 2003 to 2009 is illustrated in Fig. 4. SWAT resulted in more concentrated snowmelt in temporal distribution; the amount of snowmelt from June to September constituted 76.84 % of the annual snowmelt. While in MIKE SHE, the period of snowmelt is longer and the temporal distribution was more dispersive, and amount of snowmelt from June to September constituted 58.43 %.

Fig. 4
figure 4

The simulated daily snowmelt from SWAT and MIKE SHE in Yarkant River basin

Full size image

Degree-day approach is employed in both SWAT and MIKE SHE model to calculate snowmelt. In SWAT, variation of snowmelt factor and snow pack temperature that reflected accumulated temperature of land surface are considered, this is possibly dominant cause for the more concentrated distribution in summer, also the simulation process would be more flexible for catching natural realities rather than only considering air temperature (Zuzel and Cox 1975) in MIKE SHE. Consequently, the temporal distribution of snowmelt in SWAT model matched better with the runoff characteristics of Yarkant River basin, in which the greatest proportion of the water source derives from the melting snow and glaciers, and amount of runoff from June to September accounts for 80 % of the total annual runoff (Chen et al. 2006, 2010).

4.3.2 Spatial Distribution of Snow Storage

The average annual depth of snow storage from SWAT and MIKE SHE were compared among different elevation bands in Fig. 5, the ones above 7100 m that have very small area were not demonstrated. From Fig. 5, the consistent tendency of snow storage in different elevation bands from two models can be found. Increasing storage occurs above 5000 m in both model and the increasing amount was prominent with the higher elevation.

Fig. 5
figure 5

The snow storage in the different elevation bands from SWAT and MIKE SHE in Yarkant River basin in 2003–2009

Full size image

The distribution of elevation bands was rough and cannot illustrate the spatial status of snowpack; hence that, the fine distribution based on model’s resolution is shown in Fig. 6. In the first column, the results are from SWAT and are distributed as a sub-basin; in the second column, the results from MIKE SHE are averaged in each sub-basin of SWAT; and in the third column, the results from MIKE SHE are displayed with original grid view. According to first and second column in Fig. 6, two models have a similar tendency for spatial change during the different seasons: till the end of March, almost the entire catchment is covered by snow; snowmelt first occurs around the outlet in April, the snow melt mainly extends to the mountain region approximately 4500 m in July; and a new increasing storage appears after October; but in the region in which elevation is greater than 5500 m, the snow cover shows a continuously increasing trend throughout the entire year. However, between the first and third column, the spatial distribution of snow cover in MIKE SHE is much more distinct than that in SWAT, and differences of the snow storage depth expressed by the grids of MIKE SHE are much more obvious than by the sub-basins of SWAT.

Fig. 6
figure 6

The simulated snowpack of SWAT (1st column), MIKE SHE in sub-basin (2nd column) and MIKE SHE (3rd column) in Yarkant River basin on 31st Mar. (1st row), 30th Jun. (2nd row), 30th Sep. (3rd row) and 31st Dec. (4th row) 2003

Full size image

For the differences of snowpack, the effect of model structure cannot be detached. In SWAT, HRU is the smallest calculated unit, but spatial relation only exists between sub-basins; therefore, an areal depletion curve was introduced when considering the unequal distribution in a sub-basin and makes a uniform spatial distribution of snow storage. As a fully distributed model, in MIKE SHE model, each grid can reflect the actual spatial change of the snow storage in that place, so spatial distribution is more distinct.

4.4 Evapotranspiration

Basin on Table 1, the differences of ETa between SWAT and MIKE SHE were greatly reflected on river/pound water evaporation and soil evaporation/transpiration, while the other values are similar. Subsequently, the annual ETa in MIKE SHE was 23.5 mm more than that in SWAT. For the difference in river/pound water evaporation, it is possible that the spatial relationship between calculated units is one cause. In SWAT, HRUs are controlled by water balance independently, but the grids in MIKE SHE are linked by hydrodynamic relation, and the flood detention area is considered, what the result is that larger area of water surface in MIKE SHE model may cause more water evaporation.

The soil evaporation in SWAT is limited by soil moisture content. When soil moisture content is under field capacity, evaporation will be reduced according to their difference. In addition to restrict the amount of soil water evaporation under dry conditions, SWAT also defines a maximum value of soil evaporation as 80 % of the plant available water on a given day. However, in MIKE SHE, without any limitation, all of available water in soil will be used to meet soil evaporation. Therefore, there could be more soil evaporation in MIKE SHE.

From the aspect of spatial distribution based on third column of Fig. 7, the great influence of land use on ETa has been clearly demonstrated. Furthermore, comparing the first and second column under the same resolution, the tendency of spatial change is similar in different seasons between SWAT and MIKE SHE: the evapotranspiration around the outlet is much higher than that of the mountain area, because of higher temperature and better plants cover. This spatial difference is much more obvious in the plant grow summer season. However, comparing the first and third column, the semi-distributed SWAT model cannot support the specific information about the spatial distribution.

Fig. 7
figure 7

The daily average evapotranspiration of SWAT (1st column), MIKE SHE in sub-basin (2nd column) and MIKE SHE (3rd column) in Yarkant River from Jan. ~ Mar. (first line), Apr. ~ Jun. (second line), Jul. ~ Sep. (third line), and Oct. ~ Dec (forth line) 2003

Full size image

5 Summary and Conclusions

Based on the statistic evaluation indices of discharges calibration on Kaqun station, both SWAT and MIKE SHE got acceptable performances in Yarkant River basin. Because there is no more observations can be used to calibrate models, it would be reasonable to cautiously draw the conclusions for model responses to hydrological processes. Therefore, the differences of multiple hydrological components from two models were intercompared based on modules’ structures and algorithms, the main consequences include the following aspects.

The phenomenon that most of ice-snowmelt water will infiltrate into the soil and contribute to the stream as subsurface runoff have been described by SWAT, with a contribution proportions of 41.4 %. While, MIKE SHE generated a more reasonable base flow with a contribution of 21.3 % also ignoring the soil lateral flow. SWAT obtained a better snowmelt process corresponding to the character of the river flow, but cannot distinctly reflect the spatial features of snowpack that could be detailedly achieved by MIKE SHE. The less ETa in SWAT has been mainly caused by less water surface area and restrictive setup of soil evaporation function. In comparison, MIKE SHE also provided more information regarding the evapotranspiration being closely related to land use.

From the careful quantifying and reasoning mentioned above, the application of the SWAT and MIKE SHE models in Yarkant River basin can agree with natural observations in some aspects but not in the entire cycle processes. Combined the outputs from two models, an improved understanding of hydrological processes can be presented: SWAT complements the subsurface lateral flow to MIKE SHE which had a better groundwater simulation in this mountainous region; MIKE SHE can supplement the spatial distribution of snowpack and ETa for SWAT’s output. Moreover, the interpretations of variations between two models stated how the structure and algorithm impact on hydrological processes and provided an inspiring reference to hydrological processes study in other catchments with unique features. In future study on arid and scarcely gauged alpine basin, the joint application of multiple hydrologic models and combined results could be an effective way to control the uncertainties from modules’ structure and algorithm, and improving understanding of hydrological processes.