Widespread decline in terrestrial water storage and its link to teleconnections across Asia and eastern Europe

Recent global changes in terrestrial water storage (TWS) and associated freshwater availability raise major concerns about the sustainability of global water resources. However, our knowledge regarding the long-term trends in TWS and its components is still not well documented. In this study, we characterize the spatiotemporal variations in TWS and its components over the Asian and eastern European regions from April 2002 to June 2017 based on Gravity Recovery and Climate Experiment (GRACE) satellite observations, land surface model simulations, and precipitation observations. The connections of TWS and global major teleconnections (TCs) are also discussed. The results indicate a widespread decline in TWS during 2002–2017, and five hotspots of TWS negative trends were identified with trends between −8.94 and−21.79 mm yr−1. TWS partitioning suggests that these negative trends are primarily attributed to the intensive over-extraction of groundwater and warmthinduced surface water loss, but the contributions of each hydrological component vary among hotspots. The results also indicate that the El Niño–Southern Oscillation, Arctic Oscillation and North Atlantic Oscillation are the three largest dominant factors controlling the variations in TWS through the covariability effect on climate variables. However, seasonal results suggest a divergent response of hydrological components to TCs among seasons and hotspots. Our findings provide insights into changes in TWS and its components over the Asian and eastern European regions, where there is a growing demand for food grains and water supplies.

However, limited works have focused on the contributions of TWS components at a large scale, particularly in water-limited and densely populated regions (Tapley et al., 2019), until two recent global scale studies significantly improved our knowledge by identifying 34 hotspots of TWS changes during 2002-2016 (Rodell et al., 2018) and the changes in global endorheic basin water storage (Wang et al., 2018). 40 The Asian and Eastern European regions, home to half of the world's population and 50% of the arid/semiarid climate areas, is undergoing intensive water threats to agriculture and domestic water needs (Huang et al., 2016) (Figure 1). Most of the countries located on the borders of Asia and Eastern Europe are experiencing water resource shortages owing to low annual precipitation (less than 400 mm/y), which consequently forms the largest amount of irrigated land in the world (Rodell et al., 2009). The increasing demand for freshwater and the limited knowledge on the available water resources in this region have 45 become the key challenges to achieving sustainability in these areas (Feng et al., 2013). Therefore, knowledge on the TWS trend and its hydrological components is important for sustainable development of water resource and food supply in this region (Rodell et al., 2018).
The large-scale mode of teleconnection (TC) is an overwhelming factor of regional water resources, modulating the location and strength of storm tracks and fluxes of heat, moisture and momentum. For example, prominent teleconnection patterns 50 such as El Niño-Southern Oscillation (ENSO) show that El Niño years are related to reduced precipitation, continental freshwater discharge, and evapotranspiration on land, and therefore, TWS variability occurs over many land areas (Gonsamo et al., 2016). Many studies have invested their efforts to address the possible causes of TWS changes by connecting TWS with TC (Ndehedehe et al., 2017;Ni et al., 2018;Phillips et al., 2012). However, these studies focused primarily on the effect of ENSO on TWS. Notably, many other global and regional climate TCs also influence the changes in TWS, which 55 has been less documented and consequently limits our understanding of a comprehensive TWS-TC correlation. Therefore, knowledge of the influence of multiple TCs on TWS is critical to improve our understanding and properly manage water resources (Ndehedehe et al., 2017;Phillips et al., 2012).
Here, we conduct a comprehensive analysis of the spatiotemporal variations in TWS across the Asian and Eastern European regions and address the contributions of each hydrological component and connection with TCs using multisource data. First, 60 we calculate the deseasonalized trend and analyse the spatiotemporal variations in TWS across Asia and Eastern Europe.
Then, we partition the components of TWS into SW, SM and groundwater by using GRACE, land surface model and lakes and glacial observation data. Finally, we calculate the cross-correlation coefficients between TCs and the remainder of the TWS time series. We aim to explore 1) the spatial pattern of long trends in TWS, 2) the contributions of water components https://doi.org/10.5194/hess-2019-281 Preprint. Discussion started: 15 August 2019 c Author(s) 2019. CC BY 4.0 License.
to TWS variations among regions, and 3) the role of TCs in the changes in TWS and its components over Asian and the 65 Eastern European regions.

