Baseflow simulation using SWAT model in an inland river basin in Tianshan Mountains , Northwest China

Baseflow is an important component in hydrological modeling. The complex streamflow recession process complicates the baseflow simulation. In order to simulate the snow and/or glacier melt dominated streamflow receding quickly during the high-flow period but very slowly during the low-flow period in rivers in arid and cold northwest China, the current one-reservoir baseflow approach in SWAT (Soil Water Assessment Tool) model was extended by adding a slowreacting reservoir and applying it to the Manas River basin in the Tianshan Mountains. Meanwhile, a digital filter program was employed to separate baseflow from streamflow records for comparisons. Results indicated that the tworeservoir method yielded much better results than the onereservoir one in reproducing streamflow processes, and the low-flow estimation was improved markedly. Nash-Sutcliff efficiency values at the calibration and validation stages are 0.68 and 0.62 for the one-reservoir case, and 0.76 and 0.69 for the two-reservoir case. The filter-based method estimated the baseflow index as 0.60, while the model-based as 0.45. The filter-based baseflow responded almost immediately to surface runoff occurrence at onset of rising limb, while the model-based responded with a delay. In consideration of watershed surface storage retention and soil freezing/thawing effects on infiltration and recharge during initial snowmelt season, a delay response is considered to be more reasonable. However, a more detailed description of freezing/thawing processes should be included in soil modules so as to determine recharge to aquifer during these processes, and thus an accurate onset point of rising limb of the simulated baseflow.


Introduction
Baseflow is a streamflow component which reacts slowly to rainfall and is usually associated with water discharged from groundwater storage (Eckhardt, 2008).Knowledge about baseflow is useful in assessing water quality, forecasting streamflow, allocating water supply, and designing hydropower plants (Tallaksen, 1995) under low-flow conditions.When, where, and how much streamflow can be attributed to groundwater discharge is thus practically important.
Baseflow is, therefore, an important component in hydrological simulation.Conceptual modeling of baseflow usually assumes that outflow from the aquifer is linearly proportional to its storage (Aizen et al., 2000;Fenicia et al., 2006;Eckhardt, 2008;Ferket et al., 2010), sometimes combined with analytical solutions of the simplified Boussinesq equation (Paniconi et al., 2003;Troch et al., 2004;Hilberts et al., 2004).Wittenberg (1999) argued that the unconfined aquifer is unlikely a linear reservoir, instead, more likely a non-linear one.However, Fenicia et al. (2006) confirmed that the linear storage-discharge relationship describes groundwater behavior best.Baseflow itself may be composed of a number of components, each of which may vary seasonally with different recession constants (Nathan and McMhon, 1990).As probably a compromise, multi-reservoir algorithms, linear, non-linear, or combined were used to generate baseflow by e.g.Tallaksen (1995), Ferket et al. (2010), and Samuel et al. (2011).

Y. Luo et al.: Baseflow simulation using SWAT model in an inland river basin in Tianshan Mountains
The Soil and Water Assessment Tool (SWAT, Arnold and Fohrer, 2005) uses a conceptual linear one-reservoir (shallow aquifer storage) approach to simulate baseflow.SWAT partitions groundwater into two aquifer systems: a shallow aquifer which contributes baseflow to streams within the watershed, and a deep aquifer which contributes baseflow to streams outside the watershed and can be considered lost from the system (Arnold et al., 1993).While the shallow aquifer-baseflow was properly reproduced by SWAT (Arnold et al., 2000;Jha et al., 2007), weaker simulation of baseflow was found as well (Kalin and Hantush, 2006;Srivastva et al., 2006).Peterson and Hamlett (1998) found that SWAT was not able to simulate baseflow due to the presence of soil fragipans.Chu and Shirmohammadi (2004) found that the baseflow was not simulated properly for an extremely wet year.Wu and Johnston (2007) found underestimated baseflow by SWAT especially during dry years in a Great Lake watershed and indicated that this is primarily due to the long temporal lag between winter snowpack accumulation and spring snow melting events.Luo et al. (this study) found underestimation of baseflow during the low-flow period in the Manas River Basin in northern Tianshan Mountains using SWAT2005.For the glaciated Oigaing River basin in western Tianshan and the Ala Archa River basin also in northern Tianshan, Aizen et al. (2000) used one linear reservoir to generate baseflow and found that the discharge was underestimated during autumn-winter as well.The steep slopes of the river basins in Tianshan Mountains, the quick recession of surface runoff, and the sluggish and stable baseflow processes might indicate a quick percolation of rainfall and snow and glacier melt waters during the summertime to an underground storage which releases slowly during the wintertime.Nathan and McMhon (1990) indicated that baseflow itself may be composed of a number of components, each of which may vary seasonally with different recession constants.Supposedly, an additional slow release pool may improve the low-flow estimation for rivers in the Tianshan Mountains, which is not present in both Aizen's (2000) model and the SWAT model.
Therefore, the main objective of this paper is to develop a two-reservoir approach for baseflow simulation in SWAT and use the model to simulate the streamflow process that is characterized by combined steep and sluggish recession stages of the receding limb in an arid and cold inland river basin in the Tianshan Mountains, northwest China.

