Abstract

The nature of the heterogeneity of terrestrial water storage (TWS) in the Eastern Qinghai-Tibet Plateau (EQTP) is poorly understood because of the lack of validated datasets and the complex topographical conditions. In this study, monthly GRACE Level 2 Release 6 (RL06) products were employed to characterize TWS changes between April 2002 and August 2016 in the EQTP. Based on the observations and hydrological model output, the dominant factors contributing to the changes in TWS in sub-basins, and areas of TWS decrease and increase were analyzed systematically. We concluded that the TWS in the EQTP showed a slight decreasing trend from 2002 to 2016 with obvious spatial heterogeneity. The decrease in TWS may be attributed to the increase in evapotranspiration, which explains approximately 59% of the variations. In the region where a substantial decrease in TWS was observed, the trend primarily depended on evapotranspiration, and was certainly affected by glacial ablation. Moreover, the expansion of lakes supplemented by glaciers was the main cause of TWS change in the areas where TWS increased. A decrease in TWS mainly occurred in summer and was mainly due to the increase in evapotranspiration because of warming, an increase in wind speed, and a decrease in relative humidity.

HIGHLIGHTS

  • TWS in the EQTP decreased from 2002 to 2016 with the obvious spatial heterogeneity, and an increase in evaporation was the dominant factor that contributed to the decrease in TWS.

  • Lake expansion and glacial melt dominated the change in TWS in the region where a substantial increase in TWS was observed.

  • The decrease in TWS mainly occurred in summer and was mainly due to the increase in evapotranspiration in summer.

Graphical Abstract

Graphical Abstract
Graphical Abstract

INTRODUCTION

The Qinghai-Tibet Plateau, the Earth's third pole, possesses water resources that supply billions of people downstream in China, India, Pakistan and other countries in the region (Shrestha et al. 2015; Pritchard 2017; Zhao et al. 2019). Temporal and spatial changes in its water storage capacity have a major impact on the social economic development and livelihood of the region. In the context of global climate change, the ablation of ice and snow in this region is accelerating (Brun et al. 2017), and the land water cycle has undergone significant changes, particularly in the Eastern Qinghai-Tibet Plateau (EQTP). In recent years, the EQTP has attracted increasing attention from various sectors due to the increased frequency and intensity of natural disasters, ecological degradation, and the decreasing social and economic development. Frequent natural disasters are mainly represented by large-scale and long-term droughts and regional floods, which are generally linked to abnormal terrestrial water dynamics (Thomas et al. 2014; Xu et al. 2019). However, the lack of in-situ observations limits our understanding of the water cycle in this region.

Fortunately, studies have reported that the terrestrial water storage (TWS) derived from Gravity Recovery and Climate Experiment (GRACE) satellite data shows both the potential for revealing changes in surface water storage and the ability to support large-scale hydrologic modeling (Klees et al. 2007; Luo et al. 2016; Yang et al. 2017; Zhao & Li 2017; Frappart & Ramillien 2018). As mentioned before, GRACE has been used in many river basins, and has displayed satisfactory results for TWS changes (Li et al. 2013, 2018b; Ni et al. 2014; Feng et al. 2017; Khaki et al. 2018; Wang et al. 2019). A promising application can also be found for demonstration of TWS variability in the Three Gorges Reservoir area, the Guanzhong Plain and some medium-scale basins.

There are few studies on characterizing TWS changes on the EQTP based on GRACE observations. The previous studies have mostly addressed the entire plateau and partly discussed the EQTP, for example, Deng et al. (2018) and Meng et al. (2019) inferred a declining trend from 2002 to 2016 in the EQTP, but neglected the spatial heterogeneity of TWS change. Moreover, existing studies have been mostly concentrated on the annual scales of TWS change; few studies have pointed out the comprehensive mechanisms of TWS change, especially in the region with large spatial relief and limited in-situ measurements like the EQTP.

Generally, the spatio-temporal characteristics of TWS and its change are the result of the combined effect of climate change and human activities (Hu et al. 2018). Precipitation, runoff, and evapotranspiration are the dominating processes of changes in TWS either in space or time (Tangdamrongsub et al. 2011; Frappart et al. 2013; Yang et al. 2015). As an active process in rivers, the depletion of glaciers and snow cover is an influential factor in TWS dynamics (Syed et al. 2008). Anthropogenic processes like irrigation, the abstraction of groundwater, land use, etc. also impact TWS in highly populated areas (de Beurs et al. 2015; Huang et al. 2015; Khandu et al. 2016).

Several attempts have been carried out to understand the possible causes of TWS changes on the EQTP and its surrounding region. According to a recent study, the dominant factor contributing to an increase in TWS in eastern India and southeastern Tibet was increased precipitation, while the decrease was mainly attributed to a decrease in precipitation and an increase in irrigation water (Rodell et al. 2018). Another study, however, pointed out that evapotranspiration dominated the TWS depletion in the Brahmaputra river and precipitation played an important role on the TWS accumulation in the upper Yangtze. Additionally, glacier mass loss was the most likely cause of TWS decrease in the Brahmaputra river (Meng et al. 2019). Two findings showed somewhat different dominating factors in TWS change in almost the same area (most of the EQTP), which leads to uncertainty of our understanding on the cause of TWS changes. Each region has its own particular local climates, relief conditions, and human conditions, and some of these factors affect the evaluation of TWS in large-scale regions (Deng et al. 2018); for instance, in the EQTP, human activities are less intense, and their impact on TWS can be ignored, but the widely distributed glaciers and snow in this area require us to consider the substantial impacts of snow and ice meltwater on TWS. Consequently, there is still a lack of consistent understanding on the causes of TWS change. Therefore, the controlling factors for variations in TWS should be further discussed with additional validated datasets and detailed geographic division in the EQTP.

Additionally, the EQTP contains the source regions of the Yangtze, Lancang/Mekong, and Nujiang/Salween rivers. Their hydrologic dynamics are key to water resources in the downstream area. Huge fluctuations in downstream discharge have been focused on in recent years, especially in international rivers. Eliminating spatio-temporal patterns and the effect of climate features on the change of TWS may largely explain the water cycle process in source areas and the non-human impact on flow changes. Therefore, the present study is of great significance for local agriculture and animal husbandry, and ecological water demand predictions and management.

