Investigating water budget dynamics in 18 river basins across the Tibetan Plateau through multiple datasets

The dynamics of basin-scale water budgets over the Tibetan Plateau (TP) are not well understood nowadays due to the lack of in situ hydro-climatic observations. In this study, we investigate the seasonal cycles and trends of water budget components (e.g. precipitation P , evapotranspiration ET and runoff Q) in 18 TP river basins during the period 1982–2011 through the use of multi-source datasets (e.g. in situ observations, satellite retrievals, reanalysis outputs and land surface model simulations). A water balancebased two-step procedure, which considers the changes in basin-scale water storage on the annual scale, is also adopted to calculate actual ET. The results indicated that precipitation (mainly snowfall from mid-autumn to next spring), which are mainly concentrated during June–October (varied among different monsoons-impacted basins), was the major contributor to the runoff in TP basins. The P , ET and Q were found to marginally increase in most TP basins during the past 30 years except for the upper Yellow River basin and some sub-basins of Yalong River, which were mainly affected by the weakening east Asian monsoon. Moreover, the aridity index (PET/P ) and runoff coefficient (Q/P ) decreased slightly in most basins, which were in agreement with the warming and moistening climate in the Tibetan Plateau. The results obtained demonstrated the usefulness of integrating multi-source datasets to hydrological applications in the data-sparse regions. More generally, such an approach might offer helpful insights into understanding the water and energy budgets and sustainability of water resource management practices of data-sparse regions in a changing environment.