Baseflow modeling in SWAT
In the SWAT model, water routed through channel system to the gauges consists of four components: direct surface runoff (Q sf ), lateral flow from unsaturated soil profiles (Q lt ), drainage from tiles (Q tl ), and baseflow from underground storage (Q b ) (Fig. 1).Modeling of the direct surface runoff, the lateral soil flow, and the tile drainage are described in detail in theoretical documents of SWAT model (Neitsch et al., 2005) and thus will not be described repeatedly here.The baseflow simulation will be focused hereafter.
SWAT differentiates the underground storage into two portions, shallow aquifer and deep aquifer.The shallow aquifer receives recharge from the unsaturated soil profile percolation.An exponential decay weighting function is utilized to account for the time delay in aquifer recharge once the water exits the soil profile (Neitsch et al., 2005).The delay function accommodates situations where the recharge from the soil zone to the aquifer is not instantaneous, i.e. 1 day or less.The recharge to aquifer on a given day is calculated as below: where W rchrg is the amount of recharge entering the aquifers (mm H 2 O day −1 ), δ gw,sh is the delay time of the overlying geologic formations (days), W seep is the total amount of water exiting the bottom of the soil profile (mm H 2 O day −1 ); subscriptions "seep" indicates seepage water exiting bottom of unsaturated soil profile, "rchrg" indicates recharge, i is the sequential number of days, and "sh" indicates the shallow aquifer storage.
A fraction of the total daily recharge can be routed to the deep aquifer.The amount of water diverted from the shallow aquifer due to percolation to the deep aquifer on a given day is given by: where β dp is a coefficient of shallow aquifer percolation to deep aquifer, and subscription "dp" indicates deep aquifer.The amount of recharge entering the shallow aquifer is: Baseflow generated from the shallow aquifer on a given day i under influence of recharge is given as below (Neitsch et al., 2005): where Q b,sh,i is the baseflow from the shallow aquifer on day i (mm H 2 O day −1 ), and "b" indicates baseflow, and t is the step time length.Daily time step is used in this study.
When only one reservoir is used, the baseflow is equal to that from the shallow aquifer.
SWAT assumes that water entering the deep aquifer is not considered in the future water budget calculations and can be considered lost from the system (Neitsch et al., 2005).This study uses the deep aquifer as a parallel reservoir generating the baseflow, which enters the channel system eventually, to improve the streamflow process simulation in the low-flow period.When the two-reservoir approach is used, baseflow from the shallow aquifer is expressed as in Eq. ( 4), and following Eqs.
where W rchrg,dp is the amount of recharge entering the deep aquifer (mm H 2 O day −1 ), δ gw,dp is the delay time or drainage time of the deep aquifer geologic formations (days), W seep,dp is the total amount of water exiting the bottom of the shallow aquifer (mm H 2 O day −1 ), Q b,dp is baseflow component from deep aquifer.The total baseflow is then given as below: When the shallow storage reservoir is used only to generate baseflow, recharge to the deep aquifer is disabled.When both aquifers are used to generate baseflow, the parameter β dp is determined through calibration.Other parameters to be calibrated for baseflow modeling include the delay time δ gw,sh , δ gw,dp , the recession constants α gw,sh and α gw,dp .

Baseflow separation using automated digital filter
In consideration of difficulties in the measurement of baseflow, a third-party approach, the digital filter-based program, is used to separate baseflow from streamflow records for comparison purposes.This baseflow separation procedure is based on a recursive digital filter commonly used in signal analysis and processing (Lyne and Hollick, 1979).
It was used by Nathan and McMahon (1990), Arnold and Allen (1999) and Szilagyi et al. (2003Szilagyi et al. ( , 2004)), among others.This technique is, in fact, arbitrary and physically unrealistic.However, it does provide a subjective and repeatable estimate of baseflow that is easily automated (Nathan and McMahon, 1990).The filter given by Lyne and Hollick (1979) is expressed as below: where Q sf and i are defined as before, Q s is the surface runoff, and λ is the filter parameter.Baseflow is calculated as below: where Q b is defined as before.