The overarching goal of this study is to characterize the spatial heterogeneity of the TWS change on the EQTP and explore a comprehensive understanding of possible causes of TWS change. We first obtain the GRACE-derived TWS in the EQTP and describe its variations by decomposing the total signal into inter-annual and seasonal variability. Second, we divide the EQTP into several parts according to its geographical conditions and present the dominant factors associated with the increase and decrease in TWS throughout the region and in its subbasins by using statistical methods and multi-source data. The main findings are presented in the Results and discussion section below.

STUDY AREA

We defined the EQTP as the region that covers the entirety of the three parallel river basins (TPRB) (90°E–101°E, 27°N–36°N) with an approximate area of 1,207,400 km2 (TPRB 399,000 km2), including the sources of the Nujiang, Lancangjiang and Jinshajiang rivers, as well as their basins in the north of Yunnan province, and part of the Big Bend Region of the Brahmaputra (Figure 1). Part of the Brahmaputra is included because there is a glacierized area and the glacier/snow runoff of the area may affect notably on the TWS change in the context of glacier retreat in high mountain Asia. The topography of the region is complicated, with an average altitude exceeding 3000 m, and the terrain fluctuates substantially as it is generally high in the north and low in the south, high in the west and low in the east, and inclines downward from the northwest to the southeast. The study area is within the influence of the southwest and the southeast monsoon. It traverses the subtropical, plateau temperate and plateau subfrozen zones from the south to the north, with substantial differences in climate. Its ecological environment and meteorological conditions are affected not only by the uplift of the Qinghai-Tibet Plateau but also by the formation and evolution of the longitudinal range of the gorge. Dry-hot valleys surrounded by a relatively humid environment with high temperature are very common in this region, such as the TPRB region (Guo et al. 2006; Luo 2009). Complex terrain and special landform, coupled with varied climatic conditions, makes the spatial and temporal variability in water storage more prominent.

Figure 1

Study area mainly represented by the TPRB and part of the Brahmaputra river basin. The locations of the weather stations are shown.

Figure 1

Study area mainly represented by the TPRB and part of the Brahmaputra river basin. The locations of the weather stations are shown.

The region has been suffering from warming and drying in recent decades (Li et al. 2018a). A rise of 0.7 °C in temperature, and a reduction of 5.5 mm in area-averaged precipitation has been noted in the last 15 years. The arid climate has caused the region to suffer from frequent environmental problems, such as water shortage, grassland degradation, land desertification and soil erosion (Li et al. 2018a). Therefore, the temporal and spatial variations in TWS must be identified as an indicator for regional water resources and their responses to climatic events.

DATA AND METHODS

Data

GRACE RL06

In this study, GRACE Level-2 (RL06) GSM (spherical harmonic coefficients of the geoid model) data produced by the Center for Space Research, University of Texas at Austin (UTCSR), were employed to calculate the TWS. The data were collected from April 2012 to August 2016, with a time resolution of one month, and the spherical harmonic coefficient of gravity was chosen to degree 60 and order 60. The GRACE monthly gravity field model had to be refined because of the error. Satellite laser ranging (SLR) C20 has higher precision and is more sensitive in reflecting seasonal characteristics than GRACE. As a result, the C20 of SLR was used to replace C20 of GRACE (Cheng & Tapley 2003; Chen & Wilson 2008). In addition, improved p3m9 decorrelation filtering was used to eliminate the north–south strip error, and 300 km fan filtering was used to reduce the high-order spherical harmonic coefficient error in the present study (Han et al. 2005; Zhen et al. 2018).

Hydrometeorological data

Soil moisture data

The Global Land Data Assimilation System (GLDAS) datasets (version 2.1), established by the Goddard Space Flight Center (NASA) and the National Centers of Environmental Prediction (NCEP) were used in the present study, which is proved to be an efficient tool for examining meteorology and hydrology (Xia et al. 2014; Yang et al. 2018; Wu et al. 2019). The GLDAS, which includes three land surface models (NOAH, CLM, and Mosaic) and one hydrological model (VIC), provides datasets (global precipitation, soil temperature, runoff, radiation flux and other hydrometeorological data) that are simulated based on multiple observation data, atmospheric assimilation products and reanalysis data for use by various researchers and relevant institutions (Rodell et al. 2004). The four-layer soil moisture data of NOAH have been widely used in the inversion of GRACE-derived TWS (Li et al. 2018b; Zhen et al. 2018; Zhong et al. 2018). The data derived from NOAH from 2002 to 2016 were employed in this study, with a spatial resolution of 1 × 1° and a temporal resolution of one month, and they were mainly used in this study to calculate the scale factor and characterize the soil moisture in the TWS.

Precipitation, evapotranspiration and runoff data

We employed GLDAS precipitation, evapotranspiration, and runoff to depict the water balance to minimize the error caused by different data sources. As the largest input to the TWS, precipitation data are crucial for identifying the regional water balance process. Thus, we evaluated the accuracy of GLDAS precipitation here. The accuracy was validated by comparing it with the station data using R-square, bias, root mean square error (RMSE) and relative RMSE (rRMSE) (Junzhi et al. 2012). As shown in Figure 2, GLDAS precipitation displayed good correlation with the station data (R2 = 0.45), which generally reflects the characteristics of regional precipitation. However, from the perspective of RMSE and rRMSE, GLDAS was not close to the actual observed value. In addition, GLDAS slightly underestimates the precipitation. Considering that we cannot adequately quantify the errors owing to different data sources, and the GLDAS-derived TWS using the water balance method which was quite consistent with GRACE derived TWS (Figure 3), we suggest that GLDAS precipitation basically satisfies the requirements for analysis, although we also tested precipitation from CN501 (Jia & Xuejie 2013), TRMM.3b43v7 (Huffman et al. 2007), and GPCP (Adler et al. 2018), which exhibited better results than GLDAS in the EQTP.

Figure 2

Comparative analysis of GLDAS precipitation and measured data.

Figure 2

Comparative analysis of GLDAS precipitation and measured data.

Figure 3

(a) GRACE derived and GLDAS derived TWS time series, and (b) trend, seasonality, (c) and residuals of GRACE derived TWS. The GLDAS derived TWS was obtained using the water balance method and the Pearson correlation coefficient between two TWS results was 0.84 (p < 0.0001).