Introduction
As the highest plateau in the globe (the average elevation is higher than 4000 m above the sea level), the Tibetan Plateau (TP, also called "the roof of the world" or "the third Pole") is regarded as one of the most vulnerable regions under a warming climate and is exposed to strong interactions among the atmosphere, hydrosphere, biosphere and cryosphere in the Published by Copernicus Publications on behalf of the European Geosciences Union.
earth system (Duan and Wu, 2006;Yao et al., 2012;Liu et al., 2016b).It also serves as the "Asian water tower" from which some major Asian rivers such as the Yellow River, Yangtze River, Brahmaputra River, Mekong River, Indus River, etc., originate, which is a vital water resource to support the livelihoods of hundreds of millions of people in China and the neighbouring Asian countries (Immerzeel et al., 2010;Zhang et al., 2013).Hence, sound knowledge of water budget and hydrological regimes in TP river basins and their responses to the changing environment would have practical relevance for achieving sustainable water resource management and environmental protection in this part of the world (Yang et al., 2014;Chen et al., 2015).
Despite the importance of the TP in this geographic region, advances in hydrological and land surface studies in this region have been limited by data scarcity (Zhang et al., 2007;Li et al., 2013;Liu et al., 2017).For instance, less than 80 observation stations (∼ 10 % of a total of ∼ 750 observation station across China) have been established in the TP by the Chinese Meteorological Administration (CMA) since the mid-20th century (Wang and Zeng, 2012).These stations are generally sparse and unevenly distributed at relatively lowelevation regions (most stations are located in the eastern TP and few of them situated in the western parts) and focus only on the meteorological variables and lack of other land surface observations such as evapotranspiration, snow water equivalent and latent heat fluxes.In addition, long-term observations of river discharge, lake depth and glacier melts in the TP are also absent (Akhta et al., 2009;Ma et al., 2016).Therefore, the water budget and hydrological regimes for each river basin in the TP and their relation with atmospheric circulations are poorly understood (Cuo et al., 2014;Xu, 2011).Whilst this shortcoming could be resolved through installation of in situ monitoring systems (Yang et al., 2013;Zhou et al., 2013;Ma et al., 2015), the overall cost, labour and technical support for running the operational sites would be substantial.Another workaround would be through modelling approach, i.e. feeding remote sensing information and meteorological forcing data into physically based land surface models (LSMs) to simulate the basin-wide water budget (Bookhagen and Burbank, 2010;Xue et al., 2013;Zhang et al., 2013;Cuo et al., 2015;Zhou et al., 2015;Wang et al., 2016).However, such an approach is not immune from the issue of data scarcity at multiple river basins (with varied sizes and/or terrain complexities) for supporting model calibration and validation purposes (Li et al., 2014).
Most recently, several global (or regional) datasets relevant to the calculation of water budget have been released.They include remote-sensing-based retrievals (Tapley et al., 2004;Zhang et al., 2010;Long et al., 2014;Zhang et al., 2016), LSM simulations (Rui, 2011), reanalysis outputs (Berrisford et al., 2011;Kobayashi et al., 2015) and gridded forcing data interpolated from the in situ observations (Harris et al., 2014).For example, there are many products related to terrestrial evapotranspiration (ET), such as GLEAM_E (Global Land surface Evaporation: the Amsterdam Methodology, Miralles et al., 2011), MTE_E (a product integrated the pointwise ET observation at FLUXNET sites with geospatial information extracted from surface meteorological observations and remote sensing in a machine-leaning algorithm, Jung et al., 2010), LSM-simulated ETs from Global Land Data Assimilation System version 2 (GLDAS-2) with different land surface schemes (Rodell et al., 2004), ETs from Japanese 55-year reanalysis (JRA55_E), the ERA-Interim global atmospheric reanalysis dataset (ERAI_E), and the National Aeronautic and Space Administration (NASA) Modern Era Retrosphective-analysis for Research and Application (MERRA) reanalysis data (Lucchesi, 2012).Moreover, there are also several global or regional LSM-based runoff simulations from GLDAS and the variable infiltration capacity (VIC) model (Zhang et al., 2014).A few attempts have been made to validate multiple datasets for certain water budget components and to explore their possible hydrological implications.For example, Li et al. (2014) and Liu et al. (2016a) evaluated multiple ET estimates against the water balance method on annual and monthly timescales.Bai et al. (2016) assessed streamflow simulations of GLDAS LSMs in five major rivers over the TP based on the discharge observations.Although uncertainties might exist among different datasets with various spatial and temporal resolutions which are calculated using different algorithms (Xia et al., 2012), they offer an opportunity to examine the general basin-wide water budgets and their uncertainties in gauge-sparse regions such as the TP considered in this study.
From the multiple-dataset perspective, this study aims to investigate the water budget in 18 TP river basins distributed across the Tibetan Plateau and evaluate seasonal cycles and annual trends of these water budget components.This paper is organized as follows: the datasets and methods applied in this study are described in Sect. 2. The results of seasonal cycles and annual trends of water budget components for the river basins are presented and discussed in Sect.3. The uncertainties arising from employing multiple datasets are also discussed in the same section.In Sect.4, we generalize our findings, which could be helpful for understanding the water balances of the river basins under constant influence of interplay between westerlies and monsoons (e.g.Indian monsoon, east Asian monsoon) in the Tibetan Plateau.(1952-2012, VIC_IGSNRR simulated) with a spatial resolution of 0.25 • and a daily temporal resolution from the Geographic Sciences and Natural Resources Research (IGSNRR), Chinese Academy of Sciences, is also used.This dataset is derived from the VIC model forced by the gridded daily observed meteorological forcing (IGSNRR_forcing; Zhang et al., 2014).A degree-day scheme was used in the model to account for the influences of snow and glaciers on hydrological processes.
In terms of precipitation (P ), we used the gridded monthly precipitation dataset available at CMA (spatial resolution of 0.5 • ; 1961-2011; interpolated from observations of 2372 national meteorological stations using the thin plate spline method; Table 1).Since the reliability of this dataset might be restricted by the relatively sparse stations and complex terrain conditions of the TP, we make an attempt to incorporate two other precipitation datasets (IGSNRR_forcing and Tropical Rainfall Measuring Mission TRMM 3B43 V7).The precipitation from IGSNRR forcing datasets (0.25 • ) was derived by interpolating gauged daily precipitation from 756 CMA stations based on the Synergraphic mapping system algorithm (Shepard, 1984;Zhang et al., 2014) and was further bias-corrected using the CMA gridded precipitation.
To get the change in terrestrial storage ( S), we used three of the latest global terrestrial water storage anomaly and water storage change datasets (available on the GRACE Tellus website: http://grace.jpl.nasa.gov/)that were retrieved from the Gravity Recovery and Climate Experiment (GRACE, Tapley et al., 2004;Landerer and Swenson, 2012;Long et al., 2014).Briefly, they were processed separately at the Jet Propulsion Laboratory (JPL), the GeoForschungsZentrum (GFZ) and the Center for Space Research at the University of Texas (CSR).To minimize the errors and uncertainty of extracted S, we averaged these GRACE retrievals (2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013) from different processing centres in this study.

Vegetation and snow-glacier parameters
To quantify the dynamics of vegetation of each river basin, we applied the normalized difference vegetation index (NDVI) and the leaf area index (LAI; Table 1).Briefly, the NDVI data were obtained from the Global Inventory Modeling and Mapping Studies (GIMMS;Turker et al., 2005; https://nex.nasa.gov/nex/projects/1349/wiki/general_data_description_and_access/) while the LAI data were collected from the Global Land Surface Satellite (GLASS) products (http://www.glcf.umd.edu/data/lai/;Liang and Xiao, 2012).Whist the change in seasonal snow cover and glacier has a significant impact on the water and energy budgets in TP river basins, it remains a technical challenge to get reliable observations due to the harsh environment (especially on the basin scale).However, recently available satellitebased or LSM-simulated products might provide adequate characterization of the variation of snow and glacier cover.To quantify the change in snow cover at each basin, we applied the daily cloud-free snow composite product from MODIS Terra-Aqua and the Interactive Multisensor Snow and Ice Mapping System for the Tibetan Plateau (Zhang et al., 2012;Yu et al., 2015), in conjunction with the snow water equivalent (SWE) retrieved from the Global Snow Monitoring for Climate Research product (GlobSnow-2, http://www.globsnow.info/)and the VIC_IGSNRR simulations (Takala et al., 2011;Zhang et al., 2014).We extracted a general distribution of glaciers in the TP from the second glacier inventory dataset of China (Guo et al., 2014).All gridded datasets used were first uniformly interpolated to a spatial resolution of 0.5 • based on the bilinear interpolation to make their intercomparison possible.The datasets were then extracted for each of TP basins.

Monsoon indices
In general, the TP climate is under the influences of the westerlies, Indian summer monsoon and east Asian summer monsoon (Yao et al., 2012).To investigate the changes in monsoon systems and their potential impacts on water budgets in the TP basins, we used three monsoon indices, namely Asian zonal circulation index (AZCI), Indian Ocean dipole mode index (IODMI) and east Asian summer monsoon index (EASMI).Briefly, the IODMI (reflects the dynamics of Indian summer monsoon) is an indicator of the east-west temperature gradient across the tropical Indian Ocean (Saji et al., 1999), which can be downloaded from the following website: http://www.jamstec.go.jp/frsgc/research/d1/iod/ dmi.html.The EASMI and AZCI (60-150 • E) reflect the dynamics of the east Asian summer monsoon (Li and Zeng, 2002) and the westerlies (represented by the Asian zonal circulation index), which can be obtained from Beijing Normal University (http://ljp.gcess.cn/dct/page/65577)and the National Climate Center of China (http://cmdp.ncc-cma.net/en/), respectively.

Study basins
In this study, we selected 18 river basins of varied sizes (range: 2832-191 235 km 2 ; see Table 2 for details) with adequate runoff data over a 30-year period .They are distributed in the northwestern, southeastern and eastern parts of the plateau with multiyear-mean and basin-averaged temperature and precipitation ranging from −5.68 to 0.97 • C and 128 to 717 mm, which are solely dominated by or under the combined influences of the westerlies, the Indian Summer monsoon and the east Asian monsoon (Yao et al., 2012).There is more glacier and snow cover in the westerliesdominant basins such as Yerqiang, respectively) and less for the east-Asian-monsoon-dominated basins such as Yellow, Yangtze and Bayin (0-0.96 and 9.42-20.05%, respectively; Table 2).

Water balance-based ET estimation
The basin-wide water balance on the monthly and annual timescales could be written as the principle of mass conservation (also known as the continuity equation, Oliveira et al., 2014) of basin-wide precipitation (P , mm), evapotranspiration (ET wb , mm), runoff (Q, mm) and terrestrial water storage change ( S, mm): (1) The terrestrial water storage ( S) in Eq. ( 1) includes the surface, subsurface and groundwater changes.It has been demonstrated that S cannot be neglected in water balance calculation over monthly and annual timescales due to snow cover change and anthropogenic interferences (e.g.reservoir operation, agricultural water withdrawal; Liu et al., 2016a).
For the period 2002-2011, we calculated basin-wide ET (ET wb ) directly using the GRACE-derived S in Eq. ( 1).Since GRACE data are absent before 2002, we calculated the monthly ET wb using the following two-step bias-correction procedure (Li et al., 2014).We defined P − Q in Eq. ( 1) as biased ET (ET biased , available from 1982 to 2011) relative to the "true" ET (ET wb = P − Q − S, available during the period 2002-2011 when the GRACE data are available).
Over the period 2002-2011, we first fitted the ET biased and Table 2. Main features of the 18 TP river basins used in this study.The precipitation and temperature statistics for each basin were calculated from the observed CMA datasets while the NDVI and LAI statistics were extracted from the GIMMS NDVI dataset and GLASS LAI product.
The GA and SC percentages represent the percentages of multiyear-mean glacier cover and snow cover in each basin, which were calculated from the second glacier inventory dataset of China and the daily TP snow composite products (2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013).
No  ET wb series separately using different gamma distributions, which has been evidenced as a proper method for modelling the probability distribution of ET (Bouraoui et al., 1999).
The monthly ET biased series (2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011) can then be biascorrected through the inverse function (F −1 ) of the gamma cumulative distribution function (CDF, F) of ET wb by matching the cumulative probabilities between two CDFs as fol-lows (Liu et al., 2016a): where α biased , β biased and α wb ,   to that of ET corrected (2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011).The second step of bias correction is to eliminate the annual bias through the ratio of annual ET biased to annual ET corrected calculated in the first step using the following method: where ET final (m) is the final monthly ET after bias correction.ET biased (a) and ET corrected (a) represent the annual biased and corrected ET while ET corrected (m) is the monthly corrected ET obtained from the first step.The procedure was then applied to correct the monthly ET biased series and calculated the monthly ET corrected during the period 1982-2001 for all TP basins.We take these results as sufficient representation of the "true" ET (ET wb ) for evaluating multiple ET products and trend analysis.

Modified Mann-Kendall test method
The Mann-Kendall (MK) test is a rank-based non-parametric approach which is less sensitive to outliers relative to other parametric statistics, but it is sometimes influenced by the serial correlation of time series.Pre-whitening is often used to eliminate the influence of lag-1 autocorrelation before the use of MK test.For example, X(X 1 , X 2 , . .., X n ) is time series data that will be replaced by (X 2 −cX 1 X 3 −cX 2 , . .., X n+1 − cX n ) in pre-whitening if the lag-1 autocorrelation coefficient (c) is larger than 0.1 (von Storch, 1995).However, significant lag-1 autocorrelation may still be detected after pre-whitening because only the lag-1 autocorrelation is considered in pre-whitening (Zhang et al., 2013).Moreover, it sometimes underestimates the trend for a given time series (Yue et al., 2002).Hamed and Rao (1998) proposed a modified version of the MK test (MMK) to consider the lag-1 autocorrelation and related robustness of the autocorrelation through the use of equivalent sample size, which has been widely used in previous studies during the last 5 decades (McVicar et al., 2012;Zhang et al., 2013;Liu and Sun, 2016).
In the MMK approach, if the lag-1 autocorrelation coefficients are significantly distinct from zero, the original variance of MK statistics will be replaced by the modified one.In this study, we used the MMK approach to quantify the trends of water budget components in 18 TP basins and the significance of trend was tested at the > 95 % confidence level.
3 Results and discussion

ET evaluation and general hydrological characteristics of 18 TP basins
We first assessed the VIC_IGSNRR-simulated runoff against the observations for each basin (for example, at Tang-naihai and Pangduo stations in Fig. 2).If the Nash-Sutcliffe efficiency (NSE) coefficient between the observation and simulation is above 0.65, the VIC_IGSNRRsimulated runoff is acceptable and could be used to replace the missing runoff values for a given basin.Moreover, the CMA precipitation is consistent with TRMM (correlation coefficient (corr.)= 0.86, root mean square error (RMSE) = 8.34 mm month −1 ) and IGSNRR forcing (corr.= 0.94, RMSE = 7.15 mm month −1 ) precipitation for multiple basins (i.e. for the smallest basin above Tongren station, Fig. 2).Moreover, the magnitudes of GRACE-derived annual mean water storage change ( S) in 18 TP basins are relatively less than those for other water balance components such as annual P , Q and ET (Tables 2 and 3).The uncertainties among GRACE-derived annual mean S from different data processing centres (CSR, GFZ and JPL) are small for 18 basins except for the basins controlled by Gadatan and Tangnaihai stations.
We then evaluated 6 ET products in 18 TP basins against our calculated ET wb at a monthly basis during the period 1983-2006 (Fig. 3).The ranges of monthly averaged ET among different basins (approximately 439 mm month −1 ) are very close for all products compared to that calculated from the ET wb (642 mm month −1 ).However, GLEAM_E (corr.= 0.85 and RMSE = 5.69 mm month −1 ) and VIC_E (corr.= 0.82 and RMSE = 6.16 mm month −1 ) perform relatively better than others.Although Zhang_E and GNoah_E were found to be closely correlated to monthly ET wb in the upper Yellow River, the upper Yangtze River, and Qiangtang and Qaidam basins (Li et al., 2014), they did not exhibit overall good performances (corr.= 0.61, RMSE = 7.97 mm month −1 for Zhang_E and corr.= 0.42, RMSE = 10.16 mm month −1 for GNoah_E) for 18 TP basin used in this study.We thus use GLEAM_E and VIC_E together with ET wb to analyse the seasonal cycles and trends of ET in 18 TP basins in the following sections.
To investigate the general hydro-climatic characteristics of river basins over the TP, we classify 18 basins into 3 categories, namely westerlies-dominated basins (Yerqiang, Yulongkashi and Kelia), Indian-monsoon-dominated basins (Brahmaputra and Salween), and east-Asian-monsoondominated basins (Yellow, Yalong and Yangtze) according to Tian et al. (2007), Yao et al. (2012) and Dong et al. (2016).Interestingly, they are clustered into three groups under Budyko framework (Budyko, 1974;Zhang et al., 2016) with a relatively lower evaporative index in Indian-monsoon-dominant basins and a higher aridity index in westerlies-dominant basins, which reveal various long-term hydro-climatologic conditions (Fig. 4).Overall, in the westerlies-dominant, Indian-monsoon-dominant and east Asian monsoon-dominant basins, the annual mean air temperature (−5.68-0.97• C) and ET (and thus runoff coefficient gradually decreases) increase while the multiyear-mean glacier area (and thus the glacier melt normalized by precipitation) gradually decreases (Fig. 4 and Table 2).More-   over, the vegetation status (NDVI range: 0.05-0.43;LAI range: 0.03-0.83)tends to be better.The R 2 between basinaveraged NDVI and ET (0.76) is much higher than that between T and NDVI (0.35), which indicates that the water availability plays a more important role than the heat stress (i.e.colder status) over such basins.The results are in line with Shen et al. (2015), which indicates that the spatial pattern of ET trend is significantly and positively correlated with NDVI trend over the TP.The dominant climate systems are discrepant overall for the three TP regions with different water-energy characteristics and sources of water vapour.For example, in the westerlies-controlled basins, more glaciers developed due to their relatively colder air temperature and special seasonality of precipitation.Therefore, there are more snowmelt contributions to total river streamflow with global warming during the period 1983-2006.It is a general picture of hydrological regimes in high-altitude and cold regions (Zhang et al., 2013;Cuo et al., 2014), which could be interpreted from the perspective of multi-source datasets in the data-sparse TP.

Seasonal cycles of basin-wide water budget components for the TP basins
The multi-year means of water budget components (i.e.P , Q, ET, snow cover and SWE) and vegetation parameters (i.e.NDVI and LAI) are calculated for each calendar month and for 18 TP river basins using multi-source datasets available from 1982 to 2011.Overall, the seasonal variations of P , Q, ET, air temperature and vegetation parameters are similar in all TP basins, with peak values occurring in May to September (Figs. 5 and 6).The seasonal cycles of snow cover and SWE are generally consistent among the basins (the peak values mainly occur from October to next April, Fig. 7).With the ascending air temperature from cold to warm months, the basin-wide precipitation increases and vegetation cover expands gradually (the basin-wide ET also increase).Meanwhile, snow cover and glaciers retreat gradually with the meltwater, supplying the river discharge together with precipitation.The inter-basin variations of hydrological regime are to a large extent linked to the climate systems that prevail over the TP.
Although the temporal patterns of hydrological components are generally analogous, they vary among the parameters, climate zones and even basins (Zhou et al., 2005).For example, relative to air temperature, the seasonal pattern of runoff is similar to precipitation, which reveals that runoff is mainly controlled by precipitation in most TP basins.It is in agreement with that summarized by Cuo et al. (2014).In the westerlies-dominated basins, the peak values of precipitation and runoff are mainly concentrated in June-August, which contribute approximately 68-82 and 67-78 % of annual totals, respectively.During this period, the runoff always exceeds precipitation, which indicates large contributions of glacier and snow meltwater to streamflow.It is consistent with the existing findings in the Tarim River (Yerqiang, Yulongkashi and Keliya rivers are the major tributaries of Tarim River), which indicated that the meltwater accounted for about half of the annual total streamflow (Fu et al., 2008).The ET (vegetation cover) values in three westerliesdominated basins are relatively smaller (scarcer) than that in other TP basins while the percentages of glacier and seasonal snow cover are higher in these basins which contribute more meltwater to river discharge (Figs. 6 and 7).Overall, the SWE in Yerqiang, Yulongkashi and Keliya rivers are higher in winter than other seasons, but they vary with basins and products which reflect considerable uncertainties in SWE estimations.
In the Indian-monsoon-and east-Asian-monsoondominated basins, the runoff is concentrated during June-September (or June-October), with precipitation being the dominant contributor of annual total runoff.For example, the peak values of precipitation and runoff occur during June-September at Zhimenda station (contributing about 80 and 74 % of the annual totals) while those occur during June-October at Tangnaihai station (contributing about 78  .The ET used in this figure is calculated from the bias-corrected water balance method.and 71 % of the annual totals, respectively).The results are quite similar to the related studies in the eastern and southern TP such as Liu (1999), Dong et al. (2007), Zhu et al. (2011), Zhang et al. (2013) and Cuo et al. (2014).The vegetation cover (ET) in most basins is denser (higher) than that in the westerlies-dominant basins.Moreover, the seasonal snow mainly covers from mid-autumn to spring and correspondingly the SWE is relatively higher in these months in all basins except for the Yellow River above Xining station, Salwee River above Jiayuqiao station and Brahmaputra River above Nuxia and Yangcun stations.

Trends of basin-wide water budget components for the TP basins
The Q, P and ET wb values increased overall under regional warming during the past 30 years in the westerlies-dominated basins (Fig. 8) except for P in the Yerqiang River basin (Kulukelangan station), but only Q in Keliya River basin (Numaitilangan station) showed a statistically significant increase at the 0.05 level.The aridity index (PET/P ), which is an indicator for the degree of dryness, slightly declined (not significantly) in all basins in the northwestern TP.Although both P and PET increased in the Keliya River basin since the 1980s (Shi et al., 2003;Yao et al., 2014), the PET/P declined due to the higher rates of the increase in P than that in PET.The climate moistening (Shi et al., 2003) in the headwaters of these inland rivers would be beneficial to the water resources and oasis agro-ecosystems in the middle and lower basins.The increase in streamflow was also found in most tributaries of the Tarim River (Sun et al., 2006;Fu et al., 2010;Mamat et al., 2010).Moreover, the westerlies, revealed by the Asian zonal circulation index (60-150 • E), slightly enhanced (linear trend: 0.21, p value: 0.26) over the period 1982-2011 (Fig. 9).With the strengthening westerlies, more water vapour may be transported and fall as rain or snow in the northwestern TP (e.g. the eastern Pamir region).Both SWE products (VIC_IGSNRR-simulated and GlobaSnow-2 product) showed a marginal increase across these basins, with rising seasonal snow and glacier cover (Yao et al., 2012).More precipitation was transformed into snow or glacier cover and the runoff coefficient (Q/P ) exhibited a decrease while precipitation obviously increased (Fig. 8).In addition, the transpiration in these basins might decrease overall with vegetation degradation (Yin et al., 2016) as revealed by the NDVI and LAI (both decrease significantly in all westerliesdominated basins except for NDVI in Yerqiang and Yulongkashi rivers), but the atmospheric evaporative demand indicated by CRU PET increased (significantly increase in the Yulongkashi and Keliya rivers) during the period 1982-2011.
In the east-Asian-monsoon-dominated basins, there are two types of change for basin-wide water budget components.For example, P and Q showed marginally decrease in the upper Yellow River (Tangnihai, Huangheyan and Jimai stations) and Yalong River (Yajiang station) but slightly increased in other basins (Zelingou, Gandatan, Xining, Tongren and Zhimenda stations) over the period of 1982-2011 (Fig. 10).Only P in Zhimenda station exhibited a statistically significant increase at the 0.05 level.The decline in Q   of water budget components in westerlies-dominated (row 1), east-Asian-monsoon-dominated (rows 2-4) and Indian-monsoon-dominated (rows 5-6) TP basins.and P in the upper Yellow and Yalong Rivers (locates at the eastern Tibetan Plateau) was consistent with that found by Cuo et al. (2013Cuo et al. ( , 2014) ) and Yang et al. (2014) and were in line with the weakening east Asian summer monsoon (linear slope: −0.01, p value: 0.56; Fig. 9).The vegetation turned green markedly while ET wb and PET increased (distinct trends were found for ET wb in basins controlled by Zelingou, Tangnaihai and Zhimenda and for PET in all rivers except for basins controlled by Xining, Tongren and Zhimenda stations) in all east-Asian-monsoon-dominated basins except for ET wb in basins above Tongren and Yajing stations with the ascending air temperature during the period 1982-2011.The aridity index (PET/P ) slightly decreased in all basins (a significant decrease was found in the upper Yangtze River basin above Zhimenda station) except for the upper Yellow River basin above Jimai station and the upper Yalong River basin above Yajiang station.Moreover, the SWE showed a slight but insignificant decrease in most east-Asianmonsoon-dominated basins (SWE_VIC exhibited a marked decline in basins above Tangnaihai, Huangheyan and Zhimenda stations while SWE_Globsnow showed a significant decrease in basins above Xining station) except for the Bayin River above Zelingou station and the upper Yellow River above Tongren station.The P , ET wb and Q increased slightly in the Indianmonsoon-dominated basins (except for ET wb in the basin above Yangcun station; Fig. 11), which are in line with the strengthening (linear trend: 0.01, p value: 0.12) of the Indian summer monsoon (revealed by the Indian Ocean dipole mode index) during the specific period 1982-2011 (Fig. 9).However, only P in basins above Nuxia and Yangcun stations, Q in Nuxia station, and ET wb in basins above Jiayuqiao, Pangduo, Tangjia and Yangcun showed statistically significant trends at the 0.05 level.The slightly increasing trend of annual streamflow at Jiayuqiao station was consistent with that examined by Yao et al. (2012) during the period 1980-2000.The vegetation status, revealed by NDVI and LAI, improved slightly (marked trends were found in NDVI in basins above Gongbujiangda stations and LAI in all Indian-monsoon-dominated basins except for one above Pangduo station) when associated with the ascending air temperature.The aridity index (PET/P ) exhibited a slight but insignificant decrease in all basins (markedly declined in basins above Nuxia and Yangcun stations) except for the Brahmaputra River above Tangjia station, which indi- cated that most basins in the Indian-monsoon-dominated regions became wetter over the period of 1982-2011.The increased PET/P in Brahmaputra River basin may be consistent with the drying moisture flux in the southeastern TP, as illustrated by Gao et al. (2014).The runoff coefficient (Q/P ) increased slightly in basins above Gongbujiangda and Nuxia stations while decreasing distinctly in basins above Jiayuqiao, Tangji and Yangcun stations.Moreover, the basinwide SWE_Globsnow exhibited a minor decrease in the upper Salween River and Brahmaputra River above Tangjia and Gongbujiangda stations while increasing significantly in Brahmaputra River above Nuxia and Yangcun stations.

Uncertainties
The results may be unavoidably associated with some uncertainties inherited from the multi-source datasets used.The primary sources of uncertainty may arise from the precipitation inputs.We compared the seasonal cycles and annual trends in different precipitation products, i.e.CMA_P, IGSNRR_P and TRMM_P (and their calculated ET wb from the water balance) during the period 2000-2011 (Figs. 12 and13).We found there are some uncertainties among different precipitation products and thus among their estimated ET wb , especially in the westerlies-dominated basins.However, for each basin, the seasonal cycles of precipitation (and their calculated ET wb ) calculated from different products are similar overall (especially for the observation-based products, CMA_P and IGSNNR_P).The signs of trends for annual CMA_P and IGSNRR_P (and their calculated ET wb ) are consistent in most river basins (i.e. 14 out of 18 basins for two precipitation products and 17 out of 18 basins for their calculated ET wb ) during the period 1982-2011.The consistency of trends between two precipitation products, to some extent, revealed that the trends in CMA_P were not obviously influenced by the changing density of rain gauges in TP basins.Although some uncertainties exist due to limited and unevenly distributed meteorological stations used in the plateau and the influences of complex terrain, CMA_P is still the best observation-based precipitation product nowadays in China, which could be applied to hydrological studies in the TP.
Although the seasonal cycles of ET wb could be captured by GLEAM_E and VIC_E, they still have considerable uncertainties at some stations (e.g.Numaitilangan, Gongbujiangda and Nuxia; Fig. 5).Compared to the annual trend of ET wb (Table 4), most ET products (including the well-performed GLEAM_E and VIC_E) could not detect the decreasing trends in 7 out of 18 basins (Kulukelangan, Tongguziluoke, Xining, Tongren, Jimai, Nuxia and Gongbujiangda) due to their different forcing data, algorithms used and varied spatial-temporal resolutions (Xue et al., 2013;Li et al., 2014;Liu et al., 2016a).In particular, it is well known that land surface models have some difficulties (e.g.parameter tuning in boundary layer schemes) when applying to the TP, even though they sometimes have good performances in different regions or basins (Xia et al., 2012;Bai et al., 2016).For example, Xue et al. (2013)    The interpolation of missing values of runoff with VIC_IGSNRR-simulated runoff, and the gridded precipitation data (which interpolated from limited gauged precipitation over the plateau) also introduced uncertainties.There are also considerable uncertainties arising from the process of extending the ET series back empirically prior to the GRACE era.However, the trends in ET wb have not been significantly affected by erroneous trends in the precipitation inputs to the bias-correction based water balance calculation.For exam-   ple, the trends in CMA_P and IGSNRR_P are opposite in few basins (nos.01, 07, 08, 13 in Fig. 13), but the trends in their calculated ET wb are both consistent for each basin.It has, to some extent, certified the effectiveness of the bias-correctionbased ET-estimate approach.With these caveats, we can interpret the general hydrological regimes and their responses to the changing climate in the TP basins solely from the perspective of multi-source datasets, which are comparable to the existing studies based on the in situ observations and complex hydrological modelling.

Summary
In this study, we investigated the seasonal cycles and trends of water budget components in 18 TP basins during the period 1982-2011, which is not well understood so far due to the lack of adequate observations in the harsh environment, through integrating the multi-source global and regional datasets such as gauge data, satellite remote sensing and land surface model simulations.By using a two-step bias-correction procedure, we calculated the annual basinwide ET wb through the water balance approach, considering the impacts of water storage change.We found that the GLEAM_E and VIC_E perform better relative to other products against the calculated ET wb .
From the Budyko framework perspective, the general water and energy budgets are different in the westerliesdominated (with higher aridity index, runoff coefficient and glacier cover), the Indian-monsoon-dominated and the east-Asian-monsoon-dominated (with higher air temperature, vegetation cover and evapotranspiration) basins.In the 18 TP basins, precipitation is the major contributor to the river runoff, which is concentrated mainly during June-October (June-August for the westerlies-dominated basins, June-September or June-October for the Indian-monsoondominated and the east-Asian-monsoon-dominated basins).The basin-wide SWE is relatively high from mid-autumn to spring for all 18 TP basins except for Keliya River and Brahmaputra River above the Nuxia and Yangcun stations.There is relatively less vegetation cover, whereas there is more snow-glacier cover in the westerlies-dominant basins compared to other basins.
During the period 1982-2011, P , Q and ET wb showed a slight but insignificant increase across most of the basins in the Tibetan Plateau with the exception of some tributaries located at the upper Yellow River and Yalong River due to the weakening east Asian monsoon.The aridity index (PET/P ) exhibited an indistinctively decreasing trend in most TP basins which corresponds to the warming and moistening climate in the TP and western China.Moreover, the runoff coefficient (Q/P ) declined marginally in most basins, which may be, to some extent, due to ET increase induced by vegetation greening and the influences of snow and glacier changes.Although there are considerable uncertainties inherited from the multi-source data used, the general hydrological regimes in the TP basins could be revealed, which are consistent with the existing results obtained from in situ observations and complex land surface modelling.It indicates the usefulness of integrating the multiple datasets (e.g. in situ observations, remote-sensing-based products, reanalysis outputs, land surface model simulations and climate model outputs) for hydrological applications.The generalization here could be helpful for understanding the hydrological cycle and supporting sustainable water resources management and ecoenvironment protection in the Tibetan Plateau.

Figure 1 .
Figure 1.Map of river basins and hydrological gauging stations (green dots) over the Tibetan Plateau (TP) used in this study.The grey shading shows the topography of the TP in metres above sea level, and the blue shading exhibits the glaciers distribution in the TP extracted from the second glacier inventory dataset of China.

Figure 2 .
Figure 2. Comparison of VIC_IGSNRR-simulated and observed monthly runoff for Tangnaihai and Panduo stations (a, b) as well as (c) basin-averaged monthly TRMM, CMA gridded and IGSNRR forcing precipitations for the smallest basin (Tongren station) over the period1982 -2011.Panel (d) .Panel (d)  shows the comparison of TRMM (blue) and IGSNRR forcing (red) precipitations against CMA gridded precipitation for 18 river basins over the TP during the period 2000-2011.

Figure 3 .
Figure 3.Comparison of different ET products against the calculated ET through the water balance method (ET wb ) on the monthly timescale for 18 TP basins during the period 1983-2006.The box plot of monthly estimates of different ET products for 18 TP basins are shown in (a), while the correlation coefficients and root mean square errors (RMSEs, mm month −1 ) for each ET product relative to ET wb are exhibited in (b).

Figure 4 .
Figure 4. General water and energy status (a, the perspective of Budyko framework) and their relationships with glaciers (b) and vegetation (c, d) for 18 TP river basins.The ET used in this figure is calculated from the bias-corrected water balance method.

Figure 7 .Figure 8 .Figure 9 .
Figure7.Seasonal cycles of snow cover and snow water equivalent (SWE) in westerlies-dominated (row 1), east-Asianmonsoon-dominated (rows 2-4) and Indian-monsoon-dominated (rows 5-6) TP basins.The snow cover was extracted from a cloud-free snow composite product during the period 2005-2013.It should also be noted that the GlobSnow data are not available for some basins.

Figure 10 .
Figure10.Similar to Fig.8but for east-Asian-monsoon-dominated TP basins.It should be noted that the GlobSnow data are not available for some basins.The double red stars showed that the trend was statistically significant at the 0.05 level.

Figure 11 .
Figure11.Similar to Fig.8but for Indian-monsoon-dominated TP basins.It should be noted that the GlobSnow data are not available for some basins.The double red stars showed that the trend was statistically significant at the 0.05 level.

Figure 12 .
Figure12.Uncertainties in seasonal cycles of ETwb calculated from three precipitation products (CMA gridded, IGSNRR_Forcing and TRMM precipitation) in 18 TP basins.The comparisons were conducted during the period 2000-2011 when TRMM data were available.

Figure 13 .
Figure 13.Uncertainties in annual trends of ET wb (b) calculated from two precipitation products (CMA gridded and IGSNRR_Forcing) (a) in 18 TP basins.The comparisons were conducted during the period 1982-2011 (TRMM data were not available for the whole period).

Table 1 .
Overview of multi-source datasets applied in this study.
Liu et al.:Water budget dynamics in the Tibetan Plateau rected and biased ET, respectively.The bias correction procedure can be flexibly applied to the period 1983-2011 by matching the CDF of ET biased wb are shape and scale parameters of gamma distributions for ET biased and ET wb .ET corrected (m) and ET biased (m) represent the monthly cor-W.
indicated that GNoah_E underestimated ET wb in the upper Yellow River and Yangtze River basins on the Tibetan Plateau, mainly due to its negative-biased precipitation forcing.We thus only Hydrol.Earth Syst.Sci., 22, 351-371, 2018 www.hydrol-earth-syst-sci.net/22/351/2018/ W. Liu et al.: Water budget dynamics in the Tibetan Plateau

Table 4 .
Non-parametric trends for different ET estimates during the period 1982-2006, detected by modified Mann-Kendall test.The bold numbers show the detected trend is statistically significant at the 0.05 level.
GlobaSnow-2 SWEs have not been validated in the TP due to the lack of snow water equivalent observations, but in some basins (e.g.Zelingou and Numaitilangan) they showed similar seasonal cycles and annual trends.