Model setup
MRB is located at the northern side of the middle Tianshan Mountains, northwest China (Fig. 2).MRB originates from the Yilianhabierga Mountain, runs 160 km to the outlet at Kenswat Hydrological Station (KHS, 85 • 57 E, 43 • 58 ), and runs further 240 km through the oasis and the desert and finally merges into Manas Lake.The catchment area of the MRB above the outlet KHS is 5163 km 2 .
Maps of a 1:250 000 DEM, a 1:100 000 land cover, and the China Glacier Inventory (CGI) were used to setup the Arc-SWAT2005.The CGI data used as the initial glacier layout were mainly derived from topographical maps (1:100 000) based on aerial photos acquired during 1962-1977(Shangguan et al., 2009)).Eventually, the watershed is delineated into 27 subbasins and 163 Hydrological Response Units (HRUs).Each subbasin is divided into ten bands with equal elevation increment for simulating the snow and glaciers.within a subbasin are significant.On average, the difference is 2561 m, with the biggest at 3876 m for the subbasin 11 and the smallest at 448 m for the subbasin 1.The land cover types include range grassland of 8.11 % elevation band 2500-3500 m a.s.l., short bushes of 39.1 % within the elevation band of 1500-2500 m a.s.l., and forest of 5.32 % within the below the elevation of 1500 m a.s.l., bare land of 33.58 %, and glaciers.Among the 163 HRUs, there are 28 glacierized ones with total glacier area of 717 km 2 .Ratios of glacier area to subbasin area range from 0.7 % for the subbasin 4 to 51.2 % for the subbasin 22, with a mean ratio of 13.9 % over the basin.Watershed glacier processes simulation was detailed in Luo et al. (this study). (

2) Soils
The main soils in the basin include alpine meadow soil, subalpine meadow soil, subalpine meadow and steppe soil, mountain chernozem soil, mountain grey cinnamon soil, mountain chestnut soil, which take account of 36 %, 11 %, 42 %, 1 %, 7 %, and 3 % of the basin area, respectively.Textures and properties of these soils were derived from the field-collected and lab-tested data of the publication "Soils in Xinjiang (technical report)". (

3) Climate
The For each subbasin, a virtual weather station (VWS) is defined.For each VWS, the temperature and precipitation data were derived from the SWS by using the temperature and precipitation lapse rates.Default value −6 • C km −1 in the SWAT model was used for the temperature lapse rate,and 45 mmkm −1 was used for the precipitation lapse rate (Luo et al., this study).

(4) Streamflow
Daily streamflow records at the KHS from 1961-1999 were used.The mean daily discharge rate is 39.3 m 3 s −1 and the average annual volume 12.15 × 10 8 m 3 .The recorded maximum annual volume is 20.08 × 10 8 m 3 in 1999 and the minimum 9.39 × 10 8 m 3 in 1983.The flow volume from June to August takes account of 70.5 % of the annual value and of 28.9 % and 25.9 % for July and August, respectively.During the seven months from October to the next April, flow volume accounts for 15.9 % of the year with the monthly ratio decreasing from 4.2 % in September to 1.2 % in February and then going up gradually to 2.0 % in April.
MRB is snow and glacier melt dominated.Snowmelt starts usually in late April or early May, glacier melts as snowpack depletes, and streamflow starts to rise consistently till peak discharge in late July.Glacier melt contribution ceases in late September.As temperature falls below 0 • C, a new snow season begins, and the direct surface runoff to streamflow ceases.Steep rising and receding streamflow curve is then followed with an almost flat low-flow line during the winter and spring seasons, while the streamflow is very stable and has quite a long duration (Fig. 3), which is a common feature for rivers in northwest China.
where Q obs i is the i-th observation for the daily flow, Q sim i is the i-th simulation value for the daily flow, mean is the mean of observed data for the daily flow, and n is the total number of the daily flow observations.NSE indicates how well the plot of observed versus simulated data fits the 1:1 line.NSE ranges between −∞ and 1.0 (1 inclusive), with NSE = 1 being the optimal value.PBIAS measures the average tendency of the simulated data to be larger or smaller than their observed counterparts.The optimal value of PBIAS is 0.0, with low-magnitude values indicating accurate model simulation (Morasi et al., 2007).