Figure 3

(a) GRACE derived and GLDAS derived TWS time series, and (b) trend, seasonality, (c) and residuals of GRACE derived TWS. The GLDAS derived TWS was obtained using the water balance method and the Pearson correlation coefficient between two TWS results was 0.84 (p < 0.0001).

Groundwater data

Observed wells data can accurately reflect groundwater reserves, but the number of wells is small and data acquisition is difficult. In recent years, data from the WaterGAP Global Hydrology Model (WGHM) simulated by IPG (University of Frankfurt) have better reflected changes in groundwater reserves. In addition, the model also reflects the changes in all forms of water, except glaciers. Currently, the results from this model have been widely used (Eicker et al. 2014; Feng et al. 2018). The groundwater data from 2002 to 2013 were used in the present study, with a spatial resolution of 0.5 × 0.5° and a temporal resolution of one month (Döll et al. 2015).

Measured data

To improve the validity and accuracy of the results, we derived the precipitation, evapotranspiration, wind speed, wind direction, temperature, humidity and other observation data recorded by 36 stations in the study area, obtained from the national meteorological center (http://data.cma.cn/) (Figure 1). After preprocessing the observation data, the time series of the data recorded by the 36 stations from 2002 to 2016 were obtained. The measured data from each station were mainly used to verify the accuracy of the data obtained from the model and evaluate the effects of meteorological factors on TWS.

Seasonal-trend decomposition procedure based on Loess (STL)

Since the data used in this study are time series data, they can be decomposed to obtain details and to clarify the specific characteristics of changes in timing data. As a commonly used time-series decomposition method, STL is widely used in studies related to meteorology and hydrology. It divides the time series Y into a trend component T, periodic component (seasonal component) S, and residual term R to refine and characterize variations in time series data at different levels (Cleveland & Cleveland 1990; Sanchez-Vazquez et al. 2012; Rojo et al. 2017). For time series Y, the STL decomposition at time i is calculated using Equation (1). More detailed information related to STL can be found in the description of Cleveland & Cleveland (1990):
formula
(1)

In the present study, STL is mainly used to decompose time series data, such as TWS, and to analyze the characteristics of each component and extract effective information.

Inversion of TWS

One of the main explanations for the change in the spherical gravity field is the dynamic changes in the terrestrial water mass (Tapley et al. 2004), which is expressed in the form of equivalent water thickness (EWT) using the equation described by Wahr et al. (1998). More information about TWS inversion is described in the first section of the Supplementary material. Errors may occur in a series of postprocessing operations of the GRACE gravity field data (Klees et al. 2007), such as spherical harmonic expansion and spatial filtering (decorrelation filtering, fan filtering, etc.), which cause leakage of the signal (Feng et al. 2017). In this study, we used double scale factors to reduce the leakage error and acquired satisfied results. A detailed description of the processes can be found in the second section of the Supplementary material.

Water balance

In the vertical direction, regional TWS includes soil moisture (SM), groundwater (GW), surface waterbodies (SW) (such as lakes and reservoirs), melting water (including snow water equivalent (SWE) and glacier melting water (G)), and canopy water content (CW) (Cao et al. 2015). Considering that the hydropower stations of Nujiang and Lancangjiang are mainly concentrated in the middle and lower reaches of the basins, reservoir filling in the study area does not exert an effect on the results. Although a plan to establish eight hydropower stations in the upper reaches of the Jinshajiang has been proposed, only three stations had been implemented as of 2016, and their scale is relatively small. For the Brahmaputra, there is no large reservoir located in the study area. Therefore, the impacts of reservoirs in the study area have been ignored. In this study, the changes in TWS were roughly measured by calculating changes in SM, GW, SW, CW, SWE, and G. Furthermore, from the perspective of the water cycle, the variation in TWS can also be denoted by changes in precipitation (), evapotranspiration ( and runoff (Yang et al. 2014; Chao & Wang 2017). Ignoring some subtle hydrological processes, was calculated using Equation (2). As Meng et al. (2019) described, the total change rate of TWS () can be expressed by the sum total of the change rate of its components () and the contribution of each variable (SM, GW, SW, SWE, G) to the total TWS changes is expressed as :
formula
(2)

RESULTS AND DISCUSSION

The overall trend of TWS in the EQTP

Using the TWS inversion methods described above under ‘Inversion of TWS’, the GRACE data were processed to obtain the water storage equivalent water thickness of the EQTP. Figure 3(a) shows the temporal variations in the regional averaged monthly TWS anomalies. Two TWS results are highly correlated, which indicates the strength of the linear relationship between the GRACE-derived TWS and GLDAS-derived TWS; the coefficient of determination R2 (hereafter R2; here R2 = 0.84 before removing the seasonal cycle, but R2 = 0.52 after removing the seasonal cycle) shows we were able to calculate a non-seasonal variation in TWS of approximately 52% in the EQTP using GLDAS. Therefore, we could evaluate controlling factors of TWS change by analyzing the changes of TWS components derived from GLDAS.

As shown in Figure 3, the TWS showed a decreasing trend from 2002 to 2016, with an approximate rate of –0.15 mm/month, but it showed an increasing trend from 2002 to 2005, with an approximate rate of 0.89 mm/month. The time series was decomposed to the trend, seasonal and residual items with STL to refine the characteristics of TWS. The seasonal variation in TWS was distinguished from the decomposed signal, which showed substantial annual losses in winter and spring and the largest annual surpluses in summer and autumn. Its amplitude was maintained at ±50 mm, and no significant inter-annual variation in the seasonal characteristics was observed. The characteristics of the changes in the trend term were basically consistent with the overall change in TWS, indicating that the change in TWS mainly depended on the trend term. The residual term of TWS reflected the characteristics of TWS fluctuation and was the most powerful evidence for drought or flood events in the region. No obvious pattern of fluctuation in TWS was observed, but some obviously abnormal signals were detected, such as in the spring of 2003 and 2004, and in the winter of 2006 and 2015, etc. We counted the annual fluctuations in TWS to capture the surplus and deficits in TWS, and the results are shown in Table 1. In general, except for certain years, the loss of water storage capacity is greater than 35 mm, resulting in a drying condition; however, the surplus of water reserves is maintained at greater than 30 mm. These findings are consistent with the characteristics of a dry-hot valley. Regarding the balance of TWS, losses were observed for seven of the years between 2003 and 2015, particularly in 2004, 2006 and 2011. According to the statistics from the China Meteorological Disaster Yearbook issued by the China Meteorological Administration, the drought in northern Yunnan and southeastern Tibet was severe in these years, and agriculture and animal husbandry were substantially affected. Therefore, we confirmed that the losses of TWS in these years were related to droughts. From the peak of TWS surplus and deficit, except for 2003, 2004, 2008, 2010 and 2015, the largest losses in other years appeared in the summer and autumn, results that were completely opposite to the seasonal variation in TWS. We speculated that, in these years, the study area experienced a certain degree of drought. More interestingly, in 2012 and 2014, the maximum deficit and maximum surplus appeared only one month apart, which may have been caused by the lagged effect of heavy precipitation on regional soil moisture supplementation.

Table 1

Characteristics of fluctuation in TWS

YearSurplus (mm)Deficit (mm)Budget (mm)Max surplus dateMax deficit date
2003 52.47 42.29 10.19 2003–06 2003–04 
2004 37.04 57.37 –20.33 2004–07 2004–03 
2005 77.00 47.87 29.13 2005–05 2005–07 
2006 59.49 92.48 –32.98 2006–12 2006–08 
2007 40.19 44.85 –4.66 2007–02 2007–07 
2008 53.39 19.75 33.64 2008–11 2008–01 
2009 39.58 48.67 –9.10 2009–08 2009–06 
2010 42.57 43.04 –0.47 2010–06 2010–02 
2011 30.49 50.38 –19.90 2011–06 2011–09 
2012 49.15 37.38 11.77 2012–08 2012–07 
2013 52.56 36.87 15.69 2013–08 2013–10 
2014 42.84 42.45 0.39 2014–07 2014–06 
2015 58.58 65.55 –6.97 2015–10 2015–12 
YearSurplus (mm)Deficit (mm)Budget (mm)Max surplus dateMax deficit date
2003 52.47 42.29 10.19 2003–06 2003–04 
2004 37.04 57.37 –20.33 2004–07 2004–03 
2005 77.00 47.87 29.13 2005–05 2005–07 
2006 59.49 92.48 –32.98 2006–12 2006–08 
2007 40.19 44.85 –4.66 2007–02 2007–07 
2008 53.39 19.75 33.64 2008–11 2008–01 
2009 39.58 48.67 –9.10 2009–08 2009–06 
2010 42.57 43.04 –0.47 2010–06 2010–02 
2011 30.49 50.38 –19.90 2011–06 2011–09 
2012 49.15 37.38 11.77 2012–08 2012–07 
2013 52.56 36.87 15.69 2013–08 2013–10 
2014 42.84 42.45 0.39 2014–07 2014–06 
2015 58.58 65.55 –6.97 2015–10 2015–12 

Spatial differences in TWS in the EQTP

To understand the spatial differentiation of TWS in the EQTP, we next performed a grid-to-grid process to obtain the long-term trend, which is defined as the monthly average rate of TWS change. As shown in Figure 4, overall, the spatial variation rate of TWS ranged from –0.94 to 0.91 mm/month (the loss in local areas exceeded 1.4 mm/month). The spatial differences were significant. The southwest direction showed a large decreasing trend (defined as R1, the blue line frame in Figure 4), and this direction is the main flow of the Brahmaputra River, and the north direction showed an increasing trend (defined as R2, the yellow line frame in Figure 4), particularly in the area recognized as the source of the Jinshajiang Basin. TWS in the central and eastern parts of the study area (defined as R3, between the blue line and the yellow line in Figure 4) changed slightly and generally maintained a balance. In addition, TWS in the middle and lower reaches of the TPRB experienced a significant loss that was largely related to the unique topography of a typical dry-hot valley. This loss has led to frequent large-scale drought events in recent years. In the TPRB, the Jinshajiang Basin showed an increasing trend overall, with a rate of 0.19 mm/month, and a weak downward trend was observed downstream of the basin. The headstream of the Lancangjiang Basin showed a slight increase in TWS; nevertheless, the TWS of the whole basin generally decreased at a rate of 0.28 mm/month. The Nujiang Basin experienced a vast deficit in TWS and its negative change rate reached 0.49 mm/month.

Figure 4

Map of the spatial distribution of changes in TWS in the EQTP. The study area is divided into three parts by the yellow line and the blue line, and the parts are labeled R2, R3 and R1 from north to south, respectively. Please refer to the online version of this paper to see this figure in colour: http://dx.doi.org/10.2166/nh.2020.039.

Figure 4

Map of the spatial distribution of changes in TWS in the EQTP. The study area is divided into three parts by the yellow line and the blue line, and the parts are labeled R2, R3 and R1 from north to south, respectively. Please refer to the online version of this paper to see this figure in colour: http://dx.doi.org/10.2166/nh.2020.039.

Partitioning components in TWS dynamics

We used soil moisture (SM), snow water equivalent (SWE) and total amount of water within canopy (CW) from the GLDAS-NOAH model and groundwater (GW) from WGHM model to further explore the characteristics of the variations in each component of TWS. As shown in Figure 5, the regional GW does not display an obvious change trend, consistent with the results reported by Feng et al. (2018), who showed that GW was not the main component contributing to the regional change in TWS. CW shows little contribution to the TWS change and tendency. The change in the SWE tends to be 0, which is mainly because the snow and ice is mainly concentrated in the upstream areas of the basins, and the area covered by ice and snow is much smaller than the total area of the EQTP. Regarding the whole area, the contribution of the SWE to changes in TWS is small, but the amount of change is generally greater than 0 (approximately 0–2 mm), particularly after 2011, indicating that the snow melt in basins is increasing. A continuous time series data representation for the change in glacier meltwater was not available. However, according to previous studies, the glacier in the EQTP has shown a negative mass balance in recent years, namely, the glacier melt water is increasing (Kääb et al. 2015; Berthier et al. 2016; Zhou et al. 2018). SM is the main component of TWS, which is a necessary condition for the survival and growth of regional plants. Drought occurs when SM is unable to meet the needs of plant growth. In 2009–2014, southwestern China, particularly the Yunnan-Guizhou Plateau, continued to experience a drought (Rong et al. 2018), and the drought was the most severe in 2010. As shown in Figure 5, the extent to which the SM decreased during this period was large. Consequently, we inferred that the drought during this period was mainly caused by insufficient soil moisture. Notably, in the southwestern drought in 2010, the reduction in TWS was large, while the SM did not fluctuate substantially, and a rapid decrease in SM occurred in early 2011. This difference may have been caused by the lag effects of P and ET on soil moisture. The decrease in the soil water content in 2011–2015 was greater than TWS. According to Zhang et al. (2015a), an increase in the surface temperature enhanced the regional ET, so the drought was strengthened during the period. At the same time, due to the increase in temperature, the melting of the glaciers was accelerated. Thus, we speculate that the acceleration of glacial ablation and the increase in ET partially contributed to the changes in TWS in the southwest of the region.

Figure 5

The changes (a), statistical characteristics (b), and trend (c) of TWS and its components including soil moisture (SM), groundwater (GW), total amount of water within canopy (CW), and snow water equivalent (SWE).

Figure 5

The changes (a), statistical characteristics (b), and trend (c) of TWS and its components including soil moisture (SM), groundwater (GW), total amount of water within canopy (CW), and snow water equivalent (SWE).

Attributes of spatial heterogeneity in TWS

In this section, we analyze the causes of spatial differences in water reserves based on changes in precipitation (P), evapotranspiration (ET), and runoff (R). The TWS in R3 remained stable and will not be discussed here. We mainly analyzed the entire study area, zones R1 and R2.

Based on the results described above, in the study area, the variability in ET, P and R was calculated as 0.0180, −0.0107, and −0.0326 mm/month, respectively (Figure 6(a)). The general increase in ET and the decrease in P agreed with the negative TWS change with high correlation coefficient of 0.77 (p < 0.0001) and 0.66 (p < 0.0001), respectively. In terms of the variability rates, ET was nearly twice as high as P, which showed that ET more sensitive to TWS change. Additionally, the R2 between ET and TWS was approximately 0.59, which indicated that ET could explain 59% of the variations of TWS. R integrated the effects of P and ET and it was mainly supplied by glacier melting water in this area, and thus it made a limited contribution to the changes in TWS, although the rate of R was even larger than that of ET. Thus, the increase in ET was a primary driver for water storage dynamics across the EQTP during 2002–2016.

Figure 6

Time series of P, ET and R for the whole study area (a), R1 (b) and R2 (c). ‘k’ (mm/month) represents the change rate from 2002 to 2016.

Figure 6

Time series of P, ET and R for the whole study area (a), R1 (b) and R2 (c). ‘k’ (mm/month) represents the change rate from 2002 to 2016.

R1 displayed a severe TWS deficit, with a loss rate of −0.71 mm/month. As shown in Figure 6(b), the P and ET exhibit higher variability in R1 than in the study region as a whole, and the P shows a positive tendency. This positive trend in P means that the loss of TWS was mostly caused by a decrease in ET. Similarly, we noted that 54% of the variations in TWS could be explained by ET by calculating the R2. Figure 6(b) shows that the regional R displays a decreasing trend, which is inconsistent with the variation of P. Thus, we can infer that the water supply from glaciers and snow meltwater was decreasing. Based on the results reported by Brun et al. (2017), the mass balance of representative glaciers in the study area was −0.62 ± 0.23 m w.e. yr–1, we concluded that the glacial melt variability in R1 was approximately −0.64 mm/month (the glacier area is 4,368.836 km2, and the water equivalent of the glacial melt is 620 mm/year). The effect of glacier change on TWS change is indeed not negligible, although we cannot quantify the contribution of glacier meltwater to R in the present study. Therefore, our results suggest that glacier mass loss exerts a significant role in the variation of TWS apart from ET in the R1.

R2 showed an increase in TWS at a rate of 0.49 mm/month. As shown in Figure 6(c), the regional P, R, and ET show changes that are inconsistent with TWS, although the extent of these changes is small, which means that P and ET are not the key factors that affect change in TWS. Figure 1 shows the presence of many large lakes in R2. According to the relevant research, the water volume of lakes in the Qinghai-Tibet Plateau has increased significantly in recent years, and many new lakes have appeared (Li et al. 2019). Therefore, the contribution of the increase in lake water storage to the change in TWS must be considered. Currently, detailed data regarding the amount of change in the water storage capacity of these lakes are not available, thus we can only speculate that the increase in TWS observed in this area is closely related to lake expansion. The recharge source of these lakes is mainly the glaciers in R2 and those located outside the study area, which indirectly shows that glacial changes contribute to the increase in TWS. In addition, in recent years, the depth of the frozen soil active layer on the Tibetan Plateau has significantly increased, and the freezing depth has decreased significantly (Zhang et al. 2015b; Peng 2018). From this perspective, the transfer of mass from the ablation of frozen soil around the study area may also affect the changes in TWS.

To obtain a better understanding of the effect of lake expansion and glacial melt on TWS, we calculated the contribution of each variable. As shown in Figure 7(a), CW and SWE exhibit a slight influence on TWS change and SM dominates the TWS change before 2012. was estimated at 0.35 mm/month from 2002 to 2013 and showed a contribution of approximately 54% to TWS change, which verifies the speculation above. It is noteworthy that GLDAS-derived TWS was quite different from GRACE-derived TWS after 2012 (Figure 7(b)) and the SM showed a continuous downward trend after 2012. The decrease in SM was consistent with the increase in ET and the decrease in R. This suggests that the contribution of to the change in TWS was increasing.

Figure 7

(a) Soil moisture (SM), groundwater (GW), total amount of water within canopy (CW), snow water equivalent (SWE) and glacier melting water (G) change rates in R2. (b) GRACE and GLDAS derived TWS in R2. GLDAS derived TWS was defined as the sum of ΔSM, ΔCW, ΔSWE and ΔGW (obtained from WGHM). The k_* (mm/month) excluded k_tws(GRACE) (2002–2016) represents the change rate from 2002 to 2012.

Figure 7

(a) Soil moisture (SM), groundwater (GW), total amount of water within canopy (CW), snow water equivalent (SWE) and glacier melting water (G) change rates in R2. (b) GRACE and GLDAS derived TWS in R2. GLDAS derived TWS was defined as the sum of ΔSM, ΔCW, ΔSWE and ΔGW (obtained from WGHM). The k_* (mm/month) excluded k_tws(GRACE) (2002–2016) represents the change rate from 2002 to 2012.

Factors affecting the decrease in TWS in summer

Based on the analysis described above, the increase in ET was the main cause of the decrease in TWS in the EQTP during the study period. Additionally, the decrease in P also contributed to the decrease to some extent. Generally, of all seasons, the highest levels of ET were observed in the northern regions during the summer (defined as June, July and August, and represented by Q3) (Su & Feng 2015), when the study area was also significantly affected by the Indian summer monsoon. Therefore, in this section, we discuss the changes in TWS in the northern region during the summer and the factors that substantially influence ET. Namely, whether the variations in air temperature, air humidity and wind speed converge with the characteristics of change in ET. We will also analyze the variations in the Indian summer monsoon and its feedback effect on P to analyze the factors that affect the change in TWS in the northern region during the summer.

As shown in Figure 8, it was observed that the water deficit in summer dominates the TWS annual trend in the northern region. Meanwhile, P shows a decreasing trend that converges with the varying trends of TWS. Although the rate of decrease in P is small, as the main component of the regional water reserves this small decreasing rate slows the loss of TWS. P varies substantially in regions with complex topographic conditions, but the long-term tendencies of P are mainly affected by atmospheric conditions, and topographic and geomorphic features are amplifiers of atmospheric conditions. During summer, the study area is mainly affected by the Indian summer monsoon (Yao et al. 2017), and thus we used the Indian Monsoon Index (IMI) (Wang & Fan 1999) to measure the relationship between the Indian summer monsoon and P during the study period. In Figure 9, the IMI exhibits a decreasing trend, indicating that the IMI is weakening and the moisture brought by the Indian summer monsoon is decreasing during this period, which is one reason for the overall decrease in regional P. After the summer of 2005, the Indian summer monsoon continues to weaken, but the P only decreases slightly with a larger interannual fluctuation in general, which is particularly large from 2011 to 2015. Thus, the weakening of the Indian summer monsoon is not the main cause of the decrease in P in the EQTP. IMI is more consistent with the change in TWS than the change in P. The correlation coefficient between the TWS change and IMI is 0.53, which indicates that the TWS loss in summer is partly affected by the Indian summer monsoon. To quantify the relationship, we calculated the amplitude of effect of the Indian summer monsoon on TWS using a method that serves to characterize the association between ENSO and TWS (Phillips et al. 2012) (see ‘Quantification of the relationship between TWS change and Indian Monsoon Index’ in the Supplementary Material). The amplitude of 0.018 mm/(IMI*year) suggested that the Indian summer monsoon exerts a small influence on the yearly change in TWS. The nonconformity between the weakness of the Indian summer monsoon and the decrease in P requires further research and we will not delve into this topic in the present study. However, we speculate that the increases in summer temperature and ET (Figure 9), together with the changes in underlying surface in the region (the decrease in vegetation coverage caused by human activities and the decrease in glacier cover area caused by warming in the region, etc.), have resulted in a decrease in the atmospheric water vapor content, which may be another explanation for the regional decrease in P.

Figure 8

Changes in TWS, precipitation and IMI in summer. The k_pre (mm/month) and k_tws (mm/(3*month)) represent the change rate of precipitation and TWS, respectively.

Figure 8

Changes in TWS, precipitation and IMI in summer. The k_pre (mm/month) and k_tws (mm/(3*month)) represent the change rate of precipitation and TWS, respectively.

Figure 9

Changes in evapotranspiration (ET), relative humidity, temperature and wind speed in summer. The k represents the change rate of each variable. The unit of k depends on the unit of the variable.

Figure 9

Changes in evapotranspiration (ET), relative humidity, temperature and wind speed in summer. The k represents the change rate of each variable. The unit of k depends on the unit of the variable.

As shown in Figure 9, ET is the main factor contributing to the decrease in TWS, and it increases significantly in summer, showing a higher trend than is observed throughout the year. Based on these findings, the decrease in TWS that is mainly observed in summer is generally controlled by the increase in ET. To clarify the reasons for the increasing ET in summer, we selected the meteorological elements that had a greater impact on ET in the region and analyzed their changes. Generally, during the study period, the air temperature increase, the wind speed increase, and the variation in air humidity are consistent with the variation in ET in the region. In terms of the rate of change, the change in relative humidity has the greatest contribution to ET, followed by the wind speed. However, the relative humidity is also affected by the changes in temperature and air pressure. Therefore, the extent of the contributions of the three factors to ET is difficult to determine. Nevertheless, in terms of fluctuation, we determined that changes in temperature were the most consistent with changes in ET, indicating that changes in temperature are the factor with the greatest direct effect on ET. As shown in Figure 9, the change in regional ET is consistent with the change in TWS, but the temperature is consistent with TWS after a lag. Regarding the points at which ET and temperature changed in the time series data, asynchronous correlations were observed between the two variables as well, indicating that the temperature change exerted a certain lag effect on ET.

CONCLUSIONS

Overall, from 2002 to 2016, the TWS in the EQTP shows a slight decreasing trend, and the spatial differences were significant. In the years when droughts occurred, TWS fluctuated abnormally. Soil moisture (SM) was the main component of TWS. The drought in the EQTP was mainly due to the continuous decrease in SM. Groundwater (GW) fluctuated substantially during the dry season, particularly in the years of persistent droughts.

The increase in evapotranspiration (ET) was the dominant factor contributing to the decrease in TWS in the EQTP (approximately 59% of the variations in TWS can be explained), while the decrease in precipitation (P) had a certain effect on TWS change as well. In the area R1, 54% of the variation in regional TWS could be explained by the variation in evapotranspiration. In addition, glacier ablation played a significant role in the loss of TWS, although its contribution could not adequately be quantified. In the area where water reserves were increasing (R2), lake expansion and glacial melt played a leading role in changes in TWS. Furthermore, frozen soil around the study area may also have affected the variation in TWS. exhibited a contribution of approximately 54% to the increase in TWS from 2002 to 2013, and its effect was increasing.

The decrease in TWS mainly occurred in summer. The increase in ET in summer, which was the result of an increase in temperature, an increase in wind speed and a decrease in relative humidity, was the main cause of the decrease in TWS. The Indian monsoon exerted a certain effect on the decrease in P in summer, but understanding the dominant factor that contributed to the decrease requires further study.

ACKNOWLEDGEMENTS

The authors thank members of the Institute of International Rivers and Eco-Security, Yunnan University, plus Fuming Xie, Yongpeng Gao, Wenfei Miao, Xinxin Qiang and Shimei Duan for their input and stimulating discussions over past few months. We also thank Dr Wei Feng at Institute of Geodesy and Geophysics, Chinese Academy of Sciences, for helpful comments on this manuscript. This work is supported by the National Natural Science Foundation of China (Nos. 41761144075 and 41661144044), the Second Tibetan Plateau Scientific Expedition and Research Program (STEP, No. 2019QZKK0208), the Research Funds for Introducing Talents of Yunnan University (No. YJRC3201702) and Innovation Fund Designated for Graduate Students of Yunnan University (Grant No. 2018Z099).

DATA AVAILABILITY STATEMENT

Data cannot be made publicly available; readers should contact the corresponding author for details.

REFERENCES

REFERENCES
Adler
R. F.
,
Sapiano
M.
,
Huffman
G. J.
,
Wang
J.
,
Gu
G.
,
Bolvin
D.
,
Chiu
L.
,
Schneider
U.
,
Becker
A.
,
Nelkin
E.
,
Xie
P.
,
Ferraro
R.
&
Shin
D. B.
2018
The global precipitation climatology project (GPCP) monthly analysis (New version 2.3) and a review of 2017 global precipitation
.
Atmosphere (Basel)
9
,
138
.
Brun
F.
,
Berthier
E.
,
Wagnon
P.
,
Kaab
A.
&
Treichler
D.
2017
A spatially resolved estimate of High Mountain Asia glacier mass balances, 2000–2016
.
Nature Geoscience
10
,
668
673
.
Chen
J. L.
&
Wilson
C. R.
2008
Low degree gravity changes from GRACE, Earth rotation, geophysical models, and satellite laser ranging
.
Journal of Geophysical Research Solid Earth
113
,
B06402
.
Cheng
M. K.
&
Tapley
B. D.
2003
Variations in the Earth's oblateness during the past 26 years
.
Journal of Geophysical Research Solid Earth
109
,
1404
1406
.
Cleveland
R. B.
&
Cleveland
W. S.
1990
STL: a seasonal-trend decomposition procedure based on loess
.
Journal of Official Statistics
6
,
3
33
.
de Beurs
K. M.
,
Henebry
G. M.
,
Owsley
B. C.
&
Sokolik
I.
2015
Using multiple remote sensing perspectives to identify and attribute land surface dynamics in Central Asia 2001–2013
.
Remote Sensing of Environment
170
,
48
61
.
Feng
W.
,
Wang
C. Q.
,
Mu
D. P.
,
Zhong
M.
,
Zhong
Y. L.
&
Xu
H. Z.
2017
Groundwater storage variations in the North China Plain from GRACE with spatial constraints
.
Chinese Journal of Geophysics
60
,
1630
1642
.
Frappart
F.
,
Ramillien
G.
&
Ronchail
J.
2013
Changes in terrestrial water storageversusrainfall and discharges in the Amazon basin
.
International Journal of Climatology
33
,
3029
3046
.
Guo
J.
,
Wang
Z.
,
Bai
B.
&
Lu
Y.
2006
The analysis of abrupt climate change caused by meteorological stations moving in Yunnan
.
Yunnan Geographic Environment Research
18
,
48
52
.
Han
S. C.
,
Shum
C. K.
,
Jekeli
C.
,
Kuo
C. Y.
,
Wilson
C.
&
Seo
K. W.
2005
Non-isotropic filtering of GRACE temporal gravity for geophysical signal enhancement
.
Geophysical Journal International
163
,
18
25
.
Hu
W.
,
Liu
H.
,
Bao
A.
&
El-Tantawi
A. M.
2018
Influences of environmental changes on water storage variations in Central Asia
.
Journal of Geographical Sciences
28
,
985
1000
.
Huffman
G. J.
,
Bolvin
D. T.
,
Nelkin
E. J.
,
Wolff
D. B.
,
Adler
R. F.
,
Gu
G.
,
Hong
Y.
,
Bowman
K. P.
&
Stocker
E. F.
2007
The TRMM Multisatellite Precipitation Analysis (TMPA): Quasi-global, multiyear, combined-sensor precipitation estimates at fine scales
.
Journal of Hydrometeorology
8
,
38
55
.
Jia
W.
&
Xuejie
G.
2013
A gridded daily observation dataset over China region and comparison with the other datasets
.
Chinese Journal of Geophysics
56
,
1102
1111
.
Khaki
M.
,
Forootan
E.
,
Kuhn
M.
,
Awange
J.
,
Dijk
V.
,
Schumacher
A. I. J. M.
&
Sharifi
M. A. M.
2018
Determining water storage depletion within Iran by assimilating GRACE data into the W3RA hydrological model
.
Advances in Water Resources
114
,
1
18
.
Klees
R.
,
Zapreeva
E. A.
,
Winsemius
H. C.
&
Savenije
H. H. G.
2007
The bias in GRACE estimates of continental water storage variations
.
Hydrology and Earth System Sciences
11
,
1227
1241
.
Li
Q.
,
Luo
Z. C.
,
Zhong
B.
&
Wang
H. H.
2013
Terrestrial water storage changes of the 2010 southwest China drought detected by GRACE temporal gravity field
.
Chinese Journal of Geophysics
56
,
1843
1849
.
Li
T.
,
Wang
B.
,
Chang
C. P.
,
Zhang
Y.
&
Li
P. T.
2018a
A theory for the Indian Ocean dipole mode
.
Journal of the Atmospheric Sciences
60
,
2119
2135
.
Li
W.
,
Wang
W.
,
Zhang
C. Y.
,
Yang
Q.
,
Feng
W.
&
Liu
Y.
2018b
Monitoring groundwater storage variations in the Guanzhong area using GRACE satellite gravity data
.
Chinese Journal of Geophysics
6
,
2237
2245
.
Luo
Y.
2009
Research on Geologic Environment in Three Parallel Rivers Region
.
China University of Geosciences
,
Beijing
.
Meng
F.
,
Su
F.
,
Li
Y.
&
Tong
K.
2019
Changes in terrestrial water storage during 2003–2014 and possible causes in Tibetan Plateau
.
Journal of Geophysical Research: Atmospheres
124
,
2909
2931
.
Ni
S.
,
Chen
J.
,
Jin
L.
,
Chao
C.
&
Liang
Q.
2014
Terrestrial water storage change in the Yangtze and Yellow River Basins from grace time-variable gravity measurements
.
Journal of Geodesy & Geodynamics
34
,
49
54
.
Peng
X.
2018
Spatial-Temporal Variations of Seasonally Frozen Ground and its Response to Climate Change in the Northern Hemisphere
.
Lanzhou University
,
Lanzhou
,
China
.
Phillips
T.
,
Nerem
R. S.
,
Fox-Kemper
B.
,
Famiglietti
J. S.
&
Rajagopalan
B.
2012
The influence of ENSO on global terrestrial water storage using GRACE
.
Geophysical Research Letters
39
,
L16705
.
Rodell
M.
,
Houser
P. R.
,
Jambor
U.
,
Gottschalck
J.
,
Mitchell
K.
,
Meng
C. J.
,
Arsenault
K.
,
Cosgrove
B.
,
Radakovich
J.
,
Bosilovich
M.
,
Entin
J. K.
,
Walker
J. P.
,
Lohmann
D.
&
Toll
D.
2004
The global land data assimilation system
.
Bulletin of the American Meteorological Society
85
,
381
394
.
Rodell
M.
,
Famiglietti
J. S.
,
Wiese
D. N.
,
Reager
J. T.
,
Beaudoing
H. K.
,
Landerer
F. W.
&
Lo
M. H.
2018
Emerging trends in global freshwater availability
.
Nature
557
,
651
659
.
Rojo
J.
,
Rivero
R.
,
Romero-Morte
J.
,
Fernandez-Gonzalez
F.
&
Perez-Badia
R.
2017
Modeling pollen time series using seasonal-trend decomposition procedure based on LOESS smoothing
.
International Journal of Biometeorology
61
,
335
348
.
Rong
Y. S.
,
Gong
L.
&
Lu
S.
2018
Analysis on characteristics and causes of persistent meteorological and hydrological drought in Yunnan from 2009 to 2014
.
Water Resources Protection
34
,
22
29
.
Shrestha
M.
,
Koike
T.
,
Hirabayashi
Y.
,
Xue
Y.
,
Wang
L.
,
Rasul
G.
&
Ahmad
B.
2015
Integrated simulation of snow and glacier melt in water and energy balance-based, distributed hydrological modeling framework at Hunza River Basin of Pakistan Karakoram region
.
Journal of Geophysical Research: Atmospheres
120
,
4889
4919
.
Syed
T. H.
,
Famiglietti
J. S.
,
Rodell
M.
,
Chen
J.
&
Wilson
C. R.
2008
Analysis of terrestrial water storage changes from GRACE and GLDAS
.
Water Resources Research
44
,
W02433
.
Tapley
B. D.
,
Bettadpur
S.
,
Watkins
M.
&
Reigber
C.
2004
The gravity recovery and climate experiment: mission overview and early results
.
Geophysical Research Letters
31
,
L09607
.
Thomas
A. C.
,
Reager
J. T.
,
Famiglietti
J. S.
&
Rodell
M.
2014
A GRACE-based water storage deficit approach for hydrological drought characterization
.
Geophysical Research Letters
41
,
1537
1545
.
Wahr
J.
,
Molenaar
M.
&
Bryan
F.
1998
Time variability of the Earth's gravity field: hydrological and oceanic effects and their possible detection using GRACE
.
Journal of Geophysical Research: Solid Earth
103
,
30205
30229
.
Wang
B.
&
Fan
Z.
1999
Choice of South Asian summer monsoon indices
.
Bulletic of the American Meteorological Society
80
,
629
638
.
Xia
L.
,
Gao
Y.
,
Wang
W.
,
Lan
Y.
,
Xu
J.
&
Kai
L.
2014
Climate change and applicability of GLDAS in the headwater of the Yellow River Basin
.
Advances in Earth Science
29
,
531
540
.
Yang
D.
,
Yang
H.
&
Lei
H.
2014
Watershed Hydrology
.
Tsinghua University Press
,
Beijing
.
Yang
T.
,
Wang
C.
,
Chen
Y.
,
Chen
X.
&
Yu
Z.
2015
Climate change and water storage variability over an arid endorheic region
.
Journal of Hydrology
529
,
330
339
.
Yao
T.
,
Shilong
P.
,
Shen
M.
,
Jing
G.
,
Wei
Y.
,
Zhang
G.
,
Lei
Y.
,
Yang
G.
,
Zhu
L.
&
Xu
B.
2017
Chained impacts on modern environment of interaction between Westerlies and Indian Monsoon on Tibetan Plateau
.
Bulletin of Chinese Academy of Sciences
32
,
976
984
.
Zhang
L.
,
Wang
J.
,
Huang
Y.
,
Hao
W. U.
&
Duan
Q. C.
2015a
Characteristics of drought based on standardized precipitation evapotranspiration index from 1961 to 2010 in Yunnan province
.
Journal of Meteorology & Environment
31
,
141
146
.
Zhang
W.
,
Yi
Y.
,
Jafarov
E. E.
,
Yang
K.
,
Kimball
J. S.
&
Song
K.
2015b
Simulation of permafrost and seasonally frozen ground conditions and their response to recent climate warming in the Tibetan Plateau
. In:
Agu Fall Meeting
. .
Zhen
L.
,
Chuanyin
Z.
,
Baogui
K.
,
Yang
L.
,
Wanqiu
L.
&
Cai
Y.
2018
North China Plain water storage variation analysis based on GRACE and seasonal influence considering
.
Acta Geodaetica Et Cartographica Sinica
47
,
940
949
.

Supplementary data