Study area
The Asian and Eastern European regions are located between latitudes 6°S and 56°N and longitudes 4°E and 109°E. These regions comprise highly populated lands that vary spatially in climate, topography and hydrology. A total of 54% of the area 70 is arid and semiarid ( Figure 1a). The topography is dominated by plains, plateaus, and mountains such as the Himalaya Mountains, the Tien Shan Mountains, the Gangetic Plains, the Loess Plateau, and Mount Kilimanjaro. More importantly, the Asian and Eastern European regions are the most densely populated regions in the world (available at http://sedac.ciesin.columbia.edu/data/ collection/gpw-v4), sustaining nearly half of the total population of the world and contain the largest, most intensively irrigated land in this area (Figure 1b). Under the combined effect of climate change and 75 anthropogenic activities, the amount of sustainable available freshwater has become increasingly dire for food and water security and hence sustainable social economies.

Data
In this work, both of the mass concentration (mascon) solutions, including the Jet Propulsion Laboratory (hereafter JPL-M) (Watkins et al., 2015) and Center for Space Research (CSR-M) (Save et al., 2016) datasets were used to detect TWS changes 80 across the Asian and Eastern European regions with spatial resolutions of 0.5° × 0.5° for the period between April 2002 and June 2017 for a total of 183 solutions. The 20 months of missing data in the original time series were filled by the linear interpolation method (Long et al., 2015;Yang et al., 2017). The GRACE anomalies reported in these mascon solutions are relative to a 2004-2009 mean baseline (Scanlon et al., 2016). For details on the data processing, readers are referred to Save et al. (2016) and Watkins et al.(2015) for the mascon solutions. 85 The Global Land Data Assimilation System (GLDAS) was used to partition the GRACE-observed TWS changes into SW (snow water equivalent, canopy water, lakes and glaciers), SM and groundwater (Rodell et al., 2004). The monthly data products from the GLDAS version 2.1 Noah model contain thirty-six parameters, including canopy water storage, snow water equivalent and SM data. Noah has a total of 4 layers of SM thickness: 0-10, 10-40, 40-100, and 100-200 cm. To compute the Noah TWS, the SM in all layers, snow water equivalent, and canopy SW are summed. In this study, the 90 monthly canopy water storage, snow water equivalent and SM data were obtained from GLDAS version 2.1 with a spatial resolution of 0.25° × 0.25° during 2002 and 2017 (available at https://disc.gsfc.nasa.gov/datasets?page=1&project =GLDAS).
Then, we resampled these data to a 0.5° × 0.5° spatial resolution using the bicubic method prior to the analysis.
The integrated total water content, which is obtained from the GLDAS output by summing the layers of canopy water storage, four SM layers and snow equivalent water, is directly comparable to that which GRACE measures over land, and users need to be aware that the GLDAS version used here does not include groundwater and separate SW components (such as rivers and lakes). Thus, deviations from the GRACE total water storage changes can be expected. The comparison between GRACE and GLDAS is shown in Figure S1. Lake level data derived from multimission altimeter observations were obtained from the Database for Hydrological Time Series of Inland Water (Schwatke et al., 2015) and Hydroweb (Crétaux et al., 2011). Glacier mass change data are available from the published literature (Brun et al., 2017;Wang et al., 2018). The 100 groundwater counterpart was estimated by deducted the estimated SW and SM changes from the GRACE-observed TWS change (Wang et al., 2018). Notably, interannual variations in biomass have been shown to be well below the detection limits of GRACE (Rodell et al., 2005;Rodell et al., 2009); therefore, the variability in SW counterpart mainly involves changes in lakes, snow and ice in this study.

Time series decomposition
The deseasonalized time series was used to calculate the TWS trend by Sen linear regression. The original TWS signal is decomposed into long-term trends, seasonality signals, and residual components by implementing the Seasonal Decomposition of Time Series by Loess (STL) approach. The residuals (total-trend-seasonality) are used to calculate the 115 cross-correlation between the TWS and TC indices to partly reduce the interference with glacier melt and anthropogenic activities, i.e., TWS and TC indices, groundwater withdrawal for irrigation, mining, and human-made reservoirs for water resource management (Humphrey et al., 2016). The STL method is a robust method for detecting nonlinear time series in trend estimates, and the equation is as follows: where Stotal is the original signal, Slong-term is the long-term trend, Sseasonality is the seasonality signal and Residual is the residual component. We also provide an example for the TWS time series and its decomposition components in Figure S2. For detailed principles and applications of STL, readers are encouraged to refer to (Bergmann et al., 2012;Cleveland et al., 1990).

Theil-Sen trend analysis 125
We calculated the linear rates of TWS and precipitation change for the Asian and Eastern European regions from April 2002 to June 2017 using the Theil-Sen trend method (Sen, 1968) and tested the significance of the trend using the Mann-Kendall statistical test (Kendall, 1955). The advantage of the Theil-Sen trend is that it is nonsensitive to outliers and therefore can be more accurate than a simple linear regression for skewed and heteroscedastic data. This method competes well against the non-robust least squares method even for normally distributed data. The TWS trend β for a particular pixel is as follows: 130 where i,j are the time series sequences and xi,xj are the TWS values at times i,j. When β>0, the TWS in the corresponding pixel reveals an increasing trend; when β<0, the TWS in this pixel reveals a decreasing trend.

Cross-correlation analysis
Contemporaneous weather conditions impact TWS, often with evident time lags. Therefore, in this paper, we employ the 135 cross-correlation method to explore the relationship between TWS and teleconnection indices. The cross-correlation measures the similarity of the two time series datasets as a function of the displacement of one set relative to the other (Oppenheim et al., 2009). The cross-correlation is defined as follows: where σ_12 (τ) is the cross-covariance function of x_1 (t) leading the lagged x_2 (t), and σ_11 and σ_22 represent the auto-140 variances of x_1 (t) and x_2 (t), respectively. Moreover, we focus on the current and historical influence of TCs on TWS changes, and hence, we constrain that the value of ρ(τ) lies between 0 and 24 (Ni et al., 2018). Higher cross-correlation values indicate a stronger influence of a TC on TWS and its components. Spatially, we calculated the cross-correlation at the residual time series (full signal -(long-term trends + seasonality cycle)) to partially remove the influence of unaccounted factors on TWS changes, assuming that the trend was mainly induced by secular climate change or direct human impacts 145 rather than interannual climate variability (Phillips et al., 2012;Wang et al., 2018).  Figure S5 show the spatial distribution of the cross-correlation coefficients, illustrating the possible relationship of TCs with interannual variability in TWS. The results indicate that ENSO, AO and NAO have a comparable area 175 percentage (more than 20% of the study area) in influencing TWS variability. Spatially, the pattern of correlation coefficients between TWS and ENSO is heterogeneous, with positive correlations mostly in southeast China and boreal regions and negative correlations in southeast Asia, India, and eastern boreal regions ( Figure S5). The second and third dominant teleconnection modes are AO and NAO, respectively. AO mainly affects TWS variations across high-latitude regions through its impact on temperature variability, and NAO has a wider footprint that is scattered across the whole study.  Figure 2d.

Figure 2 and
Notably, the spatial pattern of the dominant TC having only a limited extent with respect to their influence on climate 190 conditions. The heterogeneous pattern highlighted the importance of focusing on the effect of multiple TCs on TWS rather than one teleconnection index.

Contributions of water storage components to TWS
TWS variations during 2002-2017 contributed differently to SW (lakes, biomass, snow, and ice), SM and groundwater ( Figure 3) among the five hotspot regions (circled in Figure 2a). Groundwater depletion dominance (-8.68±2.89 mm y -1 ) 195 compared to TWS loss is observed in region 1 compared to the other two components (-0.27±2.97 mm y -1 for SM, 0.02±0.11 mm y -1 for SW) during the study period. Similarly, a more prominent of groundwater depletion than TWS loss is also observed in northwest India (region 3 with -21.35 mm y -1 ), which experienced the most extensive land irrigation worldwide ( Figure 1b). The contributions of the SW and soil water in the above two regions are very small and neglected. The rapid glacial mass loss in the Tien Shan Mountain range induced a 41% water loss (Brun et al., 2017;Jacob et al., 2012). The melt 200 water could replenish the soil water, leading to an increase in soil water components (2.18±2.82 mm y -1 ) in northwest China (region 2) ( Figure 3). Notably, irrigated agriculture is the major style of agricultural development in this region (Yang et al., 2017), which contributes to more than half (-11.61±13.02 mm y -1 ) of the total water loss. Both groundwater depletion and glacial melt may enhance evapotranspiration by pumping water from aquifers and mountains to the surface, which induces a decline in TWS (Thomas et al., 2014;Wang et al., 2018). The contributions of the fourth region are 2.48±1.81 mm y -1 , -205 3.22±7.24 mm y -1 and -11.00±10.43 mm y -1 for SW, SM and groundwater, respectively, suggesting that groundwater withdrawal is the primary reason for TWS depletion. This region is also spatially coherent with precipitation deficits, as depicted in Figure 2b, particularly in the eastern part of this hotspot, accentuating the impact of meteorological drought on TWS. The prominent SW loss (-3.83±1.68 mm y -1 ) in the fifth region can be attributed to the sharply shrinkage in the Caspian Sea with a decrease in sea level of -73.2 mm y -1 , which is comparable to -68 mm y -1 during 2002-2016 (Wang et al., 210 2018). The dramatic decline in the Caspian Sea level contributes to a third of the total TWS loss in this region. The large decrease in TWS in these areas can also be attributed to the heavy reliance on groundwater for irrigation and domestic needs due to the construction of dams upstream (Voss et al., 2013;Rodell et al., 2018). The slight decreasing trend in precipitation and associated droughts will significantly exacerbate the TWS depletion, as the recharge of TWS relies largely on precipitation in those water-limited regions. 215

Divergent response of water storage components to TCs
Our results indicate that the water storage components are simultaneously influenced by several teleconnection.  (Table S2). For each component, TWS and groundwater are generally 220 less sensitive to TC signals compared to the surface and SM counterparts (Table S3). A possible explanation may be that TWS is a synthesis signal that its trend will be offset by its components in different ways, and the groundwater intensively interferes with anthropogenic activities such as irrigation by groundwater withdrawal, which indicates less correlation with TCs.
Further seasonal analysis indicates that the response of water storage to TCs is seasonally different from one region to 225 another. For example, TC signals have a dominant control on TWS and component variability in spring and summer for region 3 and region 1, respectively, whereas TC signals control more of the changes in SM in region 5 in autumn and winter.
Notably, although it has been well documented that the dramatic decline in TWS in northwest India can be attributed to the overexploitation of groundwater (Rodell et al., 2009), our seasonal response of water components to TCs suggests that the SM in this region is highly correlated with spring ENSO signals (Figure 4, r = 0.71). These results highlight the importance 230 of understanding the seasonal responses to TCs to improve the understanding and prediction of changes in TWS and associated water storage components.

Comparison of our results to previous studies
We investigated the spatiotemporal trend of TWS and its components over Asia and Eastern Europe during 2002-2017. The 235 overall pattern and trend of TWS over the study area are consistent with those of previous studies (Humphrey et al., 2016;Scanlon et al., 2016). In region 1, our estimate is within that of previous studies in this region (22±3 mm y-1 during 2003-2010) (Feng et al., 2013). Due to a long-term warm and dry climate and intensive anthropogenic activities (agriculture, industry, and urbanization), the groundwater in the region has been overexploited since the 1970s, and more than 70% of the groundwater exploitation is used for regional irrigation (Wang et al., 2007). The glaciers of Tien Shan Mountain are melting 240 rapidly (-0.28±0.20 m w.e.yr-1) (Jacob et al., 2012), accelerating an increase in the loss of water resources since the glacial meltwater will provide additional water that was lost to rivers or evaporation. Therefore, the negative trend in TWS indicates that water demand is larger than supply in region 2, which can be attributed to both continuous withdrawal of groundwater and extensive evaporation in the endorheic basin (Rodell et al., 2018). The increase in precipitation is expected to offset a certain portion of water losses in this region. In region 3, the rate of groundwater loss was also reported by a previous study 245 (approximately 40±10 mm y-1 from August 2002 to October 2008) (Rodell et al., 2009). As Indian agriculture leads the world in total irrigated land by consuming ~85% of the utilizable water resources (Panda et al., 2016;Salmon et al., 2015), a concluding consensus has been reached that the dramatic decline in TWS is mainly due to the overexploitation (extraction exceeding recharge) of groundwater for irrigation (Shamsudduha et al., 2019). Although precipitation in this region shows an increasing trend during the GRACE period, the rapid depletion of TWS in northwest India induced by unsustainable consumption of groundwater for irrigation and other anthropogenic uses has attracted worldwide attention because it is a major threat to India's sustainability (Panda et al., 2016;Rodell et al., 2009). Region 4 in our study is also heavily irrigated (Figure 1), so intensive irrigation is likely to induce groundwater decline. The increase in SW induced by melt water from mountains (Brun et al., 2017) was offset by the decrease in soil water that may be related to the decrease in precipitation ( Figure 2b). In addition to the water body loss, the overreliance on groundwater for domestic and agricultural needs was also 255 identified due to human-made dams influencing region 5 (Joodaki et al., 2014;Rodell et al., 2018;Voss et al., 2013), which is consistent with our results.
Notably, our study focused primarily on a water storage deficit hotspot analysis, which is free of bias compared to that obtained by the whole basin-based evaluation, where basins that experience both drying and wetting trends in the different sub-basins contribute to the calculation of the basin average TWS (Sun et al., 2018). Additionally, previous studies mainly 260 focused on the influence of ENSO on TWS (Ni et al., 2018;Phillips et al., 2012), and we provide a comprehensive analysis of 12 commonly used teleconnection indices and suggest that dominant TCs vary considerably for each component and region. Any conclusions or predictions of TWS and component variability based solely on TCs should be approached with caution, and multiple TCs are strongly recommended to explain the changes in TWS and its components.

Possible mechanisms of TC influence on TWS variability 265
Periodic variability in the climate system can strongly influence regional meteorological patterns and the associated TWS.
Unlike a single meteorological variable, teleconnection patterns control heat, moisture, and momentum balances through their effects on temperature, precipitation, and solar radiation reaching the Earth's surface. Therefore, the inherent mechanisms of the TC's influence on TWS variations are related to the combined simultaneous effects of TCs on regional climate factors (precipitation, temperature, and radiation); the changes in climate factors will significantly affect the recharge 270 (precipitation) and loss (evapotranspiration) of regional water resources, which eventually influences the changes in TWS.
This mechanism is also documented in former works (Ni et al., 2018;Zhu et al., 2017). We have identified several dominant TCs that influence the variability in TWS and its components. For example, during positive ENSO phases, warmer and drier conditions over India and southeast Asia are associated with a decrease in TWS. When the AO index is positive and the vortex is intense, the winds tighten like a noose around the North Pole, locking cold air in place and contributing to the 275 unusual warmth over the Northern Hemisphere land masses. The positive phase of the NAO, which is highly correlated with AO (r = 0.64, p < 0.01), leads to an intensification of the west wind drift due to a reinforcement of the Iceland low and the Azores high pressure systems, which particularly apparent in the boreal winter months (Wallace et al., 1981). In turn, this results in positive precipitation anomalies in central Europe and negative anomalies in southern Europe and on the Norwegian coast, which is reflected in the water storage variations. IOD is similar to ENSO and often co-occurring with 280 ENSO (Du et al., 2009). During the positive IOD phase, anomalous cool (warm) waters appear in the eastern (western) Indian Ocean in association with large-scale circulation changes that bring anomalous dry conditions to Indonesia, while East Africa experiences above-average rainfall (Webster et al., 1999). These results imply that climate variability may exert important influences on the TWS, which is generally consistent with the existing literature (Ni et al., 2018;Phillips et al., 2012;Zhang et al., 2015). However, challenges remain in separating the long-term relative roles of natural climatic variation 285 and anthropogenic forcings on TWS changes.

Implications for future hydrological studies
Our results indicate that climate variability could explain the increase in TWS and the decrease in most remote and sparsely populated regions. To a certain extent, climate variability could also indirectly explain glacial melt-induced changes in TWS, such as warming-induced glacial retreat. However, climate variability will interact with human activities, i.e., groundwater 290 abstraction, in regions with intensive human activities. Although we could obtain the contributions of different water storage components (SW, SM, and groundwater) through TWS partitioning, each component also influences both climate variability and human activities, which makes it extremely complicated to elicit the contribution of climate change and other processes.
Thus, well-designed experiments and coupled human-natural system models are still needed to clarify the quantitative contributions of each influencing factor on TWS and its component variability. However, our study provides a new view of 295 teleconnections to better understand the changes in TWS, enriching the limited knowledge pool of TWS dynamics attributions. Notably, we claim that climate variability-induced extremes, such as drought and heatwaves, will exacerbate the TWS loss through increased consumption of water resources from groundwater for irrigation and human water demand in these hotspots, rather than climate variability alone causing the observed TWS loss.
There are several recommendations for future hydrological studies. First, withdrawal of freshwater from groundwater in 300 water-limited regions is important for sustainable development of water resource and food supply (Rodell et al., 2018).
However, groundwater drought is a distinct class of drought resulting from a decrease in groundwater storage (Thomas et al., 2017). Understanding groundwater drought is important in water-limited regions where the interplay between groundwater recharge and abstraction results in variable groundwater stress conditions. GRACE has the unique potential to obtain data on groundwater storage by introducing subsidiary datasets. Second, glacier mass loss in mountainous areas can not only relieve 305 drought stress in drought years (Huss et al., 2018) but also results in hydroclimatic extremes, e.g., floods. Neither of these phenomena can be detected using only precipitation datasets, such as those commonly used in monitoring drought and flood events (Sherwood et al., 2014), which highlights the importance of TWS-related hydroclimatic extremes. With the release of the GRACE follow-up satellite , consecutive prolonged data records could provide valuable solution for evaluating hydrological conditions from a long-term perspective and would lead to considerable improvements in our 310 knowledge of TWS-related hydroclimatic extremes . Third, a recent study found that the CO2 growth rate is strongly sensitive to observed changes in TWS (Humphrey et al., 2018), and the coupling between the water and carbon cycles highlights the need for stronger interactions between the hydrological and biogeochemical research communities to achieve sustainable development of the Earth system.

Conclusions 315
In this study, we characterize the spatiotemporal variations in TWS and its components and connect these variations with TCs over the Asian and Eastern European regions from April 2002 to June 2017 using multiple sources of data. The results indicate a widespread decline in TWS during 2002-2017, and five hotspots of TWS negative trends were identified with trends between -8.94 mm y-1 and -21.79 mm y-1. Partitioning of TWS suggests that these negative trends are mainly attributed to intensive groundwater extraction and warming-induced SW loss, but the contributions of each hydrological 320 component vary among regions. The results also indicate that ENSO, AO and NAO are the largest three dominant factors in controlling variations in TWS through their covariability effects on climate variables. However, seasonal results suggest a divergent response of hydrological components to TCs among seasons and regions, highlighting the importance of knowledge of the seasonal responses to TCs to improve the understanding and prediction of changes in TWS and associated water storage components. Our study provides a comprehensive analysis of TWS variability and its connection to TCs across 325 the Asian and Eastern European regions, facilitating the target strategy of water resource management.

Data availability
The data and code generated in this study are available from the authors upon request (liuxianfeng7987@163.com).

Author contributions
XL and XF conceived and designed the research, XL conducted the experiments and analysed the results, XL wrote the 330 manuscript with contributions from XF, CP and BF.