Climate change impact on water resource extremes in a headwater region of the Tarim basin in China

The Tarim river basin in China is a huge inland arid basin, which is expected to be highly vulnerable to climatic changes, given that most water resources originate from the upper mountainous headwater regions. This paper focuses on one of these headwaters: the Kaidu river subbasin. The climate change impact on the surface and ground water resources of that basin and more specifically on the hydrological extremes were studied by using both lumped and spatially distributed hydrological models, after simulation of the IPCC SRES greenhouse gas scenarios till the 2050s. The models include processes of snow and glacier melting. The climate change signals were extracted from the grid-based results of general circulation models (GCMs) and applied on the station-based, observed historical data using a perturbation approach. For precipitation, the time series perturbation involves both a wet-day frequency perturbation and a quantile perturbation to the wet-day rainfall intensities. For temperature and potential evapotranspiration, the climate change signals only involve quantile based changes. The perturbed series were input into the hydrological models and the impacts on the surface and ground water resources studied. The range of impact results (after considering 36 GCM runs) were summarized in high, mean, and low results. It was found that due to increasing precipitation in winter, snow accumulation increases in the upper mountainous areas. Due to temperature rise, snow melting rates increase and the snow melting periods are pushed forward in time. Correspondence to: P. Willems (patrick.willems@bwk.kuleuven.be) Although the qualitive impact results are highly consistent among the different GCM runs considered, the precise quantitative impact results varied significantly depending on the GCM run and the hydrological model.