Results
Parameters calibrated for baseflow components are listed in Table 1 for reference.In case one-reservoir was used only, it was assumed that water exiting the bottom of the unsaturated soil profiles recharged the shallow aquifer only.When the two-reservoir method was employed, it was found that 40 % of recharge to the deep aquifer was proper to match the measured streamflow during low-flow period.The deep aquifer has a much longer delay time for recharge and a much slower recession rate than the shallow aquifer.Parameters for the baseflow filter are also listed in Table 1.It is interesting to notice that using the baseflow days given by the filter program (Arnold et al., 1995;Arnold and Allen, 1999) as the recharge delay time of the deep aquifer in the two-reservoir approach simulated the stream flow during low-flow period very well.Figure 4 demonstrates the comparison of stream flows among those measured and simulated.A five-year clip was taken from the whole simulation period of 39 yr to give a clearer picture.The two-reservoir method improved the lowflow simulation remarkably in vision.Statistical indices NSE and PBIAS indicate that the two-reservoir method yielded "good" or even "very good" results in the sense of either NSE or PBIAS following the rating rules given by Moriasi et al. (2007), which are better than the one-reservoir method, Table 2.
Table 3 lists the summary statistics of measured and simulated streamflow volumes.Annual flow volumes simulated by SWAT using one-reservoir and two-reservoir methods are approximated with only minor differences.The simulated mean annual flow volumes are slightly larger than those measured, and the simulation was rated as very good, Table 2.However, significant differences were found between the simulated and the measured maximum flow volume.The maximum flow volume was observed in 1999, while simulated in 1988.The difference might be due to the uncertainty of meteorological input in mountain areas, which was derived from records of the base station at the foot of the mountain using a single precipitation lapse rate.Moriasi et al. (2007).

Discussion
The annually averaged baseflow volumes given by different approaches were listed in Table 4.The one-reservoir and two-reservoir approaches produced similar results in annual baseflow volumes.The digital filter program gave a much larger baseflow volume than the model-based methods.The model-based baseflow volume accounted for 45 % of the annual flow volume, while the filter-based 60 %.Among the model-based baseflow, shallow aquifer contributed 58 %, and the deep aquifer 42 %.According to the recharge partition coefficient, these should be 60 % and 40 %, respectively.The minor difference was attributed to a portion of the shallow aquifer storage depleted during the simulation period.For the deep aquifer, its storage fluctuated seasonally while equilibrium was maintained, as the simulation revealed.
The observed streamflow flattened out with delayed flow supply from deeper subsurface stores and eventually became nearly constant, which is sustained by outflow from groundwater storage (Fig. 3). Figure 4 illustrates the baseflow processes.When only one-reservoir was used, the base flow was underestimated as shown in Fig. 4. The observed mean daily flow rate was 9.8 m 3 s −1 and the volume of 1.88 × 10 8 m 3 during the period from October to April, while the simulated value was 3.0 m 3 s −1 and the volume of 0.58 × 10 8 m 3 for the same period.Similar results were also found for Oigaing River basin in western Tianshan and the Ala Archa River basin in northern Tianshan when a onereservoir baseflow approach was used (Aizen et al., 2000).The steep slopes of the MRB (Luo et al., this study), the quick recession of surface runoff, and the sluggish and stable baseflow processes during late autumn and late spring (Fig. 1) indicate a quick percolation of rainfall and snow and glacier melt waters during the summertime to an underground storage which releases slowly during the wintertime, as found by isotopic measurements in the Wind River Range of Wyoming of US (Cable et al., 2011).When the deep aquifer was employed as an additional slow release pool, the baseflow simulation was improved remarkably, as demonstrated in Figs. 5 and 6a and b.
The model-and filter-based daily baseflow processes were averaged over the period from 1966 to 1999 (Fig. 6).The simulated streamflow peak time seemed to shift a little earlier.Nevertheless, the simulated rising receding limbs match the measured ones well.
During the low-flow period (from November to April), it was noticed that the one-reservoir method gave a serious  underestimation of streamflow, while the two-reservoir method reproduced the streamflow properly.Combination of the quick release reservoir and the slow release reservoir matched both the quick ad sluggish receding stages of recession limb of the streamflow very well.During the quick receding stage, the quick-release pool played a more important part than the slow-release pool and vice-versa during the sluggish stage (Figs. 3 and 7).
The filter-based baseflow started to rise earlier, reached its peak later, and turned to low-flow stage earlier again than the model-based (Fig. 6).The earlier peak time of the modelbased baseflow might be attributed to streamflow peak time shift.
Onset points of rising limbs are worth noting (Figs. 4  and 6).The filter-based baseflow responds to runoff occurrence immediately at the onset of rising limb.The model-based baseflow responds with a delay.The Manas River basin is snow and glacier melt dominated.Snowmelt starts usually in the middle of April.Infiltration in frozen soils is affected by soil permeability, water content, repeated thawing and refreezing, and many other factors and their complex interactions (French and Binley, 2004;Stahli, 2005).Experimental (Hayashi et al., 2003;Iwata et al., 2010) and mathematical (Flerchinger and Saxton, 1989;Zhao and Gray, 1999) investigations revealed the impeding effects of frozen soil layer to snowmelt infiltration, and hence the potential recharge to aquifers.As a matter of fact, baseflow should respond with a delay to the snowmelt, other than immediately, as given by the filter-based approach (Fig. 6).This is achieved by the SWAT model through a simple assumption of no water flow during the frozen season.An issue remaining to be addressed is infiltration to and recharge from the soil profile during freezing and thawing, and eventual determination of the onset point of rising limb (Eqs.4 and 7).This needs more detailed description of soil freezing/thawing processes, which are described insufficiently in most watershed hydrological models.Both the model-based and the filter-based approaches captured the reflection point of recession limb in late September, when direct surface runoff usually ceases to discharge the aquifers.
Compared to the one-reservoir method, the two-reservoir method requires three extra parameters to calibrate.This case study revealed that the extra parameters will not provide too much extra work for calibration if proper steps are followed.Firstly, calibrate the parameters with recharge to deep aquifer disabled, and optimize parameters with the quick rising and receding limbs as target.Then, activate the recharge to the deep aquifer and optimize the three parameters for the slow release pool.An important reminder is that the slow-release pool has a longer recharge delay time and smaller recession constant than the quick-release pool.And, as aforementioned, the baseflow days given by the filter-based program can be a tip for calibrating the delay time.The recession constant given by the filter program can also be an initial estimation of the recession constant of the slow-release pool.

