Projections for headwater catchments of the Tarim River reveal glacier retreat and decreasing surface water availability but uncertainties are large

In the Tarim River Basin, water resources from the mountain areas play a key role due to the extremely arid climate of the lowlands. This study presents an analysis of future climate change impacts on glaciers and surface water availability for headwater catchments of the Aksu River, the most important tributary to the Tarim River. We applied a glacio-hydrological model that underwent a comprehensive multivariable and multiobjective model calibration and evaluation, based on daily and interannual discharge variations and glacier mass changes. Transient glacier geometry changes are simulated using the Δh-approach. For the ensemble-based projections, we considered three different emission scenarios, nine global climate models (GCMs) and two regional climate models, and different hydrological model parameters derived from the multiobjective calibration. The results show a decline in glacier area of −90% to −32% until 2099 (reference ∼2008) (based on the 5–95 percentile range of the ensemble). Glacier melt is anticipated to further increase or stay at a high level during the first decades of the 21st century, but then declines because of decreased glacier extents. Overall discharge in the Aksu headwaters is expected to be increased in the period 2010–2039 (reference 1971–2000), but decreased in 2070–2099. Seasonally, projections show an increase in discharge in spring and early summer throughout the 21st century. Discharge changes in mid to late summer are more variable, with increases or decreases depending on the considered period and GCM. Uncertainties are largely caused by differences between the different GCMs, with further important contributions from different emission scenarios in the second half of the 21st century. Contributions from the hydrological model parameters to the ensemble uncertainty were generally found to be small.


Updated tables
The amendments affect table 1 and supplementary  table S1. Lines 13-16 in table 1 accidently reported the mean and the 33 and 66 percentiles in the original paper. They have now been amended to correctly report the median and the 5 and 95 percentiles.
Supplementary Table S1. Projected reductions in glacier area for studies in Central Asia.

Updated figures
The amendments affect figures 3-10 and supplementary figures S5-S8.      , and temperature (c), (f), shown as coefficient of variation (CV) for discharge and precipitation, and as standard deviation (std) for temperature. To exclude influences of systematic discharge changes, a 20 year moving average filter was applied before calculating the CV or standard deviation (see text).
Supplementary Figure S5. Projections of the cumulative glacier volume loss for the Sari-Djaz (a) and Kakshaal catchment (b).The black line and the grey area represent the median and the 5-95 percentile range of the ensemble (Sari-Djaz N=162 and Kakshaal: N=756). Coloured lines represent median values for three selected GCMs representing cold (light blue), warm-wet (yellow), and warm-dry (red) projections. The grey dashed and dotted lines show the median of the projections based on the RCMs CCLM and REMO.
Supplementary Figure S6. Projected changes in the discharge regime in percentage for the '2020s', '2050s', and '2080s' as compared to the control period 1971-2000 for the Sari-Djaz (top) and the Kakshaal catchment (bottom). The black line and the grey area represent the median and the 5-95 percentile range of the ensemble (Sari-Djaz N=162 and Kakshaal: N=756). Coloured lines represent median values for three selected GCMs representing cold (light blue), warm-wet (yellow), and warm-dry (red) projections. The grey dashed and dotted lines show the median of the projections based on the RCMs CCLM and REMO.
Supplementary Figure S7. Projection for the evolution of (a) rain and snowmelt, (b) glacier melt, (c) AET, and (d) the water balance (calculated as the sum of rain, snowmelt, and glacier melt minus AET) over the 21st century for the Sari-Djaz catchment (smoothed curves using a 10 year moving average). The black line and the grey area represent the median and the 5-95 percentile range of the ensemble (Sari-Djaz N=162 and Kakshaal: N=756). Coloured lines represent median values for three selected GCMs representing cold (light blue), warm-wet (yellow), and warm-dry (red) projections. The grey dashed and dotted lines show the median of the projections based on the RCMs CCLM and REMO.