Introduction
Because the ecological situation of the Tarim River basin in China has seriously declined since the early 1950's and this is caused by over-consumption of irrigation water and unreasonable hydrological conditions (Han and Luo, 2010;Tong et al., 2006), many hydrological studies have been done to properly understand the processes of that watershed, to restore the ecological system and to support the more sustainable development of the basin's water resources (Huang et al., 2010;Jiang et al., 2005;Liu et al., 2006;Zhang et al., 2007).In the meantime, the global climate is changing due to increasing greenhouse gas emissions.Recent research on historical trends of hydro-climatic variables based on historical measurement records have confirmed that temperature at different climate stations of the Tarim River basin has increased consistently with about 1 degree in the past 50 yr (Chen et al., 2007(Chen et al., , 2008;;Li et al., 2003;Shi et al., 2006;Wang et al., 2003).The precipitation tends to have increasing trends since the 1970s and the total water resources are slightly increasing.Some authors reveal that the shrinking of the glacier areas due to climate change in the past decades is the main cause (Jiang and Xia, 2007;Li and Williams, 2008;Li et al., 2007;Xu et al., 2007).The impacts of climate Published by Copernicus Publications on behalf of the European Geosciences Union.
T. Liu et al.: Climate change impact on water resource extremes in China change on the hydrologic systems have already been observed and further changes are expected.Hereby, there is the demand to detect the possible impacts of climate change on water resources at the regional scale to provide important scientific information in support of the basin's water resources management.The future availability of water resources will among other factors, depend on the future greenhouse gas (GHG) emissions.
This study has quantified the impact of scenarios on future changes in GHG emissions on the surface and ground water resources in one of the most important headwater regions of the Tarim river basin: the Kaidu river basin.Use is made of two hydrological models developed for the basin in previous studies (Liu et al., 2011;Liu and Willems, 2011): a lumped conceptual model (VHM) and a more detailed, physically based and spatially distributed model (MIKE-SHE).The strong difference in model structure between these models allowed study of the sensitivity the climate change impact results on the way the hydrological system of the Kaidu basin has been conceptualized.
The climate change signal was obtained for precipitation and temperature for the period 2046-2065 from simulations with general circulation models (GCMs).All GCM runs available for the study region in the AR4 Archive of the IPCC were considered.These are 36 runs obtained based on the IPCC SRES greenhouse gas emission scenarios A1B, A2 and B1 (Nakicenovic, 2000).The observed data for the period 1979-1998 at the single meteorological station site were considered as the baseline data for the climate change impact analysis.An advanced perturbation technique was applied to transfer the climate change signals from the grid-based GCMs to the point-based meteorological station data series.
This paper is organized as follows: Sect. 2 gives a brief description of the study area; Sect. 3 summarizes the hydroclimatic and hydrometric data considered; Sect. 4 gives a description of the two hydrological models considered; Sect. 5 describes the climate change signals obtained from the GCM runs and the downscaling method applied; thereafter the climate change impact results based on the hydrological models are reported and discussed in Sect.6; and final conclusions are obtained in Sect.7.

Study area
The Kaidu river basin is located in the upper region of the Tarim river basin in China in the area ranging from 82 • 58 to 86 • 55 E longitudes and from 41 • 47 to 43 • 21 N latitudes (Fig. 1).It has a complex topography including grassland, marsh and surrounding mountainous alpine areas.The basin area is around 31 750 km 2 .In that area, the elevation of the Kaidu river basin ranges from 5298 m down to 1905 m.
The Kaidu subbasin area contains one meteorological station (at BYBLK) and two flow gauging stations (at BY-BLK and DSK).The BYBLK meteorological and flow gauging stations are located in the upper mountainous region at 2465 m elevation.The DSK flow gauging station is located downstream in the basin.Based on the available observation period (1953 ∼ 2006) at BYBLK, 270 mm of average annual precipitation was observed.The highest and lowest daily average temperatures are respectively 18.8 • C and −42.6 • C. The frozen period lasts (on average) for approximately 180 days yr −1 in the mountainous area and 140 days yr −1 in the grassland area.
The water resources in the Kaidu basin mainly originate from the rainfall in the summer season, the snowfall in the winter season, the snow melting in spring and the perennial glacier melting in summer (Chen et al., 2006;Yang and Cui, 2005).The precipitation in the snow season (winter) is stored in solid form (snow accumulation or frozen soil) and released by increasing temperatures in spring.As shown in Fig. 2, the first flood peaks caused by snow melting appear in April.The latter flood peaks in June-July-August are generated by precipitation.The precipitation between June and August contributed to 45.2 % of the total river flow volume at DSK in the observation period.This total river flow from the Kaidu basin supplied 47 % of the total flow volume in the lower Tarim river (based on the data of 2003) (Wang et al., 2005).

Data overview
The data required to set up the hydrological models and as a baseline for the climate change impact study were provided by the Xinjiang Meteorological Administration (XMB; for the meteorological data); by the Tarim Water Resources Management Bureau (TMB; for the river gauging and other hydrological data); and by the Xinjiang Institute of Ecology and Geography (XIEG) of the Chinese Academy of Sciences (for the geographic data).These station based and geographic data were complemented with remotely sensed data products for the spatially distributed hydrological modelling.These products were derived by the Flemish Institute for Technological Research (VITO) from Moderate Resolution Imaging Spectroradiometer (MODIS) and Fengyun satellite data.
For the lumped conceptual model, the period from 1 July 1980 to 31 December 2002 was selected for model calibration and the period from 1 July 1958 to 30 June 1962 for model validation.For the distributed model, the period from 2 January 2000 to 30 December 2001 was used for model calibration and the period from 2 January 2005 to 30 December 2005 applied for model validation.These periods were selected so that they included several dry and wet periods, following the recommendation of Klemes (1986).
Hydrol.Earth Syst.Sci., 15, 3511-3527, 2011 www.hydrol-earth-syst-sci.net/15/3511/2011/In support of the climate change impact investigations, GCM simulation results were extracted from the database of the World Climate Research Programme's (WCRP) Coupled Model Intercomparison Project phase 3 (CMIP3).Scenario runs from 2046 to 2065 and control runs from 1979 to 1998 were used.Both types of runs have a standard length of 20 yr.The control runs are climate simulation results for the end of the 20th century (20CM3) and were used in this study to evaluate the performance of the GCM simulations.The simulations use the A2, A1B and B1 IPCC SRES GHG emission scenarios.In total, 36 runs with 18 different GCMs were available for the study region.The spatial resolutions of these GCMs range from 1.4 to 5 degrees, which means that only 4 ∼ 64 grid points cover the study area, depending on the GCM's resolution.
Daily precipitation and maximum and minimum surface temperature results were directly obtained from the GCM outputs.The hydrological models, however, also require data on potential evapotranspiration (PET).This PET was calculated from the daily maximum and minimum surface temperature values by the Modified FAO Penman Eq. ( 1).More advanced methods exist for the estimation of PET, but these require data on additional variables like relative humidity, wind speed, number of sunshine hours, etc.In this study, only temperature is used to calculate the PET for two reasons: one is that the additional variables typically have higher uncertainties in the GCM outputs; the second reason is that long-term historical series for the BYBLK meteorological station were only available for temperature.
The R n net radiation at the surface was obtained from the FAO manual (Allen et al., 1998), u 2 was taken equal to the yearly average wind velocity, e s was calculated as the mean between the saturation vapour pressure at both the daily maximum and minimum air temperature e o (T max ), e o (T min ) by Eq. ( 2), e a based on the dew point temperature T dew , which in this study was approximated by the daily minimum air temperature T min (Eq.3), and based on the mean temperature (Eq.4).The constant γ was derived in function of the altitude (Eq.5).Hereby, the entire equation could be simplified with only two variables: maximum and minimum surface temperature.(5) where z [m] is the elevation above sea level.

Lumped conceptual model
The lumped conceptual model has been built (model equations identified) and its parameters calibrated using a databased approach.It starts from a generalized, lumped, conceptual model structure which is detailed in a step-wise way, using information extracted from the available meteorological and flow data records by means of a number of time series techniques (Willems, 2009(Willems, , 2010)).These include the separation of the river flow series in subflows, the splitting of the flow series in nearly independent, quick and slow flow hydrograph periods, and the extraction of nearly independent peak and low flows.The final model structure obtained is shown in Fig. 3.The precipitation contributes to the soil storage, baseflow, interflow and overland runoff subflows.Soil storage is emptied by evaporation and recharged by precipitation, where the recharge speed is determined by the soil moisture status of the previous time step.Subsequently, precipitation contributes to the overland flow and interflow according to the soil moisture status and/or precipitation intensity.The baseflow is the total precipitation residual of soil water storage, overland flow, and interflow.Furthermore, the surface discharge is also recharged by the seasonal snow melting in spring and the permanent glacier melting in summer.The temperature is the only meteorological variable considered to quantify the thermal input which controls the melting.The melted water contributes to the different runoff subflows as rainfall does.More details on the lumped conceptual model and its identification and calibration process are provided by Liu and Willems (2011).The Nash-Sutcliffe model efficiency coefficient (EF; Nash and Sutcliffe, 1970) of the calibrated lumped conceptual model is 0.49 for the daily runoff flows at BYBLK river gauging station and 0.629 at DSK station.For the validation period, the EF is 0.44 at BYBLK station and 0.767 at DSK station.Following the guidelines by Moriasi et al. (2007), these simulation results can be considered acceptable.

Spatially distributed model
The more detailed and more physically based spatially distributed hydrological model was implemented in the MIKE-SHE software of DHI Water and Environment (DHI, 2008a, b).The model includes five model components: overland flow, channel flow, evapotranspiration, unsaturated flow and saturated flow (DHI, 2008a, b).A simulation resolution of 5 km is used.The spatial variability of precipitation was assessed based on a lapse rate of 25 % per 100 m when below 2500 m altitude and 6 % per 100 m when above 2500 m altitude.The temperature correction lapse rate was set as −0.5 • C per 100 m, which was calculated from MODIS remote sensing-based land-surface temperature data.A crop coefficient of 0.2 was used to transform pan evaporation to potential evapotranspiration.The melting temperature on the southern mountainous slopes was set as 0.5 • C and 1.5 • C on the northern slopes.The melting was modelled based on the degree-day temperature method, where the degree-day coef-ficient is set 2.5 mm • C −1 day −1 .The initial total snow storage was derived from MODIS-based snow cover data and the snow storage variations were based on the cumulative precipitation in the early winter period.The leaf area index (LAI) input to the model was calculated based on the MODIS product of the Normalized Difference Vegetation Index (NDVI).Soil physical profiles were described by the Van Genuchten equations for both retention curve and unsaturated hydraulic conductivity.The vertical unsaturated soil profile was discretized in 5-50 cm computational layers.The hydraulic conductivity in the saturated zone was calculated from the soil grade.A computational distance of 5 km and a time interval of 15 s were used.The MIKE-SHE catchment model was bi-directionally linked to a full hydrodynamic model of the river, implemented in the MIKE11 software of DHI Water and Environment, to account for the surface-ground water interactions.For full details on the spatially distributed model and its calibration and validation results, the reader is referred to Liu et al. (2011).The EF of the calibrated MIKE-SHE model is 0.50 for the daily runoff at BYBLK station and 0.72 at DSK station; the validated EF is −0.69 at BYBLK station and 0.27 at DSK station.

Calibration and validation results
Next to the EF coefficient calculated based on the total flows, the goodness-of-fit of the hydrological models was www.hydrol-earth-syst-sci.net/15/3511/2011/ Hydrol.Earth Syst.Sci., 15, 3511-3527, 2011 evaluated based on the following statistics: the correlation coefficient (CF) between simulated and observed flows (Rodgers and Nicewander, 1988) and the cumulative volume error (VF) in the simulated flows.Table 1 lists the calibration and validation results at BYBLK and DSK gauging stations.Next to the evaluation of the total flows (TF), which are the combined result of the different subflows, also results on the baseflow or slow flow (SF), interflow (IF) and quick flow (QF) were tested.Given that this study includes a focus on high and low flow extremes and on climate change impact analysis, extreme wet and dry conditions also have to be considered explicitly in the model evaluation.This was done based on nearly independent peak and low flows extracted from the daily time series, using hydrological independence criteria.Following the method proposed by Willems (2009), two successive events are considered independent if (i) the time length of the decreasing flank of the first event exceeds a given value (to be based on the recession constant of the quick runoff subflow), (ii) the minimum flow in between the two events drops down to a given fraction of the peak flow, and (iii) the peak flow increment (above the lowest flow in between the two events) has a minimum height.
5 Downscaling of the hydro-climatic data

Consistency check
A consistency and accuracy check has been carried out on the GCM simulation results to filter certain GCM runs which are not suitable for the particular study region (Baguis et al., 2010).This accuracy check was based on certain criteria comparing the results for the control runs of the selected GCMs with the observations at the ground meteorological station at BYBLK.The criteria considered in this study area are the Root Mean Square Relative Error (RMSRE; Eq. 6) and the relative bias (Eq.7).
where k is the total number of (daily) time steps, y i the climate model result and x i the observed value at time step i.
The data period ranging from 1979 to 1998 has been considered for the GCM control runs (CN) as well as for the observed data at BYBLK station.The bias and RMSRE statistics have been calculated for the average annual and the monthly volume of precipitation and mean monthly temperature.It seems that the GCMs do not precisely consider the geographic features of the Xinjiang Province, given the consistent positive bias in mean annual temperature.If the GCMs' temperature has consistent bias but high correlation with observations, the error can be considered as a systematic error and does not affect the model accuracy (Wilby and Wigley, 1997).The limited RMSRE values (less than 6 %; see Fig. 4) show that the GCMs have good correlation with the observed annual average temperature.When the monthly variations are analyzed (Fig. 5), higher similarities are found in the relative monthly variations than in the absolute values.The differences in absolute values can be partly due to the differences in hydro-climatic variables at the local point scale of the station considered versus the grid-averaged GCM results.
For precipitation (PR), the bias ranges from −40 % to positive values up to 100 % for the yearly scale or higher for the monthly scale.The RMSRE is higher than 20 %.At the monthly scale, Figs. 4 and 5 show that the climate models are less consistent and reliable in describing precipitation variations than they are in predicting temperature variations.

Climate change and time series perturbation of precipitation
When the changes in hydro-climatic variables are analyzed (from the control period 1979-1998 till the future period 2046-2065), large differences are found between different GCM runs.This uncertainty has been accounted for in this study by making use of an ensemble modelling approach.The entire range of climate change projections from all GCM runs have been clustered into a more limited number of scenarios (mean, high and low) before input into the hydrological models.The mean scenario corresponds to the average value of the change factor for all considered GCM runs, while the high/low scenarios correspond to the maximum/minimum values of the change factors of all considered GCM runs.
The climate change signals were for each scenario extracted from the coarse gridded GCM results and applied to the historical observations at the BYBLK gauging station by using a statistical downscaling method.New time series of daily precipitation for the period 2046-2065 were generated after changing (perturbing) the number of wet days in the historical series (1979)(1980)(1981)(1982)(1983)(1984)(1985)(1986)(1987)(1988)(1989)(1990)(1991)(1992)(1993)(1994)(1995)(1996)(1997)(1998) and after applying wet day perturbation factors to the daily precipitation intensities of the wet days.The latter change in wet day intensity changes was done in a quantile-based way (depending on the return period of the precipitation intensity), following the approach by Ntegeka et al. (2011) and Willems and Vrac (2011).
In the approach, the absolute amounts of precipitation were derived from the GCM scenario runs (SC) and the GCM control runs (CN), and applied to the observation records (OBS) of daily precipitation.This is based on the wet day perturbation factors (for both the wet day frequency and the wet day intensity quantiles) for each of the three (mean, high, low) scenarios.The perturbations were applied in a two-step way.In a first step, the number of wet days was changed.When the wet day frequency increased, wet days Hydrol.Earth Syst.Sci., 15,[3511][3512][3513][3514][3515][3516][3517][3518][3519][3520][3521][3522][3523][3524][3525][3526][3527]2011 www.hydrol-earth-syst-sci.net/15/3511/2011/  were randomly selected from the time series and replaced by zero precipitation values.For months where the wet day frequency decreased, dry days were randomly selected for that month and replaced by a random wet day intensity (randomly selected from the wet day intensities in that month).
In a second step, the wet day intensities were perturbed in a quantile-based way.This means that for each month, all wet day intensities were ranked from the highest to the lowest intensity value and empirical probabilities linked to each intensity value.The same was done for the GCM based results: intensity perturbation factors were calculated as a function of the probability.To each wet day in the input se-ries of the hydrological model (after changes in the wet day frequency), an intensity perturbation factor was applied, depending on the intensity probability of that day.See Ntegeka et al. (2011) and Willems and Vrac (2011) for more details on this method.
The approach is summarized in Fig. 6.For each of the 12 months (January to December), the following steps were followed: -Count the number of wet days (WD; precipitation more than 0.1 mm day −1 ) in SC and CN according to each calendar month; these are denoted WD SC and WD CN; -The perturbation factor (PF) is calculated for each empirical probability based on the ratio of SC over the CN PR intensity for that probability (linear interpolation is applied between the PR intensities when different sets of empirical probabilities are obtained for SC and CN; -The PF is applied to each wet day in OBS (depending on the empirical probability of that day) to generate the perturbed observed series NOBS which includes the climate change signal/trend; -The previous steps are repeated 5 times and one of the 5 generations selected (the one that leads to results closest to the mean monthly PR value of the 5 generations); this step avoids that a random generation is considered containing strong outliers.
Since the wet day perturbation method was originally developed for humid regions where the precipitation events are relatively frequent (Ntegeka et al., 2011), its application for dry regions might be problematic.For the extremely dry conditions of the Tarim basin, the number of dry days might significantly enlarge.Therefore, the results of the wet day perturbation method have been compared in this study with another (simpler) approach, based on monthly volume correction factors.The correction factor was for each month derived from the ratio of average monthly volume of generated precipitation (2046-2065) to average monthly volume of observed precipitation (1979)(1980)(1981)(1982)(1983)(1984)(1985)(1986)(1987)(1988)(1989)(1990)(1991)(1992)(1993)(1994)(1995)(1996)(1997)(1998).It is clear that this change is due to the combined effect of the wet day frequency and quantile perturbations.The two methods are hereafter referred to as the wet day perturbation (WDP) method and the monthly volume (MV) change method.

Wet day frequency and quantile perturbation factors
Figure 7 shows the wet day perturbation factors for the precipitation intensity extremes.The perturbation factor of 1 means that the precipitation in the future period will have no change in comparison with the control period.
In the winter months, precipitation quantile perturbation factors higher than 1 are found for most GCM runs, while higher uncertainty is noticed for the summer months.The precipitation quantile perturbation factors in summer indeed range from positive to negative values, although negative mean values are found for August and September.
Figure 8 shows the mean of the projected future change in the number of wet days per month compared with the control runs.This mean has been calculated based on all GCM runs for each of the three GHG emission scenarios considered.Results indicate that the number of wet days in winter will increase but in summer it might slightly decrease.These changes are similar to the changes found for the precipitation quantiles: the winters become wetter, the summers might become dryer.When the changes are compared for the three GHG emission scenarios considered, it is noticed that the A1B and A2 scenarios have a stronger change between winter and summer than the B1 scenario.

Monthly volume change factors
Figure 9 shows the mean monthly precipitation volume changes (for all GCM runs of the A1B emission scenario).The changes directly obtained from the GCM runs (MV method) have been compared with the changes after application of the WDP method (changes in wet day frequency and quantiles in the input series of the hydrological model).Based on this comparison the validity of the WDP method can be partly tested.
Similar monthly volume changes (by both WDP and MV method) are found in Fig. 9 (increases up to 20 %-30 % in winter months, changes close to zero for summer months at BYBLK station).The monthly volume change factors by the WDP method are slightly higher than by the MV approach.The differences are, however, minor and small in comparison with the maximum and minimum MV changes among all GCM runs considered.Given these minor differences, and given that the WFD method accounts for the dependency of the changes in the precipitation extremes and the precipitation probability, the WDP method was selected for further use in the hydrological impact analysis of climate change (including the impact on hydrological extremes).

Climate change and time series perturbation of temperature and evapotranspiration
Statistical downscaling of temperature and PET has been performed by applying a simpler, more straightforward approach than for precipitation.Given that temperature and PET temporal variations are not organized in events or  wet/dry days, obviously no wet day frequency perturbations were considered.Only mean monthly changes (MV approach) were implemented.The changes are for temperature defined as absolute differences (between future and control periods), while they are defined as relative changes (ratios of the values during future versus control periods) for PET.
Again, high, mean, and low levels of changes were considered, to account for the entire ensemble range of GCM runs considered.Figure 10 shows the mean levels of monthly temperature change for the A2, A1B and B1 emission scenarios.
Hydrol.Earth Syst.Sci., 15, 3511-3527, 2011 www.hydrol-earth-syst-sci.net/15/3511/2011/The temperature change is positive and ranges (for the yearly averaged temperature) from 1 to 3 degrees at BYBLK station.The monthly changes range from 1.5 to 3 degrees at that station.The winter and summer months have higher increases in temperature than the months in spring or autumn.When the emission scenarios are compared, it can be seen that the B1 scenario has the lowest temperature increase.The A1B and A2 scenarios have similar changes.The temperature changes for the A1B and A2 scenarios are about 1 degree higher than for the B1 scenario.
Correspondingly, also the PET (Fig. 11) has mostly positive changes, which are again higher for the A1B and A2 emission scenarios.The highest changes are in the months September and October, with increases of up to about 10 % at BYBLK station.

Climate change impact analysis on hydrological variables
The lumped conceptual model only requires precipitation, temperature and PET as input data (considered at specific locations for model calibration reasons).The distributed model demands (for the same hydro-climatic variables) spatially variable input data.Given that the spatial input variability in the MIKE-SHE model was based on interpolation of the available station data, the same locations were considered for input after climate changes.The historical input series considered on the basis of model calibration and validation were perturbed based on the procedure outlined above (WDP method for precipitation, MV method for temperature and PET).
Next to the combined impact of all hydro-climatic input variables, also the separate effects of each of the perturbed inputs were simulated in order to have a thorough understanding of the contribution of each variable to the total impact.For the lumped conceptual model (due to its model structure and concept), hydrological impacts were analyzed only for the total flow results at the flow gauging stations BYBLK and DSK.For the MIKE-SHE, more detailed investigations could be done on more hydrological variables and including spatial variability.Because the hydrological models need a warming up time (to correct assumed initial catchment conditions), the period of the first three months was excluded from the simulation results.

Combined impact of all hydro-climatic input variables
Figure 12 (1st row) shows the mean level of impact results on the mean monthly flow at the BYBLK and DSK gauging stations.These show the combined effect of the changes in precipitation, temperature and PET.The impact results from single climate variables (precipitation, temperature, ET) aiming to provide insight in the impact results for all climate variables combined.
For BYBLK station, the climate change impact on the total catchment runoff and river flow is positive during the entire year (average changes in monthly volumes in the range 1.3 % ∼ 12.8 %).For the VHM results at DSK station, the flow changes are strongly positive during the snow melting period April-May (17.7 % ∼ 29.7 %).The snow melting period becomes longer at both the BYBLK and DSK stations.The impact results of the MIKE-SHE model are, however, highly different, mainly at DSK station.At BYBLK station, a positive impact (0.4 % ∼ 25.7 %) is found during July-August, while at DSK station the MIKE-SHE discharges reduce to −6.2 % ∼ 3.7 %.During the snow melting period, a negative impact (−9.5 % ∼ −5.4 %) occurs in April and a positive impact (2.3 % ∼ 3.2 %) in May.Discharge increases are found in winter (4 % ∼ 29.2 %).Interestingly, for both hydrological models and both stations, the three emission scenarios (A1B, A2, B1) give similar impact results.
For the VHM model, the total flow in the snow melting period increases due to more snow accumulation in winter.The melting period becomes longer at DSK station in comparison with BYBLK station.After the snow melting period, the flow at BYBLK station increases while it decreases at DSK station.Next to the increase in snow melting volume, there is also a shift in time of the snow melting period due to the temperature rise.This contributes to the flow increase in April-May and decrease in the following months.
For the MIKE-SHE model, higher temporal variability in the impact results are found when compared with the impact results of the lumped conceptual model.Considerable discharge decreases are found in April and July at DSK station, while consistently positive impacts up to about 25 % are found during the entire year at BYBLK station.
These changes by the VHM and MIKE-SHE models are explained based on the impact results for the individual hydro-climatic input variables; see Sect.6.2.
The MIKE-SHE also provides information on the spatial variability of the hydrological outputs.Figure 13 shows the mean level of impact results on the groundwater levels based on the simulation period 2000 ∼ 2001.The groundwater level impacts averaged over the Kaidu river basin are shown.The average groundwater change over the entire study area ranges from −0.08 cm to 0.41 cm. Figure 14 shows that the northern mountainous areas have positive impacts.This is because of the increase in rainfall and related infiltration and percolation.The groundwater recharge increases in the snow melting period (March to May), but decreases in the rainy season (May to September) and continuously feeds the groundwater zone even in winter (September to March).The groundwater table changes are more positive in the higher elevation areas and less positive or even negative in the southern plain region.
While Fig. 12 compared the mean level of impact for the two hydrological models, the two stations and the three emission scenarios, Fig. 15 shows the high, mean, and low levels of impact for the VHM results at BYBLK discharge station.It can be seen in Fig. 15 that the different GCM runs give different levels of impact.The range of impact results covered by the high, mean and low climate change scenarios are indeed very wide, and are mainly due to the differences in climate changes projected by the different GCM runs considered, hence due to uncertainties in the GCMs.

Separate impacts of each hydro-climatic input variable
Similar impact investigations were carried out, but after changes to each input variable, separately.Figure 12 also summarizes these results.The changes in precipitation strongly contribute to the changes in total flows.They increase the monthly flows in the snow melting period (April and May) in the VHM model up to 29.2 % at DSK and up to 15.7 % at BYBLK station due to more snow accumulation in winter.After the snow melting period, the total flow increases in the VHM model at BYBLK station and  The temperature rise increases the discharge that originates from the mountainous area in winter.Due to the increased temperature, the time span that the temperature is over the CT becomes longer; subsequently more precipitation remains in liquid form.The temperature increases also shifts the snow melting period earlier in time by around 1 month.These two phenomena contribute to a flow increase of 7.3 % ∼ 9.6 % at BYBLK station and 1.7 %∼ 2.4 % at DSK station in March after the VHM modelling.The temporal shift of the snow periods causes a decrease in runoff flow during the original snow period April-May (−4.5 %∼ −3.9 % at BYBLK and −6.4 % ∼ −6.1 % at DSK).BYBLK station is more sensitive to the temperature rise due to its location in the upper mountainous region.The temperature increase barely influences the discharge after the snow melting at DSK station since the runoff from the mountainous region only covers a small portion of total flow at this station, which is at the outlet of the Kaidu River basin.In the MIKE-SHE model, the increased temperature also brings the snow melting period forward, which leads to a positive change in March (6.7 % ∼ 8.6 % at BY-BLK and 7.7 % ∼ 8.4 % at DSK) and a negative change in April (−2.7 % ∼ −1 % at BYBLK and −11.2 % ∼ −9.4 % at DSK).The increased temperature also postpones the frozen season and results in more liquid rainfall that transforms into discharge.Both stations respond to temperature in a similar way but with different magnitudes.The mountainous region (BYBLK station) is more sensitive to the temperature change for both models.
The impact results of the changes in PET are also interesting.They do not simply reflect the changes in temperature, but also make clear how the differences in the soil moisture model implementation affect the results.The snow melting submodel in VHM assumes that the soil moisture content keeps constant in the winter season and keeps its level at saturation during the snow melting period.This is because the frozen soil does not change so much by PET.Moreover, the snow evaporation has not been accounted for due to lack of data to estimate/calibrate that process.The PET change leads in the VHM model to small negative flow impacts all over the year.The impact magnitudes are similar at BYBLK and DSK stations and range from −1.1 % to −0.2 % (annual changes).These small changes are due to the soil moisture modelling assumptions in VHM, as explained above.In the MIKE-SHE model, the PET changes also lead to negative flow impacts (−2.1 % ∼ −0.3 % annual change) but with higher temporal variations.The increased PET has stronger effects on the mountainous runoff, while it has stronger effects on the discharges downstream of the catchment during the snow melting period.This is due to the higher water availability in that period.
Only for the VHM results at DSK station, stronger (negative) changes are found.This is due to the stronger influence of the precipitation changes at that station.The control area of BYBLK station (mainly mountainous region) is only 20.4 % of the catchment area at DSK station.The runoff from the mountainous region therefore only has a limited contribution to the total flow at DSK station, while the effect of precipitation is stronger.The reduced precipitation in July and August (44.1 % of annual precipitation) strongly influences the runoff variations at DSK station.
Figure 17 shows how the impact on the peak flow extremes depends on the return period for the VHM results.It is con- The results show that the distributed model leads to stronger changes in mean and peak flows than the lumped conceptual model.For DSK station, strongly different changes are found for the lumped conceptual and distributed models.For the precipitation and combined impacts, the VHM model projects negative mean, peak, and low flow impacts, while the MIKE-SHE model projects clear positive impacts.Section 6.3 discusses the physical reasons for these differences.

Discussion on VHM versus MIKE-SHE results
In previous impact results, the spatially distributed model was found to be more sensitive to scenario precipitation than the lumped conceptual model.The reason is that the lumped conceptual model treats the increased scenario precipitation in the frozen season as snow storage, which does not contribute to runoff after the temperature declines below the critical level.The increased flow with the lumped conceptual model in winter is due to the increased baseflow in summer.For the spatially distributed model, however, even when the temperature at the gauging station declines to the critical level, part of the catchment is not yet cooler than the critical temperature and the precipitation can still contribute to the runoff.This is different for the lumped conceptual model, for which the temperature input is only based on the temperature gauging station located at mid altitude in the catchment.When the gauged temperature reaches the critical level, the whole catchment starts melting at the same time in the lumped conceptual model.The differences in temperature scenario impacts (snow melting differences) are thus mainly due to the difference in spatial detail and related variability across the river basin.Also in the summer season, the results are significantly different.Despite the decline in summer precipitation, the glacier submodel in the spatially distributed model stores more precipitation in winter and consequently provides higher runoff in summer.This is also because the glacier is assumed to have enough ice storage in the spatially distributed model.Both types of models show similar impact results of the scenario PET.Because the PET has a higher influence on baseflow than on overland flow, the range of change for total flow due to change in PET is rather small.It can be seen that the integrated impacts from all climate variables are almost fully explained by the impacts of the changes in precipitation and temperature.This is trivial given that precipitation and temperature are dominant variables in control of the local surface water resources.

Conclusions
Climate change impacts on the surface and ground water resources of one of the main headwater regions of the Tarim basin have been studied.They were based on high, mean and low climate change scenarios, after statistical analysis and downscaling of a large set of GCM runs available for the region.Wet day frequency and quantile perturbations were applied to the inputs of the hydrological models for precipitation, while a simplified procedure (based on changes in mean monthly volumes) was applied for temperature and PET.
Although the climate change scenarios are by far the most important influencing factor, the climate change impact results are also strongly affected by the hydrological model structure.As shown in Figs. 12 and 16, the impact results of the lumped conceptual model and the spatially distributed model are similar in order of magnitude, but the impact temporal distributions are quite different.The different emission scenarios, however, result in a similar impact magnitude.
The spatially distributed model was found to be much more sensitive to scenario precipitation than the lumped conceptual model due to the ignorance of the spatial variability by the latter model.In the lumped model, when the temperature input reaches the critical level, the whole catchment starts/stops melting at the same time while spatial variability in temperature, PET and temperature/PET changes are accounted for in the distributed model.This higher variability leads to more gradual starting/ending of melting/runoff in the MIKE-SHE model.It was seen that the hydrological climate change impacts were mainly controlled by the changes in precipitation and temperature.
We expect from the analysis of the results and their discussion that the impact results simulated by the spatially distributed model are more realistic than by the lumped conceptual model.The distributed model accounts for the spatial variations, which have shown to be important in accurately describing the temporal variations in snow melting and runoff.On the other hand, the lumped conceptual model can more easily simulate the effect of climate change scenarios based on 20 yr input series.Such long-term simulation is needed if one wishes to study the impacts on extreme events.The long computational time of the spatially distributed mod-els does not allow such long continuous time series simulation to be carried out in a reasonable time.Simulations with that model have to be limited to a few years.Given that such a period contains a limited number of extreme events, it is difficult to assess the impact on extreme events with such a model.Extreme events (e.g.floods and water shortage) are, however, important events in support of water management.One solution to overcome the limitation of the spatially distributed model is to carry out model simulations based on multiple statistical categories, e.g. for wet years, dry years, hot years, cool years or normal years, and to attach frequencies to the categories.Another solution would be to look for a compromise solution between the lumped conceptual and detailed spatially distributed models, thus to add some (limited, well-thought out) details to the lumped model.

Fig. 1 .
Fig. 1.Location of the Kaidu river basin in the Xinjiang province, the climate stations and discharge stations (based on the map of Xinjiang Uygur Autonomous Region, certificate No. GS(2007)268).The area simulated by the lumped conceptual models is represented by the blue polygon, while the area simulated by the spatially distributed model is surrounded by the black lines.

Fig. 2 .
Fig. 2. Mean daily hydro-climatic variations throughout the year at BYBLK meteorological station inside the Kaidu basin.

Fig. 5 .--
Fig. 5. Precipitation monthly variations (a), monthly average temperature (b), monthly ETo (c) and return period of annual total precipitation (d) for GCM based results and observations for 1979-1998 at BYBLK meteorological station.

Fig. 6 .
Fig. 6.Steps in the quantile perturbation approach to precipitation data.

Fig. 7 .
Fig. 7. Perturbation factors for wet day precipitation intensity extremes: mean factors for precipitation quantiles higher than the quantile with exceedence probability 0.1 in each month (top) and perturbation factor versus exceedance probability for February (bottom).

Fig. 10 .
Fig. 10.Mean level of perturbation factor of monthly temperature for A1B, A2 and B1 emission scenarios at BYBLK meteorological station.

Fig. 12 .
Fig. 12.Mean level of impact results on mean monthly flows for climate variables (precipitation, 2nd row; temperature, 3rd row; ETo, 4th row; all 3, 1st row) at BYBLK (blue line) and DSK (red line) gauging stations after VHM (left column) and MIKE-SHE (right column) simulations (lines connect the A1B results).The horizontal lines represent the impact results of the A1B scenario, while the small bars represent the max and min impact results of the 3 scenarios (A1B, A2 and B1).

Fig. 13 .
Fig. 13.Mean level of impact results on groundwater (average over the simulation; period 2000 ∼ 2001 and averaged over the entire Kaidu river basin; after the MIKE-SHE simulation (after changes to all or individual hydro-climatic variables).

Fig. 14 .
Fig. 14.Mean level of impact results (for A1B scenario) on the average groundwater level for the period 2000 ∼ 2001 in the entire study area after MIKE-SHE simulation.

Fig. 15 .
Fig. 15.The high, mean, and low level of impact results on the mean monthly flows at BYBLK after VHM simulations.

Fig. 16 .
Fig. 16.Mean level of impact results on mean, peak, and low flows at BYBLK and DSK gauging stations after VHM and MIKE-SHE simulations (after changes to all or individual hydro-climatic variables).

Fig. 17 .
Fig. 17.Changes in empirical peak flow extreme value distribution at BYBLK station after VHM.

Table 1 .
Model performance evaluation statistics for BYBLK and DSK flow gauging stations, based on daily runoff results after model calibration and validation.