Conclusions
In this study we presented a methodology to simulate baseflow processes by adding a slow-reacting linear reservoir to the available quick-reacting reservoir of baseflow generation in the SWAT model.The baseflow-enhanced SWAT model was used to simulate the streamflow process in the snow and glacier melt-dominated Manas River basin in the Tianshan Mountains, where the streamflow process is featured with steeper rising and receding limbs from May to September, and a quite flatter recession from October to April.The following conclusions were achieved.Combination of two linear reservoirs lead to the best results in reproducing the streamflow processes.The Nash-Sutcliff efficiency values at the calibration and validation stages are 0.68 and 0.62 for the one-reservoir case, and 0.76 and 0.69 for the two-reservoir case; the percent bias for both cases is better than −4 %.The two-reservoir approach improved the streamflow flow remarkably, especially the lowflow.
The filter-based approach responds immediately to surface runoff occurrence at the onset of rising limb, while the model-based approaches respond with a delay.In consideration of retention of surface storage and impeding effects of frozen soil layers on infiltration during the initial snowmelt stage, a delay response is believed more reasonable.Meanwhile, it is suggested that freezing/thawing processes be included in soil modules to determine recharge to the aquifers and hence, the proper timing of baseflow response.The filterand model-based methods gave similar surface runoff cessation points which are in late September.
The filter-based method estimated the baseflow index as 0.60, while the model-based as 0.45 in the Manas River basin. it cannot be decided which is more representative due to the unavailability of baseflow measurements.

( 1 )Fig. 2 .
Fig. 2. The Manas River basin and subbasin delineation map (the upper-left small figure is the sketch map of China indicating the location of the study area).
Shihezi Weather Station (SWS, 43 • 29 N, 87 • 06 E) is located below the outlet with an elevation of 444 m a.s.l.The daily meteorological data include maximum and minimum temperatures, wind speed at 10 m height, relative humidity, precipitation, and 20 cm-pan evaporation from 1961 to 1999.This area displays an alpine climate, very cold winter and moderate summer temperatures.The mean high temperature is 39.6 • C, the low −31.7 • C, and the daily average 7.0 • C. The mean annual precipitation is 196 mm and the pan evaporation 1714 mm.

Fig. 5 Fig. 5 .Fig. 6
Fig. 5 Clip illustration of the seasonal baseflow processes generated by model-and filter-based approaches, the whole period covered 1961-1999,

Fig. 6 Fig. 6 .
Fig. 6 Comparison of the simulated and measured monthly averaged low-flow discharge rates.(a), the one-reservoir approach was used; (b), the two-reservoir approach was used.

Table 1 .
Baseflow parameter values for one reservoir and two reservoir approaches in SWAT and the automated baseflow filter program for the Manas River basin, Tianshan, China.

Table 2 .
The NSE and PBIAS for the simulated discharge by SWAT using the one-reservoir and two-reservoir baseflow approaches in the Manas River basin, Tianshan, China.
NSE Rating *PBIAS(%) Rating * * The rating is based on rules given by

Table 3 .
Statistical analysis for the baseflow components and index for Manas River basin, Tianshan, northwest China.