Introduction
Mountain areas influenced by snow and glacier melt play a key role in regional water supply (Viviroli et al 2007), yet these regions are particularly sensitive to temperature changes (Barnett et al 2005). Due to its semi-arid to arid lowlands, Central Asia is a prominent example for a region that strongly relies on mountain water resources (Viviroli et al 2007). Irrigation agriculture and hydropower generation are of significant economic importance. The overuse of water resources has become a severe problem since the 1960s, with the consequences of widespread desertification, the Aral Sea shrinkage (Micklin 2007), and drying-out of parts of the lower Tarim River and its previous terminal lakes (Hao et al 2009).
Observations over the last 50 years show significant temperature increase and glacier shrinkage Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. (Sorg et al 2012, Unger-Shayesteh et al 2013. Glacier volume of the Tien Shan is estimated to have decreased by 27±15% over the period 1961-2012 . With continuing climate change, increasing challenges are anticipated for the water management in Central Asia (Siegfried et al 2012), and projections of future water availability are therefore of key significance.
However, so far only a limited number of studies investigated the implications of future climate change on glaciers and water resources for this region. At the regional scale, climate change impacts on glacier extents for the upstream parts of the Amu Darya and Syr Darya have been analysed (Lutz et al 2013), and several studies investigated the response of glaciers and the hydrology to climate change for individual catchments (Hagg et al 2013, Sorg et al 2014, Gan et al 2015, Ma et al 2015. However, for the Tarim River Basin, state-of-the-art studies of the climate change impact on water availability are still missing. Earlier studies on the climate change impact on water resources in this region lacked changes in glacier extents (Liu et al 2010(Liu et al , 2011(Liu et al , 2013, or only focused on glacier runoff . The Tarim River Basin is primarily located in the Xinjiang Province in China and has a population of more than ten million people (Zhou et al 2012). Water from the mountain areas has an immense importance due to extremely low precipitation in the plains.
There are large uncertainties in climate impact analyses due to our limited knowledge and understanding of future external forcings, the climate response to them, and the response of glaciers and hydrology to the projected changes in climate. Uncertainties are further caused by the necessary simplifications in models and deficiencies of the input data. Several studies looked into the role of different uncertainty components in glacierized catchments (Schaefli et al 2007, Stahl et al 2008  Here, we analyse impacts of future climate change on glaciers and water availability for the Aksu, the most important tributary of the Tarim River Basin. The applied glacio-hydrological model includes transient glacier geometry changes, and it was calibrated thoroughly using multiple criteria based on discharge and glacier mass balances. The climate projections are based on three emission scenarios, nine global climate models (GCMs), and additionally two regional climate models (RCMs). The objectives of this study are (1) to investigate the impacts of future climate change on glacier and hydrologic variables, as well as the change mechanisms, and (2) to analyse uncertainties from three sources: the climate projections, the emission scenarios, and the hydrological model parameters.

Study area
The Tarim River forms at the confluence of the rivers Aksu, Hotan, and Yarkand. The Aksu River is the most important tributary with a discharge contribution of ∼80%, (e.g., Duethmann et al 2015). As the downstream part of the Aksu is strongly influenced by water management and river discharge does not further increase , we focused on the two mountainous headwater catchments, the Sari-Djaz (or Kumarik) River and the Kakshaal (or Toshkan) River with the gauges Shaliguilanke and Xiehela, respectively (figure 1). The catchment areas are 12 950 km 2 for the Sari-Djaz, and 18 410 km 2 for the Kakshaal Basin, of which 20% and 4% are glacierized, respectively (reference year ∼2008). The average annual discharge volume of the two catchments sums up to 7.6 km 3 a −1 , equivalent to 240 mm a −1 (average over 1957-2004) (Wang 2006).
Observed climate change signals over the recent decades include temperature increases and changes in precipitation

Climate data and bias correction
For model calibration and evaluation over 1957-2004, precipitation and temperature were interpolated from observed station measurements. Radiation and humidity data were retrieved from the Watch Forcing Data, based on ERA-40 and ERA-Interim (Uppala et al 2005, Dee et al 2011, Weedon et al 2011. The methods used to interpolate the meteorological data are reported by Duethmann et al (2015).
For the scenario analyses, we used daily climate projections from nine CMIP5 GCMs forced by three emission scenarios (RCP2.6, RCP4.5, and RCP8.5). The GCMs were selected to represent the range of changes projected by the CMIP5 ensemble in the study area (supplementary figure S1). With a spatial resolution of about 1.5°-3°, the GCMs cannot capture the small-scale topographic features. A better spatial resolution can be achieved through downscaling by RCMs. Data were available from two RCMs: the RCM REMO forced by the CMIP3 GCM ECHAM5 under the A1B emission scenario with a resolution of 0.166° (Mannig et al 2013), and the RCM CCLM forced by the CMIP5 GCM MPI-ESM-LR under the emission scenarios RCP2.6, RCP4.5, and RCP8.5 with a resolution of 0.44° (Wang et al 2013). The RCM based scenarios were used for impact assessment but to keep a consistent ensemble they were not included for uncertainty evaluation.
All climate model data were summarised to subcatchments and bias corrected using an empirical quantile mapping approach (Gudmundsson et al 2012), which allowed us to also consider changes in the variability. The bias correction was applied separately for each month.

Hydrological modelling
Discharge and glacier geometry changes were simulated with the hydrological model WASA (Güntner and Bronstert 2004), which has also previously been applied in Central Asia . The implementation, calibration, and evaluation of the WASA model to the Aksu headwater catchments is presented in detail in Duethmann et al (2015) and summarised in the following.
The model is applied at a daily time step, and the spatial discretization is based on hydrologic response units Snow and glacier melt is simulated using a temperature-index approach with seasonally varying melt factors. For each glacier, area and thickness is updated annually based on the simulated glacier mass balance using the Δh-approach (Huss et al 2010) (for details and model equations, refer to the supplement of Duethmann et al 2015). For model calibration, glacier area changes were based on the glacier inventories for ∼1975 and ∼2008 assuming linear decrease rates for each glacier.
A comprehensive multi-variable approach was used for model calibration and evaluation. Daily discharge variations were calibrated using an average of the Nash-Sutcliffe efficiency of linear and logarithmic discharge values to achieve a balanced evaluation of high and low discharges (equation (1) 1). The criterion was represented by a function that reaches an optimal value of one for the average mass balance estimate, a value of zero outside the uncertainty range, and linearly increasing/decreasing values in between. The temporal variation of the simulated glacier mass balances was calibrated by maximising the correlation to an in situ glacier mass balance series of Karabatkak Glacier (Dyurgerov andMeier 2005, WGMS 2012), located close to the study area (figure 1).
The model was automatically calibrated to the two discharge and the two glacier mass balance criteria with a multiobjective calibration algorithm (ε-NSGAII; Kollat and Reed 2006), using 24 years (1976-1999 for model calibration, and 24 years (1957-1975 and 2000-2004) for model validation. From the results of the multiobjective optimisation, solutions were selected if they showed a good performance for daily and interannual streamflow variations were inside the uncertainty range of the geodetic glacier mass balance estimate, and represented observed discharge trends in the past. Since the calibration period was regarded as too short for the calculation of trends, we calculated trends in seasonal discharge time series over the entire historical period 1957-2004. Trend evaluation was constricted to trend direction and trend significance. As a result, 6 and 28 parameter sets were selected for the Sari-Djaz and Kakshaal catchment, respectively.
Nash-Sutcliffe efficiencies for daily (monthly) discharge are within 0.77-0.87 (0.81-0.94) in the calibration and 0.67-0.84 (0.80-0.92) in the validation period (ranges over both catchments and the selected parameter sets). The model shows an acceptable performance with respect to the interannual variations of seasonal discharge and represents the observed discharge trends. Trends over 1957-2004 are 19-21 mm decade −1 (simulated) and 21 mm decade −1 (observed) for the Sari-Djaz catchment (mean annual runoff 382 mm a −1 ), and 10-11 mm decade −1 (simulated) and 10 mm decade −1 (observed) for the Kakshaal catchment (mean annual runoff 151 mm a −1 ) (trends evaluated with Sen's slope estimator Sen 1968). In agreement with the observed geodetic glacier mass balance estimates, the simulated glacier mass balance estimates over 1976-1999 are in a range of −0.83 to −0.51 m w.e. a −1 (−0.85 to −0.25 m w.e. a −1 ) for the Sari-Djaz (Kakshaal) catchment.

Uncertainty estimation
In this study, we analysed uncertainties from three sources, including three emission scenarios, nine GCMs, and six (28) parameter sets of the hydrological model for the Sari-Djaz (Kakshaal) catchment. This led to 162 ensemble members for the Sari-Djaz catchment and 756 for the Kakshaal catchment. Uncertainties are characterised by the 5-95 percentile range of the ensemble. For consistency, only GCMbased simulations were included in the ensemble. Contributions of individual uncertainty sources were analysed by their contribution to the total variance (see supplement).

Projected changes in climate variables
The analysis of changes in climate and glacio-hydrologic response was carried out by comparing the near future 2010-2039 ('2020s'), mid future 2040-2069 ('2050s'), and far future period 2070-2099 ('2080s') to the control period 1971-2000. Increased air temperatures are projected by all climate models under all emission scenarios (figure 2). Differences in the projected temperature increase between the different RCPs are negligible for the '2020s' but substantial for the '2080s'. Temperature increases at a similar rate over the whole 21st century under RCP8.5, while for RCP4.5, the rate of temperature increase is projected to slow down toward the end of the 21st century, and in RCP2.6, temperature is projected to stabilise or decrease about the end of the century. The resulting temperature increase in the '2080s' compared to 1971-2000 is 3.2°C-7.2°C for RCP8.5, 1.6°C-4.5°C for RCP4.5 and 0.7°C-2.8°C for RCP2.6. Seasonally, the strongest temperature increases are expected for winter (supplementary figure S2).
Precipitation is projected to be about 10% higher for the '2020s' as compared to 1971-2000 for all RCPs, but the variability between GCMs is high. Over the 21st century, most GCMs project further small increases in precipitation. In contrast to temperature, differences between the RCPs remain small.
Compared to the median of the GCMs, the two RCMs show a tendency toward lower precipitation and temperature increases, but are mostly within the interquartile range of the GCMs. From the ensemble, we depicted three GCMs representing warm-wet (MIROC-ESM), warm-dry (IPSL-CM5A-LR), and cold (GFDL-ESM2M) conditions to better comprehend the change mechanisms of the glacier and hydrological variables described in the next sections.
Changes in humidity and radiation are mostly small and within ±5%. Larger decreases in humidity are linked to the warm-dry climate projections (supplementary figure S4).

Changes in glacier extent
The projected changes in temperature and precipitation are simulated to cause predominantly negative mass balances and significant losses in glacier extent (figure 3) and volume (supplementary figure S5). By the end of the century, projections suggest a reduction in glacier area by −66(−89 to −28)% for the Sari-Djaz and −78(−94 to −47)% for the Kakshaal catchment (uncertainties refer to the 5-95 percentile range). The variability of the GCM ensemble contributes most to the overall uncertainty (figure 4). As expected, the strongest glacier area reduction is projected for the warm-dry GCM (figure 3). The influence of the emission scenario increases from the 2050s, and uncertainties from the GCMs and emission scenarios are similar by the end of the century. The pattern at the beginning of the scenario period in figure 4 relates to very small total uncertainties (figure 3) and is therefore not analysed.
The glacier area reduction by the CCLM simulations is slightly lower than the ensemble median, in accordance with the lower than average temperature increase ( figure 3). The REMO simulations suggest a strong glacier reduction from the 2050s onwards, which is in line with the continuous temperature increase under the A1B scenario.

Changes in the runoff regime
In the '2020s', discharge is generally projected to be higher than during 1971-2000 in all seasons (figure 5). In both catchments and over the whole 21st century, discharge in spring and early summer (April-June) is higher than during the control period for nearly all ensemble members. The largest changes in relative terms are observed during spring, where discharge during the reference period is low (supplementary figure S6). Discharge changes in late summer are more heterogeneous. In the Sari-Djaz catchment, August discharge is increased in the '2020s' but decreased in the '2080s'. For the Kakshaal catchment, both increases and decreases are observed during all scenario periods, and the direction of change strongly depends on the choice of GCM, which contributes most to overall uncertainties (figure 6). By the end of the century, emission scenarios are an important uncertainty source during spring and early summer, while the contribution of the hydrological model parameters generally remains small.

Mechanism of changes in the runoff regime
Snowmelt is shifted toward earlier in the season leading to higher snowmelt in spring and lower snowmelt in summer as shown for the '2080s' (figure 7). Increase in spring rainfall further contributes to higher discharge during spring. During the '2080s', lower summer discharge is caused by the shift in snowmelt, decreases in ice melt, and increases in actual evapotranspiration (AET). Summer rainfall shows a large variability in the GCM projections. A decrease in summer rainfall as in the warm-dry GCM may further exacerbate the discharge decreases, while rainfall increases as projected under the warm-wet GCM partly compensate the changes in snowmelt, ice melt, and AET.

Peak glacier melt and evolution of annual runoff
Glacier melt (defined as ice melt from glaciers) is simulated to change substantially over the 21st century in the Sari-Djaz catchment ( figure 8(a)). After increases at the beginning of the 21st century, glacier melt drastically decreases during 2040-2060. Peak glacier melt is projected to occur roughly around 2030 ( figure 8(a)). In the Kakshaal catchment, glacier melt is generally much lower and peak melt is passed earlier ( figure 8(b)).  The changes in glacier melt have a significant impact on the evolution of total annual runoff of the Sari-Djaz River, resulting in the reduction of discharge by −22(−51 to +3)% in the '2080s' compared to the control period ( figure 9(a), table 1, supplementary  figure S7). In contrast, the impact of changes in glacier melt is smaller for the Kakshaal River due to its lower glacier melt contribution to discharge ( figure 9(b),  table 1, supplementary figure S8).
The projected evolution of glacier melt for the RCM-based simulations is similar to those of the GCM-based simulations, while the projected discharge toward the end of the century is slightly lower than the median of the GCM-based  simulations. This is in line with a lower precipitation increase than projected by the median of the GCMs.

Changes in the interannual variability of discharge
Changes in the interannual variability of discharge were evaluated using the coefficient of variation (CV). The CV was calculated after applying a 20 year moving average filter to exclude influences of systematic discharge changes. In both catchments, a tendency toward increased interannual variability of discharge can be seen; however, considering the large variability for different ensemble members, the changes are only small (figures 10(a) and (d)). A moderate glacier cover of a catchment is linked with low interannual discharge variability: in years with low precipitation, glacier melt begins earlier and can therefore (partly) compensate lower snowmelt (Jansson et al 2003, Hock et al 2005. The increase of interannual variability may therefore be caused by decreasing glacier area and glacier melt. In addition, changes in the interannual variability of climate variables might also play a role; there is a slight increase in the interannual variability of temperature (figures 10(c) and (f)), but changes in interannual variability of precipitation are small (figures 10(b) and (e)).   table S1). The projected percentage glacier area loss is larger in the Kakshaal than in the Sari-Djaz catchment, where glaciers are generally larger and have a higher thickness.
Changes in the runoff regime with lower runoff in August and higher runoff in May/June are typical for mountain catchments influenced by snow and glacier melt and have been projected for many regions including other catchments in Central Asia (Hagg et al 2007, Sorg et al 2014, Gan et al 2015, Ma et al 2015. While this study projects decreases in summer discharge in the '2080s' (figure 5), an earlier study for headwater catchments of the Aksu, Hotan, and Yarkand projected only little changes (±5%) for summer discharge in the period 2081-2100 compared to 1981-2000(Liu et al 2013. These differences are likely due to the static (and thus overestimated) glacier areas in the model of Liu et al (2013), which lead to an overestimation of glacier melt.
Estimates of the timing of peak melt vary by catchment and climate projection. Averaged over entire Central Asia, Bliss et al (2014) found a steady decline in glacier runoff over the 21st century, indicating that many glaciers have already passed peak melt. For the Chon-Kemin catchment in the northern Tien Shan, peak melt is projected to occur in the 2020s under warm climate scenarios but glacier runoff may also remain constant until 2100 under cooler climate scenarios (Sorg et al 2014). Only small decreases in glacier runoff until 2050 were projected for the Tanimas Basin in the Pamir (Hagg et al 2013).
For estimating water availability in the mainstream of the Tarim River, one needs to consider changes in land and water management as well as  changes in other tributaries, which provide only a small contribution to the Tarim today but possibly a higher contribution in the future. Due to lacking meteorological and glacier mass balance data in the Hotan and Yarkand Basins, glacio-hydrologic simulations and climate change impact analyses in these basins are, however, highly uncertain.

Uncertainties
In agreement with many other impact studies, our evaluation of uncertainties showed that the largest uncertainties stem from the climate models, while emission scenarios become a further important uncertainty source toward the end of the 21st century (Wilby and Harris 2006, Kay et al 2009, Prudhomme and Davies 2009, Ott et al 2013, Addor et al 2014. However, in data scarce mountain regions, uncertainties of hydrological models are large, and parameter uncertainties of hydrological models can become the largest source of uncertainties (Finger et al 2012, Ragettli et al 2013. The availability of glacier mass balance data in our study likely contributed to relatively small uncertainties from the hydrological parameters. In our case, parameter uncertainties were derived from different solutions of the multiobjective calibration. One would expect more diverse parameter sets if they had been derived by independent optimisation runs or Monte Carlo simulations and further parameter sets with a lower model performance would have been accepted. However, the simulated glacier mass changes (−0.83 to −0.51 m w.e. a −1 and −0.85 to −0.25 m w.e. a −1 ) cover large parts of the range given by the estimate from the observations (−0.85 to −0.19 m w.e. a −1 ) and thus indicate that the parameterizations cover a good range of uncertainties. It is to note that uncertainties of the hydrological model are larger than the analysed parameter uncertainties due to uncertainties in model input and structure, as discussed below.
Further uncertainties in the climate change signal result from the coarse resolution of the climate models, which only partly resolve the complex topography of this region, and from the assumption of the stability of the bias correction under a changed climate. With respect to the glacio-hydrological modelling, further uncertainties stem from uncertainties in the calibration and input data, as well as the model structure. For example, while for the Sari-Djaz Basin, the glacier mass balance estimate used for model calibration largely overlaps with a more recent geodetic estimate derived specifically for this catchment , an estimate specifically for the Kakshaal catchment is still missing. The applied ice-thickness estimates were based on a distributed ice-thickness model. While this type of thickness estimation method is seen as more reliable than approaches based on volume-area scaling, which tend to overestimate ice volume for large and steep glaciers (Frey et al 2014), uncertainties are still estimated to be in the range of 30% and contribute to uncertainties in the glacier and runoff evolution (Gabbi et al 2012. Uncertainties in the model structure result from simplifications of the described processes and the assumption of stable model parameters under a changed climate. For example, the application of a temperature-index approach implicates that influences on the energy balance by changes in radiation or the wind field cannot be considered. However, comprehensive model testing showed that the model represents observed discharge trends and glacier mass loss in the Figure 10. Interannual variability of discharge (a, d), precipitation (b, c), and temperature (c, f), shown as coefficient of variation (CV) for discharge and precipitation, and as standard deviation (std) for temperature. To exclude influences of systematic discharge changes, a 20 year moving average filter was applied before calculating the CV or standard deviation (see text). past, which increases our confidence in the model. While the average influence of debris cover on glacier melt is implicitly taken into account, influences on glacier retreat by possible changes of the debris-cover extent cannot be considered by the current model. Our simulations did also not consider feedback effects of the vegetation, such as changes in the length of the growing period, the vegetation types and the water use efficiency. These aspects should be addressed by future research.

Conclusions
In this study, we analysed climate change impacts on water resources for the most important headwater catchments of the Tarim River. The results revealed a strong decline of the glacier area over the next decades. However, uncertainty ranges are large. Glacier melt is projected to be high in the beginning of the 21st century and to decrease rapidly during 2040-2060. As a result of reduced glacier melt and increased AET, which are projected to be only partly compensated by precipitation increases, overall runoff in the '2080s' is simulated to be lower than in the control period. Seasonally, runoff during spring and early summer is anticipated to be higher throughout the whole 21st century, while summer runoff is projected to be at a high level during the '2020s' but expected to decline afterwards.
In summary, it can be expected that the currently increasing discharge trends are likely not sustained in the long term, which needs to be accounted for by land use planning and water resources management. Increases in irrigated area and water use as observed in the recent decades , Feike et al 2015 already impair ecological conditions and livelihoods today, despite currently increasing water availability in the headwaters. The region will face enormous challenges if the water availability from the Aksu headwaters declines.