# Terrestrial water storage changes over Xinjiang extracted by combining Gaussian filter and multichannel singular spectrum analysis from GRACE

Terrestrial water storage changes over Xinjiang extracted by combining Gaussian filter and... Summary Water resource management is crucial for the economic and social development of Xinjiang, an arid area located in the Northwest China. In this paper, the time variations of gravity recovery and climate experiment (GRACE)-derived monthly gravity field models from 2003 January to 2013 December are analysed to study the terrestrial water storage (TWS) changes in Xinjiang using the multichannel singular spectrum analysis (MSSA) with a Gaussian smoothing radius of 400 km. As an extended singular spectrum analysis (SSA), MSSA is more flexible to deal with multivariate time-series in terms of estimating periodic components and trend, reducing noise and identifying patterns of similar spatiotemporal behaviour thanks to the data-adaptive nature of the base functions. Combining MSSA and Gaussian filter can not only obviously remove the north–south striping errors in the GRACE solutions but also reduce the leakage errors, which can increase the signal-to-noise ratio by comparing with the traditional procedure, that is, empirical decorrelation method followed with the Gaussian filtering. The spatiotemporal characteristics of TWS changes in Xinjiang were validated against the Global Land Dynamics Assimilation System, the Climate Prediction Center and in-situ precipitation data. The water storage in Xinjiang shows the relatively large fluctuation from 2003 January to 2013 December, with a drop from 2006 January to 2008 December due to the drought event and an obvious rise from 2009 January to 2010 December because of the high precipitation. Spatially, the TWS has been increasing in the south Xinjiang, but decreasing in the north Xinjiang. The minimum rate of water storage change is −4.4 mm yr−1 occurring in the central Tianshan Mountain. Hydrogeophysics, Satellite gravity, Time variable gravity, Inverse theory, Spatial analysis 1 INTRODUCTION Xinjiang is located in one inland arid area in the Northwest China (35°25΄N-49°10΄N and 73°40΄E- 96°23΄E, see Fig. 1). The scarce water resource severely restricts the local economic and social development. For example, agriculture, as the most important part of Xinjiang economy, highly relies on the water resource management. Xinjiang is divided into the southern and the northern parts by Tianshan Mountain. There are more than 570 rivers in its territory. The total water resources are about 882 × 108 m3 including the surface water resources of about 794 × 108 m3 and the ground water resources of about 88 × 108 m3. The water production per unit land area is only 5.30 × 104 m3 km−2 which accounts for only 18 per cent of the average water production in China (Zhang 2001). Study on terrestrial water storage (TWS) changes can help to scientifically optimize the allocation of water resources, protect the ecological environment and promote the sustainable development of water resources in Xinjiang. Figure 1. View largeDownload slide Geographical distribution of Xinjiang and meteorological stations. Figure 1. View largeDownload slide Geographical distribution of Xinjiang and meteorological stations. The gravity recovery and climate experiment (GRACE) satellite was successfully launched in 2002 March. It is jointly sponsored by the National Aeronautics and Space Administration (NASA) of USA and the Deutsches Zentrum für Luft- und Raumfahrt. GRACE can probe the temporal variations of Earth gravity field from which the Earth surface mass balance can be extracted (Wang et al. 2015). The Stokes coefficients estimated from GRACE mission are widely used in geophysical, glaciological, oceanographic and hydrological studies due to its unprecedented accuracy. For example, GRACE data are used extensively to estimate the TWS changes and groundwater variations over the large-scale areas (Luo et al. 2012; Seoane et al. 2013; Lu et al. 2015; Guo et al. 2016; Hassan & Jin 2016; Deng & Chen 2017), monitor local drought and flood events (Houborg et al. 2012; Li et al. 2013; Cao et al. 2015) and determine mass changes of the Antarctic ice sheet (Chen et al. 2009; Ju et al.2013, 2014; Gao et al. 2015; Mu et al. 2017). Many other applications such as the sea level change (Tamisiea 2011) and the glacial isostatic adjustment (GIA, Steffen et al. 2008) are also widely studied using GRACE data. GRACE-derived high-degree Stokes coefficients are contaminated by random noises which make the GRACE TWS changes contain the obvious north–south striping errors in the spatial domain and the high-frequency errors in the frequency domain. These errors are usually be removed or reduced by the Gaussian filtering technique (Wahr et al. 1998), the decorrelation procedure (Swenson & Wahr 2006), the DDK filter (Kusche 2007), the Wiener filter (Klees et al. 2008), the principal component analysis (PCA, Rangelova et al. 2007; Mu et al. 2014), the singular spectrum analysis (SSA, Wang et al. 2011) or the independent component analysis (Guo et al. 2014), etc. These methods have achieved good results, but there are still residual errors in the water storage series. In this paper, we combine the Gaussian filter and multichannel singular spectrum analysis (MSSA) to denoise the north–south striping errors in GRACE solutions. MSSA as one data-adaptive method without the prior information can divide periodic components, separate trend, filter out noise and distinguish patterns of similar spatiotemporal behaviours (Zotov 2012). We apply MSSA to the GRACE Stokes coefficients and then implement the Gaussian filtering technique to better separate noises and real signals. The effectiveness of this procedure is demonstrated by comparing to the results of traditional method, that is, the empirical decorrelation method (EDM, Duan et al. 2009) combined with the Gaussian filter. Using this alternative method, we extract TWS changes in Xinjiang which are then verified against the hydrological models like the Global Land Dynamics Assimilation System (GLDAS), the Climate Prediction Center (CPC) and in-situ precipitation data, and then the temporal and spatial changes of TWS in Xinjiang are analysed. 2 DATA AND METHODOLOGY 2.1 GRACE products and water storage estimates Stokes coefficients determined from GRACE data released by the Center for Space Research, the Texas University at Austin, were used in the study. The monthly GRACE level-2 products are downloaded from http://isdc.gfz-potsdam.de/grace through January 2003 to December 2013 with the exceptions of 2003 June, 2011 January and June, 2012 May and October, and 2013 March, August and September due to lack of data. These losing data are interpolated from the neighbouring monthly solutions with the cubic spline interpolation method. The GRACE Level-2 RL05 data are truncated at maximum degree and order 60 without the effects of atmosphere, ocean and tide (Bettadpur 2012). The degree-1 coefficients are not restored as their impact on this regional study is rather limited. Because the C20 coefficient measured by the satellite laser ranging (SLR) is much better than that in RL05, C20 in RL05 was replaced with those derived from the SLR data (Cheng & Tapley 2004). The terrestrial mass balance inversed from GRACE can effectively indicate TWS change with the equivalent water height (EWH, Chao et al. 2015). The EWH determined from GRACE products is   \begin{eqnarray} \Delta h(\theta ,\phi ) &=& \frac{{a{\rho _{{\rm{ave}}}}}}{{3{\rho _w}}}\sum\limits_{l = 0}^\infty \sum\limits_{m = 0}^l {{\overline p }_{lm}}(\cos \theta )\frac{{2l + 1}}{{1 + {k_l}}}\nonumber\\ &&\times \,(\Delta {C_{lm}} \cos m\phi+\Delta {S_{lm}}\sin m\phi) \end{eqnarray} (1)where a is the averaged Earth radius, ρave is the mean Earth density (5517 kg m−3), ρw = 1000 kg m−3, θ and ϕ are the colatitude and longitude, respectively, $${\overline p _{lm}}(\cos \theta )$$ is the fully normalized associated Legendre function with degree l and order m, kl is the load Love number of degree l and ΔClm and ΔSlm are the fully normalized Stokes coefficient changes determined by GRACE with respect to the mean gravity from 2003 January to 2013 December. Combining the Gaussian filter and MSSA (GMSSA) can weaken and remove the effects of satellite orbit errors and Stokes coefficient errors in GRACE products. The Gaussian smoothing kernel function (Jekeli 1981; Wahr et al. 1998) is introduced to eq. (1) into get the filtered EWH as   \begin{eqnarray} \Delta h(\theta ,\phi ) &=& \frac{{a{\rho _{{\rm{ave}}}}}}{{3{\rho _w}}}\sum\limits_{l = 0}^\infty \sum\limits_{m = 0}^l {{\overline p }_{lm}}(\cos \theta )\frac{{2l + 1}}{{1 + {k_l}}}{w_l}\nonumber\\ &&\times \,(\Delta {{\overline C }_{lm}} \cos m\phi +\Delta {\overline S _{lm}}\sin m\phi ) \end{eqnarray} (2)where wl is the Gaussian smoothing kernel function for degree l. $$\Delta {\overline C _{lm}}$$ and $$\Delta {\overline S _{lm}}$$ are the spherical harmonic coefficient changes reconstructed by MSSA. The TWS changes in Xinjiang are calculated from eq. (2) on grid of 1° × 1°, and the results are further improved using the ICE-5G model (Peltier 2004) to remove the effect of GIA which is small in the study region. The water storage changes over the whole Xinjiang are then estimated. 2.2 Hydrological models Two hydrological models, GLDAS and CPC, are used in the study. GLDAS (Rodell et al. 2004) is based on the Noah terrestrial ground model with the time resolution of 1 month and the spatial resolution of 0.25° × 0.25° from 2003 to 2013. The model is released by the Goddard Space Flight Center of NASA and the National Center of Environmental Prediction, USA. The water storage change can be computed from the soil moisture variations up to four layers (0–0.1, 0.1–0.4, 0.4–1 and 1–2 m), snow water change and canopy water change with the hydrological model. CPC model (Fan & van den Dool 2004) is released by the CPC of National Oceanic and Atmospheric Administration, USA. The time and spatial resolutions of the model are one month and 0.5° × 0.5°, respectively. The model can provide the soil moisture up to depth of 0–1.6 m caused by water balance changes. GLDAS and CPC data are all converted to spherical harmonic coefficients up to degree and order 60 and then filtered using the same strategy applied to the GRACE data. 2.3 Precipitation The monthly precipitations collected by 54 weather stations (see Fig. 1) located in Xinjiang from 2003 to 2013 are used in the study. The Kriging method is used to interpolate the precipitation data with grid of 1° × 1° over Xinjiang. The latitude cosine weighting method is used to calculate the mean precipitation for the whole Xinjiang, which can be compared with the GRACE-derived results. 2.4 MSSA method MSSA as an extended SSA can be used to extract the detailed temporal and spectral information from the multidimensional time-series (Plaut & Vautard 1994; Shen et al. 2018). It is also often assumed to be synonymous with the extended Empirical Orthogonal Function (EOF). Instead of using the sample correlation matrix in EOF/PCA, MSSA forms the trajectory matrix with a time lag to be analysed, and different channels can help each other to capture spatiotemporal correlation patterns. MSSA is more flexible to recognize trend, modulate oscillations of different periods and denoise multidimensional time-series (Zotov & Shum 2010; Oropeza & Sacchi 2011). It will be useful for filtering out strips and high-frequency noise in the GRACE data (Rangelova et al. 2010, 2012). 3 MSSA ON GRACE SPHERICAL HARMONICS There is one data set x (xn,l:n=1,N;l=1,L), which consists of N observations (time epochs) each with L variables (channels). It can be processed by MSSA (Oropeza & Sacchi 2011) with the following procedure. (1) The data set x is the time-series of Stokes coefficient {ΔCij(tn), ΔSij(tn)}, n = 1, …, N, and we select the lag parameter M (a lag window) which satisfies $$1 \le M \le \frac{N}{2}$$, where N = 132 and M = 60 in the GRACE data processing. So the trajectory matrix of ΔCij channel is   $${{\bf X}_{\Delta {C_{ij}}}} = \left( {\begin{array}{@{}*{4}{c}@{}} {\Delta {C_{ij}}({t_1})}&{\Delta {C_{ij}}({t_2})}& \quad\cdots &\quad{\Delta {C_{ij}}({t_K})}\\ {\Delta {C_{ij}}({t_2})}&{\Delta {C_{ij}}({t_3})}& \quad\cdots &\quad{\Delta {C_{ij}}({t_{K + 1}})}\\ \vdots & \vdots & \quad\cdots & \vdots \\ {\Delta {C_{ij}}({t_M})}&\quad{\Delta {C_{ij}}({t_{M + 1}})}& \quad\cdots &\quad{\Delta {C_{ij}}({t_N})} \end{array}} \right)$$ (3)where K = N − M + 1. For GRACE data, the trajectory matrix X composed of such blocks for every channel ΔCij and ΔSij is   $${\bf X} = {\left[ {{{\bf X}_{\Delta {C_{20}}}},{{\bf X}_{\Delta {S_{20}}}}, \ldots ,{{\bf X}_{\Delta {C_{ij}}}},{{\bf X}_{\Delta {S_{ij}}}}, \ldots {{\bf X}_{\Delta {C_{60}}}},{{\bf X}_{\Delta {S_{60}}}}} \right]^{\rm{T}}}$$ (4)in which X is one LM × K matrix. (2) Matrix X can be decomposed by the singular value decomposition, that is, X = USVT, in which S is one diagonal matrix with D diagonal elements. Here, D = min (LM, K). These elements are singular values which are sorted in descending order (see Fig. 2). Figure 2. View largeDownload slide Singular values for MSSA on GRACE products. Figure 2. View largeDownload slide Singular values for MSSA on GRACE products. The big singular value stands for the part with large power spectrum and the small corresponds to noises. The trajectory matrix X can be synthesized by the elementary matrices as   $${\bf X} = {{\bf T}_{\rm{1}}} + {{\bf T}_{\rm{2}}} + \cdots + {{\bf T}_{\rm{D}}}$$ (5)where the ith component corresponds to the matrix $${{\bf T}_{i}} = {s_{i}}{{\bf u}_{i}}{\bf v}_{i}^{T}$$ in which ui and vi are the left singular vector and the right singular vector for the singular value si of matrix X. (3) Since GRACE data time-series includes different sorts of signals such as long-term, annual and semi-annual terms and noises. Ti representing the same signal can be classified into one group using the ω-correlation method (Hassani 2007). The ω-correlation is   $$\rho _{i,j}^w = \frac{{({Y^{(i)}},{Y^{(j)}})}}{{{{\left\| {{Y^i}} \right\|}_w}{{\left\| {{Y^j}} \right\|}_w}}}$$ (6)where Yi is the time-series reconstructed from seriesTi using eq. (8), $${\| {{Y^i}} \|_w} = \sqrt {({Y^{(i)}},{Y^{(i)}})}$$ and $$({Y^{(i)}},{Y^{(j)}}) = \sum_{k = 1}^N {{w_k}y_k^iy_k^j}$$ in which wk = min (k, M, N − k). The large absolute $$\rho _{i,j}^w$$ indicates Ti and Tj may be correspond to the same signal so that these two components should be grouped together. Fig. 3 shows the absolute $$\rho _{i,j}^w$$ for the first 30 Ti in processing GRACE products in the manuscript. From Fig. 3, we can find when i, j > 12 the absolute $$\rho _{i,j}^w$$ are all relatively high, which means these series cannot be well separated and maybe contain noises. So these first 12 singular values and their corresponding singular vectors are used to reconstruct the interesting signals. The first two components correspond to the annual cycle, the eighth and ninth components are composed to denote the semi-annual signal, and the rest of the first 12 modes represent the long-term single, other useful signals and some noises. This process is equivalent to denoising and filtering GRACE products. The trajectory matrix reconstructed from the first 12 components is   \begin{eqnarray} {{\bf X}^{\boldsymbol{'}}} &=& {{\bf T}_1} + {{\bf T}_2} + \cdots + {{\bf T}_{12}}{\rm} \nonumber\\ &=& {\left[ {{{{\bf X'}}_{\Delta {C_{20}}}},{{{\bf X'}}_{\Delta {S_{20}}}}, \ldots ,{{{\bf X'}}_{\Delta {C_{ij}}}},{{{\bf X'}}_{\Delta {S_{ij}}}}, \ldots {{{\bf X'}}_{\Delta {C_{60}}}},{{{\bf X'}}_{\Delta {S_{60}}}}} \right]^{\rm{T}}} \end{eqnarray} (7) Figure 3. View largeDownload slide ω-correlation for the first 30 reconstructed components. Figure 3. View largeDownload slide ω-correlation for the first 30 reconstructed components. (4) Denoised data can be reconstructed from the trajectory matrix X΄ (Golyandina et al. 2001). Let $${\bf Z} = {\bf X}_i^\prime$$(i = ΔC20, ΔS20, …ΔS60), then the reconstructed series dk(ΔCij,ΔSij) is   $${d_k} = \left\{ \begin{array}{@{}l@{}} \frac{1}{k}\sum\limits_{n = 1}^k {{z_{n,k - n + 1}}} \\ \frac{1}{{{M^*}}}\sum\limits_{n = 1}^{{M^*}} {{z_{n,k - n + 1}}} \\ \frac{1}{{N - k + 1}}\sum\limits_{n = k - {K^*} + 1}^{N - {K^*} + 1} {{z_{n,k - n + 1}}} \end{array} \right.\quad \begin{array}{@{}*{1}{l}@{}} {1 \le k < {M^*}}\\ {}\\ {{M^*} \le k \le {K^*}}\\ {}\\ {{K^*} < k \le N}\\ {} \end{array}$$ (8)where M* = min (M, K), K* = max (M, K). 4 RESULTS AND ANALYSIS 4.1 Method validation Fig. 4 illustrates the effect of different methods on the GRACE monthly EWH anomaly in 2007 April. The top subfigures show the original fields with different Gaussian smoothing radii. Without smoothing the map is dominated by the prominent north–south stripes and little geophysical signal can be observed. Using the smoothing radius of 300 km, hydrological signals over large basins (e.g. Amazon Valley, Congo River Basin and Greenland) start to stand out, but are still corrupted with the stripes. When the radius increases to 500 km, most of the stripes over land are removed and the oceans appear mostly stripe free when using a smoothing radius of 700 km. While increasing the smoothing radius, the magnitude of signals (e.g. Amazon Valley) is reduced and the leakage errors (e.g. Antarctica) increase. The second-row subfigures display the fields after applying MSSA combining with different Gaussian filtering radius. The third-row subfigures display the fields after applying EDM combining with different Gaussian filtering radius (GEDM for short) to the Stokes coefficients. The bottom subfigures show the differences between the figures in the second and third rows. Without the Gaussian filtering, MSSA and EDM both can significantly reduce the stripes, and the geophysical signal over large basins can be easily observed, although some noisy patterns remain exist. For MSSA filtering only, the stripes mainly distributed in low/mid-latitudes and slight in high latitudes. However, in the EDM result the stripes are mainly concentrated in low/mid-latitudes, and the decorrelation affects geophysical signals, especially at high latitudes (e.g. Antarctica and Greenland) where the magnitudes of signals are evidently smaller than that in the MSSA result. Moreover, with the same Gaussian smoothing radius (e.g. 300 km), GMSSA can preserve significantly more signals and show less pronounced leakage effects. Examples of areas where this effect is clearly visible are Antarctica, the Amazon River Basin, the Yukon and Fraser Basins (located in Alaska and western Canada), Madagaskar et al. For instance, due to less leakage effects, the GMASS results better coincide with coastlines in the western Canada than those in the GEDM. Figure 4. View largeDownload slide Global EWHs in 2007 April with various Gaussian smoothing radii. (a) 0 km, (b) 300 km, (c) 500 km and (d) 700 km. The upper figures show the unfiltered fields, the middle figures represent the MSSA and EDM method respectively, and the bottom figures display the difference between MSSA and EDM. Units are in mm. Figure 4. View largeDownload slide Global EWHs in 2007 April with various Gaussian smoothing radii. (a) 0 km, (b) 300 km, (c) 500 km and (d) 700 km. The upper figures show the unfiltered fields, the middle figures represent the MSSA and EDM method respectively, and the bottom figures display the difference between MSSA and EDM. Units are in mm. We also use GMASS and GEDM with 400 km smoothing radius to extract the global EWHs in 2007 April and October because the land EWHs of these two months are violent and converse over low/mid-latitude. Fig. 5 shows the inversed results. The upper figures represent the results of GEDM and GMSSA in 2007 April, respectively, and the bottom figures are the results of 2007 October correspondingly processed by GEDM and GMSSA. From Fig. 5, we can find that these two processing strategies give consistent results without obvious north–south striping noises. Table 1 lists the statistical results. The root mean squares (rms) of EWH estimated by these two strategies are almost identical. The maximum EWH from GEDM is only 286.3 mm which is less than that from GMSSA, but the minimum values are basically identical in 2007 April. The maximum EWH from GEDM is only 216.3 mm which is greater than that from GMSSA, and the minimum value is less than that determined by GMSSA in 2007 October. Moreover, the mean EWH from GMSSA is larger than that from GEDM in 2007 April, but it appears small comparing to the GEDM in 2007 October. The idea of EDM is to fit a quadratic polynomial in a moving window to the Stokes coefficients of even and odd degrees for a particular order and remove this from the original Stokes coefficients. MSSA uses the Stokes coefficients to construct a trajectory matrix which is decomposed and reconstructed so that the signals representing significant geophysical can be extracted. The different processing strategies lead to different results listed in Table 1. Figs 4 and 5, and Table 1 can validate the feasibility and effectiveness of GMSSA by comparing with GEDM. Figure 5. View largeDownload slide Global equivalent water storage changes in 2007 April and October processed by (a) GEDM, and (b) GMSSA. Figure 5. View largeDownload slide Global equivalent water storage changes in 2007 April and October processed by (a) GEDM, and (b) GMSSA. Table 1. Statistics of water storage changes from GEDM and GMSSA with filtering radius of 400 km (unit in mm). Month  Method  Max  Min  Mean  rms  2007 April  GEDM  286.3  −205.4  2.6  44.1    GMSSA  327.4  −204.4  4.3  44.3  2007 October  GEDM  216.3  −425.1  1.0  43.2    GMSSA  196.0  −409.4  0.0  40.9  Month  Method  Max  Min  Mean  rms  2007 April  GEDM  286.3  −205.4  2.6  44.1    GMSSA  327.4  −204.4  4.3  44.3  2007 October  GEDM  216.3  −425.1  1.0  43.2    GMSSA  196.0  −409.4  0.0  40.9  View Large As another tool to analyse the effectiveness of GMSSA, we use the method suggest by Chen et al. (2006). The idea is based on the fact that GRACE measurement errors are approximately at the same level over both land and ocean, but the surface mass variability is stronger in the continents than that in the oceans. The optimal filter should maximize the quotient of the latitude weighted rms of the continental and oceanic signals which represent the signal-to-noise ratio. Additionally, in order to reduce the leakage of the continental signal only the ocean points farther than 400 km from the coast are included. As a result, the rms ratio from GEDM is 1.45 which is lower than that from GMSSA, that is, 1.72. So GMSSA is able to obtain greater rms ratio and retain more interesting geophysical signals which can improve the precision and reliability of inversed water storage. It is inevitable that some useful geophysical signals will be filtered out together with noises in the filtering procedure. We will evaluate to which extent the filter might remove true geophysical signal. Here, we address this problem based on the study of rms of the linear trend of surface mass changes from GRACE and noise-free synthetic data. The trend is determined by least-squares fitting the GRACE and synthetic monthly Stokes coefficients during 2003 January to 2013 December with an offset, a linear term that represents the trend, an annual term and a half-yearly term as Rangelova et al. (2012). In addition, the synthetic data consist of global grids with ocean bottom pressure over the oceans and with GLDAS over land for the same period as GRACE, and a surface mass change rate corresponding to the gravity variation of a GIA model from Paulson et al. (2007) have been added back to all grids. As GLDAS simulated values over ice covered regions (e.g. Greenland and Antarctica) are not realistic, these regions are avoided in the comparison with GRACE. The synthetic data are converted to spherical harmonic coefficients up to degree 60, and the degree-1 coefficients are set to zero. Due to the influence of the noises in GRACE data, the rms of trend is higher in the unfiltered data but lower after using a filtering technique. However, an excessive reduction of rms would mean an undesirable attenuation of the geophysical signals. The optimal filter should minimize the rms of the global linear trend and maximize it in the noise-free synthetic data. Besides, the filters attenuate the signal in both data, but noises reduce only in real GRACE data. As the GLDAS data used here does not include groundwater and separate surface water components (such as rivers and lakes) which make the linear trend rms from synthetic data lower than those from GRACE data. The optimal filter also should show the difference rms (Drms) of linear trend between GRACE and synthetic data minimum (Belda et al. 2015). We present in Fig. 6 the trend result computed by no-filtering, GEDM and GMSSA, respectively, The GEDM and GMSSA both use Gaussian radii 400 km. As the result in Fig. 6(a) is computed without filtering, it shows obvious north–south stripes noises and a little more geophysical signals over large basins. Comparing to the GEDM results that present in Fig. 6(b), the GMASS results (Fig. 6c) are better coincide with coastlines which can be easily found in the Greenland, Antarctica, Scandinavia, etc. which indicates that GMSSA can more effectively reduce the signal leakage errors than GEDM. Over Antarctic, there are some differences between the results of GMSSA and GEDM, for example, we can see the signals of the apparent highs and lows at the Amundsen Sea Embayment in the unfiltering and GMSSA results, but it does not appear in the GEDM results. Table 2 lists rms of the trend of global mass changes except Greenland and Antarctic applying GEDM and GMSSA filtering, respectively. We see that different approaches do not lead to significantly different rms for the GRACE data. When applying the GEDM and GMSSA to the noise-free synthetic data, the rms value of later is larger than that of former, which means that the GMSSA can preserve more signals. Furthermore, the Drms value between the result of GRACE and the synthetic result from GMSSA is 6.46, which is slightly smaller than that from GEDM (6.79). The results of Fig. 6 and Table 2 verify the reliability of the new method. Figure 6. View largeDownload slide GRACE global EWH rate during 2003–2013 determined from (a) unfilter, (b) GEDM and (c) GMSSA. Figure 6. View largeDownload slide GRACE global EWH rate during 2003–2013 determined from (a) unfilter, (b) GEDM and (c) GMSSA. Table 2. Rms of the trend of global mass changes except Greenland and Antarctic (unit in mm). Data  Unfilter  GEDM  GMSSA  GRACE  15.22  8.84  8.83  Synthetic  4.15  2.05  2.37  Data  Unfilter  GEDM  GMSSA  GRACE  15.22  8.84  8.83  Synthetic  4.15  2.05  2.37  View Large 4.2 Analysis on TWS changes in Xinjiang The mean water storage change series over Xinjiang from GRACE, GLDAS and CPC are shown in Fig. 7. We can find that these three series have the obvious seasonal variations with the maximum in summer (June–August) and the minimum in winter (December–February). GRACE-derived water storage changes have the great anomalies in 2008 and 2009 which are significantly smaller than other years, and the anomalies can also be found in GLDAS and CPC results. This phenomenon is mainly caused by the continuous high temperature and less precipitation from 2008 April to 2009 December in Xinjiang (Cao et al. 2015). After a continuous drought event taking place from summer to winter in 2008, Xinjiang again experienced severe drought in 2009, which led to the minimum annual precipitation since 1961. However, from 2009 December to 2010 February, the strongest snowstorm in the last 60 yr caused the whole 2010 year's annual precipitation exceeding the historical extremes (Chen et al. 2011). Those results in the water storage rapidly increase in the study area that can be obviously found in the result of GRACE, GLDAS and CPC in Figs 7 and 8. This indicates that GRACE has a good prospect in detecting drought and flood events. We also find that there was more obviously amplitude anomaly in 2008 November determined from GRACE than that from GLDAS and CPC, which can further indicate that GRACE has the ability to monitor the regional water storage changes. Figure 7. View largeDownload slide TWS changes in Xinjiang from GRACE, GLDAS and CPC. Figure 7. View largeDownload slide TWS changes in Xinjiang from GRACE, GLDAS and CPC. Figure 8. View largeDownload slide Interannual variations of water storage after removing the annual and semi-annual signals. Figure 8. View largeDownload slide Interannual variations of water storage after removing the annual and semi-annual signals. Over the course of the research time, the water storage change amplitude in Xinjiang, estimated by GRACE data, drifted from −43.5 to 31.2 mm. The maximum was present in 2011 May and the minimum took place in 2008 November. In particular, the amplitude derived from GRACE is obviously greater than those from GLDAS and CPC. These differences among the estimates can either be attributed to model deficiencies, such as missing surface and groundwater components in GLDAS and inadequate snow water equivalent in CPC, or due to uncertainties related to the GRACE data, like aliasing, data processing, leakage error, etc. Moreover, comparing GRACE results with hydrological models shows a phase lag of around 1–2 month which is mainly caused by the soil water interception and the regulation and storage capacity of lake, glacier and groundwater. The correlation coefficient between GRACE and GLDAS is 0.86 which is greater than that between GRACE and CPC of 0.72 within a 95 per cent confidence interval. This may be attributed to variances in the model structure, the forcing data choices and the components of the TWS changes. Although the hydrological models are not a perfect reproduction of TWS changes, it is sufficiently enough to capture the magnitude and variability of terrestrial hydrology in the study area. The strong correlations between GRACE and these two models can verify the reliability of GRACE-inversed results. The annual (first and second components) and semi-annual (eighth and ninth components) signals are removed in the GRACE MSSA filtering process in order to study the interannual variation of water storage in Xinjiang (Fig. 8). From Figs 7 and 8, we can find that the water storage in Xinjiang has the relatively large fluctuation form 2003 to 2013. To better illustrate the TWS changes in Xinjiang, we divided the research time span into four short intervals according to our visual detection of the variations in the TWS. The first interval spans from 2003 January to 2005 December during which the water storage presented the rising trend. The second interval covers the period from 2006 January to 2008 December during which the water storage dropped quickly. The third interval lasts from 2009 January to 2010 December during which the water storage rose obviously. The last interval starts from 2011 January to 2013 December during which the water storage first decreased and then swelled. GRACE-derived interannual change agreed relatively well with those from hydrological models, although there are small deviations in some time quantum, such as the interval spans from 2003 January to 2005 December. However, The GRACE-derived results are much notable, for example, the drought event in 2008 and the high precipitation in 2010 are all obviously detected by GRACE. The water storage changes determined from GRACE excluding the annual and semi-annual signals are fitted with the least-squares method to estimate the water storage rates in Xinjiang from 2003 to 2013 whose spatial distribution is shown in Fig. 9. The results show that the water storage overall appears the rising trend in the south Xinjiang with the maximum rising rate of 9.5 mm yr−1 and the descending trend in the north Xinjiang with the minimum rate of −4.4 mm yr−1. The water storage shows the gradual increase from the north to the south in the south Xinjiang. Figure 9. View largeDownload slide TWS change rates in Xinjiang from GRACE. Figure 9. View largeDownload slide TWS change rates in Xinjiang from GRACE. Contrary to the south Xinjiang, the north Xinjiang wholly has the descending trend of water storage. The descending rates of water storage are bigger close to Tianshan. This may be certainly relative to the great groundwater exploitation for the agriculture and living water demand over the piedmont region. But the water storage has the rising trend on the north Altai area. The water storage wholly appears the descending trend over the Tianshan except slightly rising in the east area. The central Tianshan has the minimum rate of −4.4 mm yr−1. Tianshan breeds the abundant water resources including glaciers and permanent snows which seriously affect the regional water storage change. So the descending water storage in Tianshan can indirectly respond to the accelerated melting of glaciers and snows (Li et al. 2011). 4.3 Comparison of GRACE-derived water storage and actual precipitation The regional water mass balance can be depicted by   $$\frac{{dW(t)}}{{dt}} = {\rm{precipitation}}\left( t \right){\rm{ }} - {\rm{ evaporation}}\left( t \right){\rm{ }} - {\rm{ runoff}}\left( t \right)$$ (9)where dW(t)/dt is the change in water content. According to eq. (10), the precipitation is an important factor causing the TWS change. However, there is no one-to-one relation between water storage and precipitation, as the evaporation and runoff can also lead to changes. In the arid land of Xinjiang, surface and groundwater resources are scarce, and precipitation is the primary water source. Thus analysing precipitation and water storage information on their overall relationship is important to interpret the GRACE results. Fig. 10 shows the actual precipitation and the GRACE-derived water storage change in Xinjiang. We can find that the water storage and the precipitation all have the obvious seasonal variations with the identical change procedure. The precipitations mainly concentrated from May to September during which the GRACE-derived water storage was high, and there are less precipitations from October to February of the following year during which low water storage was detected by GRACE. In 2008 and 2009, the precipitations were significantly less than the same period. Correspondingly, the GRACE-derived water storages were also very less at the same period. In particular, we can see that precipitation is larger in 2007 than that in 2006. But from 2006, the GRACE-derived water storage is descending. As the evaporation and runoff both play an important role in regional water allocations. Snow/ice melt water were also a crucial factor to supply the local water storage besides precipitation in Xinjiang. The GRACE results reflect the TWS changes that caused by all factors. Despite less precipitation, the joint effects of all factors result in a higher TWS in 2006 than that of 2007. Over the course of the research time, the GRACE TWS changes showed a relatively good agreement with precipitation, which could be explained that the precipitation is a main factor to control the water storage in Xinjiang. Figure 10. View largeDownload slide TWS and precipitation changes in Xinjiang from 2003 to 2013. Figure 10. View largeDownload slide TWS and precipitation changes in Xinjiang from 2003 to 2013. 5 CONCLUSIONS In this paper, we combine MSSA with Gaussian filter to process the GRACE products to extract the water storage signals over Xinjiang, a northwest arid area of China, from 2003 to 2013. GMSSA can efficiently denoise the high-degree Stokes coefficients and decorrelate GRACE data which reduce the north–south striping in the GRACE EWH fields. In addition, this method can also reduce the leakage errors thanks to the enhanced signal-to-noise ratio by comparing with GEDM, so it can retain more useful geophysical signals from GRACE products. However, the method can also inevitably remove some interesting signals. How to reasonably determine the window M and the number of singular values in MSSA and optimally select the filtering window in the Gaussian filter should be further studied. TWS change has obvious seasonal and periodic variations with the maximum in summer and the minimum in winter. The GRACE-derived water storage is basically consistent with those from these two hydrological models and actual precipitation, but the amplitudes of water storage change from GRACE are greater than those from GLDAS and CPC, which can indicate that the glacier, snow, groundwater, river and lake can have a great effect on TWS change in Xinjiang. The interannual variation of water storage in Xinjiang has large fluctuations. The water storage rapidly dropped in 2006–2008 and quickly ascended in 2009–2010, which is identical to the extreme weathers in the period. Therefore, we conclude that GRACE can effectively monitor drought and flood events. The spatial pattern of water storage change in Xinjiang is also determined from GRACE. The water storage mostly shows a descending trend in the northern Xinjiang and a rising trend in the southern Xinjiang. Tianshan Mountain has the largest decreasing rate of water storage with the minimum rate of −4.4 mm yr−1 in the middle Tianshan. This study can help to scientifically manage water resources which can promote the ecological environment protection, and the social and economic development in Xinjiang. GRACE can monitor the global and local water storage change from space observation. With the extension of observation time and subsequent implementation of higher precision gravimetry satellites, the spatial and temporal resolution of water storage inversion will be further improved. This can provide more scientific data to optimize the allocation of water resources and promote the sustainable development of water resources in the arid area. Acknowledgements We are very grateful to the anonymous reviewers for their insightful comments and suggestions. We thank ISDC for providing GRACE data. This study is supported by the National Basic Research Program of China (grant no. 2013CB733302), the National Natural Science Foundation of China (grant nos 41374009 and 41774001), the Special Project of Basic Science and Technology of China (grant no. 2015FY310200) and the SDUST Research Fund (grant no. 2014TDJH101). REFERENCES Belda S., García-García D., Ferrándiz J.M., 2015. On the decorrelation filtering of RL05 GRACE data for global applications, Geophys. J. Int. , 200( 1), 173– 184. https://doi.org/10.1093/gji/ggu386 Google Scholar CrossRef Search ADS   Bettadpur S., 2012. UTCSR Level-2 Processing Standards Document for Level-2 Product Release 05, Rev.4.0 , Center for Space Research, University of Texas at Austin. Cao Y., Nan Z., Cheng G., 2015. GRACE gravity satellite observations of terrestrial water storage changes for drought characterization in the arid land of Northwestern China, Remote Sens. , 7( 1), 1021– 1047. https://doi.org/10.3390/rs70101021 Google Scholar CrossRef Search ADS   Chao N.F., Wang Z.T., Sun J., 2015. The inversion of terrestrial water changes by non-isotropic combination filter, Acta Geod. Cartograph. Sin. , 44( 2), 174– 182. Chen J.L., Wilson C.R., Seo K.W., 2006. Optimized smoothing of Gravity Recovery and Climate Experiment (GRACE) time-variable gravity observations, J. geophys. Res. , 111( B6), 115- 139. https://doi.org/10.1029/2005JB004064 Google Scholar CrossRef Search ADS   Chen J.L., Wilson C.R., Blankenship D., Tapley B.D., 2009. Accelerated Antarctic ice loss from satellite gravity measurements, Nature Geosci. , 2( 12), 859– 862. https://doi.org/10.1038/ngeo694 Google Scholar CrossRef Search ADS   Chen Y., Jiang Y.A., Mao W.Y., Wang S., 2011. 2009/2010 winter snow disaster in Northern Xinjiang under global warming, Adv. Clim. Chang. Res. , 7( 2), 104– 109. Cheng M., Tapley B.D., 2004. Variations in the Earth's oblateness during the past 28 years, J. geophys. Res. , 109( B9), B09402, doi:10.1029/2004JB003028. https://doi.org/10.1029/2004JB003028 Deng H., Chen Y., 2017. Influences of recent climate change and human activities on water storage variations in Central Asia, J. Hydrol. , 544, 46– 57. https://doi.org/10.1016/j.jhydrol.2016.11.006 Google Scholar CrossRef Search ADS   Duan X., Guo J., Shum C.K., van der Wal W., 2009. On the postprocessing removal of correlated errors in GRACE temporal gravity field solutions, J. Geod. , 83( 11), 1095– 1106. https://doi.org/10.1007/s00190-009-0327-0 Google Scholar CrossRef Search ADS   Fan Y., van den Dool H., 2004. Climate Prediction Center global monthly soil moisture data set at 0.5° resolution for 1948 to present, J. geophys. Res. , 109( D10), D10102, doi:10.1029/2003JD004345. https://doi.org/10.1029/2003JD004345 Google Scholar CrossRef Search ADS   Gao C.C., Lu Y., Zhang Z.Z., Shi H.L., Zhu C.D., 2015. Ice sheet mass balance in Antarctica measured by GRACE and its uncertainty, Chin. J. Geophys. , 58( 3), 780– 792. Golyandina N., Nekrutkin V., Zhigljavsky A.A., 2001. Analysis of Time Series Structure: SSA and Related Techniques , CRC Press. Google Scholar CrossRef Search ADS   Guo J., Mu D., Liu X., Yan H., Dai H., 2014. Equivalent water height extracted from GRACE gravity field model with robust independent component analysis, Acta Geophys. , 62( 4), 953– 972. https://doi.org/10.2478/s11600-014-0210-0 Google Scholar CrossRef Search ADS   Guo J., Mu D., Liu X., Yan H., Sun Z., Guo B., 2016. Water storage changes over the Tibetan Plateau revealed by GRACE mission, Acta Geophys. , 64( 2), 463– 476. https://doi.org/10.1515/acgeo-2016-0003 Google Scholar CrossRef Search ADS   Hassan A., Jin S., 2016. Water storage changes and balances in Africa observed by GRACE and hydrologic models, Geod. Geodyn. , 7( 1), 39– 49. https://doi.org/10.1016/j.geog.2016.03.002 Google Scholar CrossRef Search ADS   Hassani H., 2007. Singular spectrum analysis: methodology and comparison, J. Data Sci. , 5( 2), 239– 257. Houborg R., Rodell M., Li B., Reichle R., Zaitchik B.F., 2012. Drought indicators based on model-assimilated Gravity Recovery and Climate Experiment (GRACE) terrestrial water storage observations, Water Resour. Res. , 48( 7), 2515– 2521. https://doi.org/10.1029/2011WR011291 Google Scholar CrossRef Search ADS   Jekeli C., 1981. Alternative Methods to Smooth the Earth's Gravity Field , Department of Geodetic Science and Surveying, Ohio State University. Ju X.L., Shen Y.Z., Zhang Z.Z., 2013. Antarctic ice mass change analysis based on GRACE RL05 data, Chin. J. Geophys. , 56( 9), 2918– 2927. Ju X., Shen Y., Zhang Z., 2014. GRACE RL05-based ice mass changes in the typical regions of Antarctica from 2004 to 2012, Geod. Geodyn. , 5( 4), 57– 67. https://doi.org/10.3724/SP.J.1246.2014.04057 Google Scholar CrossRef Search ADS   Klees R., Revtova E.A., Gunter B.C., Ditmar P., Oudman E., Winsemius H.C., Savenije H.H.G., 2008. The design of an optimal filter for monthly GRACE gravity models, Geophys. J. Int. , 175( 2), 417– 432. https://doi.org/10.1111/j.1365-246X.2008.03922.x Google Scholar CrossRef Search ADS   Kusche J., 2007. Approximate decorrelation and non-isotropic smoothing of time-variable GRACE-type gravity field models, J. Geod. , 81( 11), 733– 749. https://doi.org/10.1007/s00190-007-0143-3 Google Scholar CrossRef Search ADS   Li K.M., Li Z.Q., Gao W.Y., Wang L., 2011. Recent glacial retreat and its effect on water resources in eastern Xinjiang, Chin. Sci. Bull. , 56( 33), 3596– 3604. https://doi.org/10.1007/s11434-011-4720-8 Google Scholar CrossRef Search ADS   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, Chin. J. Geophys. , 56( 6), 1843– 1849. Lu F., You W., Fan D.M., Huang Q., 2015. Chinese continental water storage and ocean water mass variations analysis in recent ten years based on GEACE RL05 data, Acta Geod. Cartograph. Sin. , 44( 2), 160– 167. Luo Z.C., Li Q., Zhong B., 2012. Water storage variations in Heihe River Basin recovered from GRACE temporal gravity field, Acta Geod. Cartograph. Sin. , 41( 5), 676– 681. Mu D.P., Guo J.Y., Sun Z.C., Kong Q.L., 2014. Equivalent water height from GRACE gravity model based on principal component analysis, Prog. Geophys. , 29( 4), 1512– 1517. Mu D., Yan H., Feng W., Peng P., 2017. GRACE leakage error correction with regularization technique: case studies in Greenland and Antarctica, Geophys. J. Int. , 208( 3), 1775– 1786. Oropeza V., Sacchi M., 2011. Simultaneous seismic data denoising and reconstruction via multichannel singular spectrum analysis, Geophysics , 76( 3), V25– V32. https://doi.org/10.1190/1.3552706 Google Scholar CrossRef Search ADS   Paulson A., Zhong S., Wahr J., 2007. Inference of mantle viscosity from GRACE and relative sea level data, Geophys. J. R. astr. Soc. , 171( 2), 497– 508. https://doi.org/10.1111/j.1365-246X.2007.03556.x Google Scholar CrossRef Search ADS   Peltier W.R., 2004. Global glacial isostasy and the surface of the ice-age Earth: the ICE-5G (VM2) model and GRACE, Annu. Rev. Earth Planet. Sci. , 32( 1), 111– 149. https://doi.org/10.1146/annurev.earth.32.082503.144359 Google Scholar CrossRef Search ADS   Plaut G., Vautard R., 1994. Spells of low-frequency oscillations and weather regimes in the Northern Hemisphere, J. Atmos. Sci. , 51( 2), 210– 236. https://doi.org/10.1175/1520-0469(1994) Google Scholar CrossRef Search ADS   Rangelova E., van der Wal W., Braun A., Sideris M.G., Wu P., 2007. Analysis of Gravity Recovery and Climate Experiment time-variable mass redistribution signals over North America by means of principal component analysis, J. geophys. Res. , 112( F3), F03002, doi:10.1029/2006JF000615. https://doi.org/10.1029/2006JF000615 Google Scholar CrossRef Search ADS   Rangelova E., van der Wal W., Sideris M.G., Wu P., 2010. Spatiotemporal analysis of the GRACE-derived mass variations in North America by means of multi-channel singular spectrum analysis, in Gravity, Geoid and Earth Observation. International Association of Geodesy Symposia , Vol. 135, no. 7, pp. 539– 546, ed. Mertikas S., Springer, Berlin, Heidelberg. https://doi.org/10.1007/978-3-642-10634-7_72 Rangelova E., Sideris M.G., Kim J.W., 2012. On the capabilities of the multi-channel singular spectrum method for extracting the main periodic and non-periodic variability from weekly GRACE data, J. Geodyn. , 54( 2), 64– 78. https://doi.org/10.1016/j.jog.2011.10.006 Google Scholar CrossRef Search ADS   Rodell M. et al.  , 2004. The global land data assimilation system, Bull. Am. Meteor. Soc. , 85( 3), 381– 394. https://doi.org/10.1175/BAMS-85-3-381 Google Scholar CrossRef Search ADS   Seoane L., Ramillien G., Frappart F., Leblanc M., 2013. Regional GRACE-based estimates of water mass variations over Australia: validation and interpretation, Hydrol. Earth Syst. Sci. , 17( 12), 4925– 4939. https://doi.org/10.5194/hess-17-4925-2013 Google Scholar CrossRef Search ADS   Shen Y., Guo J., Liu X., Kong Q., Guo L., Li W., 2018. Long-term prediction of polar motion using a combined SSA and ARMA model, J. Geod. , 2018, doi:10.1007/s00190-017-1065-3. Steffen H., Denker H., Müller J., 2008. Glacial isostatic adjustment in Fennoscandia from GRACE data and comparison with geodynamical models, J. Geodyn. , 46( 3–5), 155– 164. https://doi.org/10.1016/j.jog.2008.03.002 Google Scholar CrossRef Search ADS   Swenson S., Wahr J., 2006. Post-processing removal of correlated errors in GRACE data, Geophys. Res. Lett. , 33( 8), L08402, doi:10.1029/2005GL025285. https://doi.org/10.1029/2005GL025285 Google Scholar CrossRef Search ADS   Tamisiea M.E., 2011. Ongoing glacial isostatic contributions to observations of sea level change, Geophys. J. Int. , 186( 3), 1036– 1044. https://doi.org/10.1111/j.1365-246X.2011.05116.x Google Scholar CrossRef Search ADS   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, J. geophys. Res. , 103( B12), 30205– 30229. https://doi.org/10.1029/98JB02844 Google Scholar CrossRef Search ADS   Wang H., Xiang L., Jia L., Wu P., Steffen H., Jiang L., Shen Q., 2015. Water storage changes in North America retrieved from GRACE gravity and GPS data, Geod. Geodyn. , 6( 4), 267– 273. https://doi.org/10.1016/j.geog.2015.07.002 Google Scholar CrossRef Search ADS   Wang X., de Linage C., Famiglietti J., Zender C.S., 2011. Gravity Recovery and Climate Experiment (GRACE) detection of water storage changes in the Three Gorges Reservoir of China and comparison with in situ measurements, Water Resour. Res. , 47( 12), W12502, doi:10.1029/2011WR010534. https://doi.org/10.1029/2011WR010534 Google Scholar CrossRef Search ADS   Zhang J.B., 2001. Strategic countermeasures for sustainable use of water resources of Xinjiang, Arid Land Geogr. , 24( 3), 217– 222. Zotov L.V., Shum C.K., 2010. Multichannel singular spectrum analysis of the gravity field data from GRACE satellite, in Proceedings of the 4th Gamow International Conference on Astrophysics and Cosmology after Gamow and the 9th Gamow Summer School “Astronomy and Beyond: Astrophysics, Cosmology, Radio Astronomy, High Energy Physics and Astrobiology”, Ukraine , Vol. 1206, no. 1, pp, 473– 479, eds Chakrabarti S.K., Bisnovatyi-Kogan G.S., Zhuk A.I., AIP Publishing. Zotov L.V., 2012. Application of multichannel singular spectrum analysis to geophysical fields and astronomical images, Adv. Astron. Space Phys. , 1( 2), 82– 84. © The Author(s) 2018. Published by Oxford University Press on behalf of The Royal Astronomical Society. http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png Geophysical Journal International Oxford University Press

# Terrestrial water storage changes over Xinjiang extracted by combining Gaussian filter and multichannel singular spectrum analysis from GRACE

, Volume 213 (1) – Apr 1, 2018
11 pages

/lp/ou_press/terrestrial-water-storage-changes-over-xinjiang-extracted-by-combining-blAKjL1v0l
Publisher
Oxford University Press
ISSN
0956-540X
eISSN
1365-246X
D.O.I.
10.1093/gji/ggy006
Publisher site
See Article on Publisher Site

### Abstract

Summary Water resource management is crucial for the economic and social development of Xinjiang, an arid area located in the Northwest China. In this paper, the time variations of gravity recovery and climate experiment (GRACE)-derived monthly gravity field models from 2003 January to 2013 December are analysed to study the terrestrial water storage (TWS) changes in Xinjiang using the multichannel singular spectrum analysis (MSSA) with a Gaussian smoothing radius of 400 km. As an extended singular spectrum analysis (SSA), MSSA is more flexible to deal with multivariate time-series in terms of estimating periodic components and trend, reducing noise and identifying patterns of similar spatiotemporal behaviour thanks to the data-adaptive nature of the base functions. Combining MSSA and Gaussian filter can not only obviously remove the north–south striping errors in the GRACE solutions but also reduce the leakage errors, which can increase the signal-to-noise ratio by comparing with the traditional procedure, that is, empirical decorrelation method followed with the Gaussian filtering. The spatiotemporal characteristics of TWS changes in Xinjiang were validated against the Global Land Dynamics Assimilation System, the Climate Prediction Center and in-situ precipitation data. The water storage in Xinjiang shows the relatively large fluctuation from 2003 January to 2013 December, with a drop from 2006 January to 2008 December due to the drought event and an obvious rise from 2009 January to 2010 December because of the high precipitation. Spatially, the TWS has been increasing in the south Xinjiang, but decreasing in the north Xinjiang. The minimum rate of water storage change is −4.4 mm yr−1 occurring in the central Tianshan Mountain. Hydrogeophysics, Satellite gravity, Time variable gravity, Inverse theory, Spatial analysis 1 INTRODUCTION Xinjiang is located in one inland arid area in the Northwest China (35°25΄N-49°10΄N and 73°40΄E- 96°23΄E, see Fig. 1). The scarce water resource severely restricts the local economic and social development. For example, agriculture, as the most important part of Xinjiang economy, highly relies on the water resource management. Xinjiang is divided into the southern and the northern parts by Tianshan Mountain. There are more than 570 rivers in its territory. The total water resources are about 882 × 108 m3 including the surface water resources of about 794 × 108 m3 and the ground water resources of about 88 × 108 m3. The water production per unit land area is only 5.30 × 104 m3 km−2 which accounts for only 18 per cent of the average water production in China (Zhang 2001). Study on terrestrial water storage (TWS) changes can help to scientifically optimize the allocation of water resources, protect the ecological environment and promote the sustainable development of water resources in Xinjiang. Figure 1. View largeDownload slide Geographical distribution of Xinjiang and meteorological stations. Figure 1. View largeDownload slide Geographical distribution of Xinjiang and meteorological stations. The gravity recovery and climate experiment (GRACE) satellite was successfully launched in 2002 March. It is jointly sponsored by the National Aeronautics and Space Administration (NASA) of USA and the Deutsches Zentrum für Luft- und Raumfahrt. GRACE can probe the temporal variations of Earth gravity field from which the Earth surface mass balance can be extracted (Wang et al. 2015). The Stokes coefficients estimated from GRACE mission are widely used in geophysical, glaciological, oceanographic and hydrological studies due to its unprecedented accuracy. For example, GRACE data are used extensively to estimate the TWS changes and groundwater variations over the large-scale areas (Luo et al. 2012; Seoane et al. 2013; Lu et al. 2015; Guo et al. 2016; Hassan & Jin 2016; Deng & Chen 2017), monitor local drought and flood events (Houborg et al. 2012; Li et al. 2013; Cao et al. 2015) and determine mass changes of the Antarctic ice sheet (Chen et al. 2009; Ju et al.2013, 2014; Gao et al. 2015; Mu et al. 2017). Many other applications such as the sea level change (Tamisiea 2011) and the glacial isostatic adjustment (GIA, Steffen et al. 2008) are also widely studied using GRACE data. GRACE-derived high-degree Stokes coefficients are contaminated by random noises which make the GRACE TWS changes contain the obvious north–south striping errors in the spatial domain and the high-frequency errors in the frequency domain. These errors are usually be removed or reduced by the Gaussian filtering technique (Wahr et al. 1998), the decorrelation procedure (Swenson & Wahr 2006), the DDK filter (Kusche 2007), the Wiener filter (Klees et al. 2008), the principal component analysis (PCA, Rangelova et al. 2007; Mu et al. 2014), the singular spectrum analysis (SSA, Wang et al. 2011) or the independent component analysis (Guo et al. 2014), etc. These methods have achieved good results, but there are still residual errors in the water storage series. In this paper, we combine the Gaussian filter and multichannel singular spectrum analysis (MSSA) to denoise the north–south striping errors in GRACE solutions. MSSA as one data-adaptive method without the prior information can divide periodic components, separate trend, filter out noise and distinguish patterns of similar spatiotemporal behaviours (Zotov 2012). We apply MSSA to the GRACE Stokes coefficients and then implement the Gaussian filtering technique to better separate noises and real signals. The effectiveness of this procedure is demonstrated by comparing to the results of traditional method, that is, the empirical decorrelation method (EDM, Duan et al. 2009) combined with the Gaussian filter. Using this alternative method, we extract TWS changes in Xinjiang which are then verified against the hydrological models like the Global Land Dynamics Assimilation System (GLDAS), the Climate Prediction Center (CPC) and in-situ precipitation data, and then the temporal and spatial changes of TWS in Xinjiang are analysed. 2 DATA AND METHODOLOGY 2.1 GRACE products and water storage estimates Stokes coefficients determined from GRACE data released by the Center for Space Research, the Texas University at Austin, were used in the study. The monthly GRACE level-2 products are downloaded from http://isdc.gfz-potsdam.de/grace through January 2003 to December 2013 with the exceptions of 2003 June, 2011 January and June, 2012 May and October, and 2013 March, August and September due to lack of data. These losing data are interpolated from the neighbouring monthly solutions with the cubic spline interpolation method. The GRACE Level-2 RL05 data are truncated at maximum degree and order 60 without the effects of atmosphere, ocean and tide (Bettadpur 2012). The degree-1 coefficients are not restored as their impact on this regional study is rather limited. Because the C20 coefficient measured by the satellite laser ranging (SLR) is much better than that in RL05, C20 in RL05 was replaced with those derived from the SLR data (Cheng & Tapley 2004). The terrestrial mass balance inversed from GRACE can effectively indicate TWS change with the equivalent water height (EWH, Chao et al. 2015). The EWH determined from GRACE products is   \begin{eqnarray} \Delta h(\theta ,\phi ) &=& \frac{{a{\rho _{{\rm{ave}}}}}}{{3{\rho _w}}}\sum\limits_{l = 0}^\infty \sum\limits_{m = 0}^l {{\overline p }_{lm}}(\cos \theta )\frac{{2l + 1}}{{1 + {k_l}}}\nonumber\\ &&\times \,(\Delta {C_{lm}} \cos m\phi+\Delta {S_{lm}}\sin m\phi) \end{eqnarray} (1)where a is the averaged Earth radius, ρave is the mean Earth density (5517 kg m−3), ρw = 1000 kg m−3, θ and ϕ are the colatitude and longitude, respectively, $${\overline p _{lm}}(\cos \theta )$$ is the fully normalized associated Legendre function with degree l and order m, kl is the load Love number of degree l and ΔClm and ΔSlm are the fully normalized Stokes coefficient changes determined by GRACE with respect to the mean gravity from 2003 January to 2013 December. Combining the Gaussian filter and MSSA (GMSSA) can weaken and remove the effects of satellite orbit errors and Stokes coefficient errors in GRACE products. The Gaussian smoothing kernel function (Jekeli 1981; Wahr et al. 1998) is introduced to eq. (1) into get the filtered EWH as   \begin{eqnarray} \Delta h(\theta ,\phi ) &=& \frac{{a{\rho _{{\rm{ave}}}}}}{{3{\rho _w}}}\sum\limits_{l = 0}^\infty \sum\limits_{m = 0}^l {{\overline p }_{lm}}(\cos \theta )\frac{{2l + 1}}{{1 + {k_l}}}{w_l}\nonumber\\ &&\times \,(\Delta {{\overline C }_{lm}} \cos m\phi +\Delta {\overline S _{lm}}\sin m\phi ) \end{eqnarray} (2)where wl is the Gaussian smoothing kernel function for degree l. $$\Delta {\overline C _{lm}}$$ and $$\Delta {\overline S _{lm}}$$ are the spherical harmonic coefficient changes reconstructed by MSSA. The TWS changes in Xinjiang are calculated from eq. (2) on grid of 1° × 1°, and the results are further improved using the ICE-5G model (Peltier 2004) to remove the effect of GIA which is small in the study region. The water storage changes over the whole Xinjiang are then estimated. 2.2 Hydrological models Two hydrological models, GLDAS and CPC, are used in the study. GLDAS (Rodell et al. 2004) is based on the Noah terrestrial ground model with the time resolution of 1 month and the spatial resolution of 0.25° × 0.25° from 2003 to 2013. The model is released by the Goddard Space Flight Center of NASA and the National Center of Environmental Prediction, USA. The water storage change can be computed from the soil moisture variations up to four layers (0–0.1, 0.1–0.4, 0.4–1 and 1–2 m), snow water change and canopy water change with the hydrological model. CPC model (Fan & van den Dool 2004) is released by the CPC of National Oceanic and Atmospheric Administration, USA. The time and spatial resolutions of the model are one month and 0.5° × 0.5°, respectively. The model can provide the soil moisture up to depth of 0–1.6 m caused by water balance changes. GLDAS and CPC data are all converted to spherical harmonic coefficients up to degree and order 60 and then filtered using the same strategy applied to the GRACE data. 2.3 Precipitation The monthly precipitations collected by 54 weather stations (see Fig. 1) located in Xinjiang from 2003 to 2013 are used in the study. The Kriging method is used to interpolate the precipitation data with grid of 1° × 1° over Xinjiang. The latitude cosine weighting method is used to calculate the mean precipitation for the whole Xinjiang, which can be compared with the GRACE-derived results. 2.4 MSSA method MSSA as an extended SSA can be used to extract the detailed temporal and spectral information from the multidimensional time-series (Plaut & Vautard 1994; Shen et al. 2018). It is also often assumed to be synonymous with the extended Empirical Orthogonal Function (EOF). Instead of using the sample correlation matrix in EOF/PCA, MSSA forms the trajectory matrix with a time lag to be analysed, and different channels can help each other to capture spatiotemporal correlation patterns. MSSA is more flexible to recognize trend, modulate oscillations of different periods and denoise multidimensional time-series (Zotov & Shum 2010; Oropeza & Sacchi 2011). It will be useful for filtering out strips and high-frequency noise in the GRACE data (Rangelova et al. 2010, 2012). 3 MSSA ON GRACE SPHERICAL HARMONICS There is one data set x (xn,l:n=1,N;l=1,L), which consists of N observations (time epochs) each with L variables (channels). It can be processed by MSSA (Oropeza & Sacchi 2011) with the following procedure. (1) The data set x is the time-series of Stokes coefficient {ΔCij(tn), ΔSij(tn)}, n = 1, …, N, and we select the lag parameter M (a lag window) which satisfies $$1 \le M \le \frac{N}{2}$$, where N = 132 and M = 60 in the GRACE data processing. So the trajectory matrix of ΔCij channel is   $${{\bf X}_{\Delta {C_{ij}}}} = \left( {\begin{array}{@{}*{4}{c}@{}} {\Delta {C_{ij}}({t_1})}&{\Delta {C_{ij}}({t_2})}& \quad\cdots &\quad{\Delta {C_{ij}}({t_K})}\\ {\Delta {C_{ij}}({t_2})}&{\Delta {C_{ij}}({t_3})}& \quad\cdots &\quad{\Delta {C_{ij}}({t_{K + 1}})}\\ \vdots & \vdots & \quad\cdots & \vdots \\ {\Delta {C_{ij}}({t_M})}&\quad{\Delta {C_{ij}}({t_{M + 1}})}& \quad\cdots &\quad{\Delta {C_{ij}}({t_N})} \end{array}} \right)$$ (3)where K = N − M + 1. For GRACE data, the trajectory matrix X composed of such blocks for every channel ΔCij and ΔSij is   $${\bf X} = {\left[ {{{\bf X}_{\Delta {C_{20}}}},{{\bf X}_{\Delta {S_{20}}}}, \ldots ,{{\bf X}_{\Delta {C_{ij}}}},{{\bf X}_{\Delta {S_{ij}}}}, \ldots {{\bf X}_{\Delta {C_{60}}}},{{\bf X}_{\Delta {S_{60}}}}} \right]^{\rm{T}}}$$ (4)in which X is one LM × K matrix. (2) Matrix X can be decomposed by the singular value decomposition, that is, X = USVT, in which S is one diagonal matrix with D diagonal elements. Here, D = min (LM, K). These elements are singular values which are sorted in descending order (see Fig. 2). Figure 2. View largeDownload slide Singular values for MSSA on GRACE products. Figure 2. View largeDownload slide Singular values for MSSA on GRACE products. The big singular value stands for the part with large power spectrum and the small corresponds to noises. The trajectory matrix X can be synthesized by the elementary matrices as   $${\bf X} = {{\bf T}_{\rm{1}}} + {{\bf T}_{\rm{2}}} + \cdots + {{\bf T}_{\rm{D}}}$$ (5)where the ith component corresponds to the matrix $${{\bf T}_{i}} = {s_{i}}{{\bf u}_{i}}{\bf v}_{i}^{T}$$ in which ui and vi are the left singular vector and the right singular vector for the singular value si of matrix X. (3) Since GRACE data time-series includes different sorts of signals such as long-term, annual and semi-annual terms and noises. Ti representing the same signal can be classified into one group using the ω-correlation method (Hassani 2007). The ω-correlation is   $$\rho _{i,j}^w = \frac{{({Y^{(i)}},{Y^{(j)}})}}{{{{\left\| {{Y^i}} \right\|}_w}{{\left\| {{Y^j}} \right\|}_w}}}$$ (6)where Yi is the time-series reconstructed from seriesTi using eq. (8), $${\| {{Y^i}} \|_w} = \sqrt {({Y^{(i)}},{Y^{(i)}})}$$ and $$({Y^{(i)}},{Y^{(j)}}) = \sum_{k = 1}^N {{w_k}y_k^iy_k^j}$$ in which wk = min (k, M, N − k). The large absolute $$\rho _{i,j}^w$$ indicates Ti and Tj may be correspond to the same signal so that these two components should be grouped together. Fig. 3 shows the absolute $$\rho _{i,j}^w$$ for the first 30 Ti in processing GRACE products in the manuscript. From Fig. 3, we can find when i, j > 12 the absolute $$\rho _{i,j}^w$$ are all relatively high, which means these series cannot be well separated and maybe contain noises. So these first 12 singular values and their corresponding singular vectors are used to reconstruct the interesting signals. The first two components correspond to the annual cycle, the eighth and ninth components are composed to denote the semi-annual signal, and the rest of the first 12 modes represent the long-term single, other useful signals and some noises. This process is equivalent to denoising and filtering GRACE products. The trajectory matrix reconstructed from the first 12 components is   \begin{eqnarray} {{\bf X}^{\boldsymbol{'}}} &=& {{\bf T}_1} + {{\bf T}_2} + \cdots + {{\bf T}_{12}}{\rm} \nonumber\\ &=& {\left[ {{{{\bf X'}}_{\Delta {C_{20}}}},{{{\bf X'}}_{\Delta {S_{20}}}}, \ldots ,{{{\bf X'}}_{\Delta {C_{ij}}}},{{{\bf X'}}_{\Delta {S_{ij}}}}, \ldots {{{\bf X'}}_{\Delta {C_{60}}}},{{{\bf X'}}_{\Delta {S_{60}}}}} \right]^{\rm{T}}} \end{eqnarray} (7) Figure 3. View largeDownload slide ω-correlation for the first 30 reconstructed components. Figure 3. View largeDownload slide ω-correlation for the first 30 reconstructed components. (4) Denoised data can be reconstructed from the trajectory matrix X΄ (Golyandina et al. 2001). Let $${\bf Z} = {\bf X}_i^\prime$$(i = ΔC20, ΔS20, …ΔS60), then the reconstructed series dk(ΔCij,ΔSij) is   $${d_k} = \left\{ \begin{array}{@{}l@{}} \frac{1}{k}\sum\limits_{n = 1}^k {{z_{n,k - n + 1}}} \\ \frac{1}{{{M^*}}}\sum\limits_{n = 1}^{{M^*}} {{z_{n,k - n + 1}}} \\ \frac{1}{{N - k + 1}}\sum\limits_{n = k - {K^*} + 1}^{N - {K^*} + 1} {{z_{n,k - n + 1}}} \end{array} \right.\quad \begin{array}{@{}*{1}{l}@{}} {1 \le k < {M^*}}\\ {}\\ {{M^*} \le k \le {K^*}}\\ {}\\ {{K^*} < k \le N}\\ {} \end{array}$$ (8)where M* = min (M, K), K* = max (M, K). 4 RESULTS AND ANALYSIS 4.1 Method validation Fig. 4 illustrates the effect of different methods on the GRACE monthly EWH anomaly in 2007 April. The top subfigures show the original fields with different Gaussian smoothing radii. Without smoothing the map is dominated by the prominent north–south stripes and little geophysical signal can be observed. Using the smoothing radius of 300 km, hydrological signals over large basins (e.g. Amazon Valley, Congo River Basin and Greenland) start to stand out, but are still corrupted with the stripes. When the radius increases to 500 km, most of the stripes over land are removed and the oceans appear mostly stripe free when using a smoothing radius of 700 km. While increasing the smoothing radius, the magnitude of signals (e.g. Amazon Valley) is reduced and the leakage errors (e.g. Antarctica) increase. The second-row subfigures display the fields after applying MSSA combining with different Gaussian filtering radius. The third-row subfigures display the fields after applying EDM combining with different Gaussian filtering radius (GEDM for short) to the Stokes coefficients. The bottom subfigures show the differences between the figures in the second and third rows. Without the Gaussian filtering, MSSA and EDM both can significantly reduce the stripes, and the geophysical signal over large basins can be easily observed, although some noisy patterns remain exist. For MSSA filtering only, the stripes mainly distributed in low/mid-latitudes and slight in high latitudes. However, in the EDM result the stripes are mainly concentrated in low/mid-latitudes, and the decorrelation affects geophysical signals, especially at high latitudes (e.g. Antarctica and Greenland) where the magnitudes of signals are evidently smaller than that in the MSSA result. Moreover, with the same Gaussian smoothing radius (e.g. 300 km), GMSSA can preserve significantly more signals and show less pronounced leakage effects. Examples of areas where this effect is clearly visible are Antarctica, the Amazon River Basin, the Yukon and Fraser Basins (located in Alaska and western Canada), Madagaskar et al. For instance, due to less leakage effects, the GMASS results better coincide with coastlines in the western Canada than those in the GEDM. Figure 4. View largeDownload slide Global EWHs in 2007 April with various Gaussian smoothing radii. (a) 0 km, (b) 300 km, (c) 500 km and (d) 700 km. The upper figures show the unfiltered fields, the middle figures represent the MSSA and EDM method respectively, and the bottom figures display the difference between MSSA and EDM. Units are in mm. Figure 4. View largeDownload slide Global EWHs in 2007 April with various Gaussian smoothing radii. (a) 0 km, (b) 300 km, (c) 500 km and (d) 700 km. The upper figures show the unfiltered fields, the middle figures represent the MSSA and EDM method respectively, and the bottom figures display the difference between MSSA and EDM. Units are in mm. We also use GMASS and GEDM with 400 km smoothing radius to extract the global EWHs in 2007 April and October because the land EWHs of these two months are violent and converse over low/mid-latitude. Fig. 5 shows the inversed results. The upper figures represent the results of GEDM and GMSSA in 2007 April, respectively, and the bottom figures are the results of 2007 October correspondingly processed by GEDM and GMSSA. From Fig. 5, we can find that these two processing strategies give consistent results without obvious north–south striping noises. Table 1 lists the statistical results. The root mean squares (rms) of EWH estimated by these two strategies are almost identical. The maximum EWH from GEDM is only 286.3 mm which is less than that from GMSSA, but the minimum values are basically identical in 2007 April. The maximum EWH from GEDM is only 216.3 mm which is greater than that from GMSSA, and the minimum value is less than that determined by GMSSA in 2007 October. Moreover, the mean EWH from GMSSA is larger than that from GEDM in 2007 April, but it appears small comparing to the GEDM in 2007 October. The idea of EDM is to fit a quadratic polynomial in a moving window to the Stokes coefficients of even and odd degrees for a particular order and remove this from the original Stokes coefficients. MSSA uses the Stokes coefficients to construct a trajectory matrix which is decomposed and reconstructed so that the signals representing significant geophysical can be extracted. The different processing strategies lead to different results listed in Table 1. Figs 4 and 5, and Table 1 can validate the feasibility and effectiveness of GMSSA by comparing with GEDM. Figure 5. View largeDownload slide Global equivalent water storage changes in 2007 April and October processed by (a) GEDM, and (b) GMSSA. Figure 5. View largeDownload slide Global equivalent water storage changes in 2007 April and October processed by (a) GEDM, and (b) GMSSA. Table 1. Statistics of water storage changes from GEDM and GMSSA with filtering radius of 400 km (unit in mm). Month  Method  Max  Min  Mean  rms  2007 April  GEDM  286.3  −205.4  2.6  44.1    GMSSA  327.4  −204.4  4.3  44.3  2007 October  GEDM  216.3  −425.1  1.0  43.2    GMSSA  196.0  −409.4  0.0  40.9  Month  Method  Max  Min  Mean  rms  2007 April  GEDM  286.3  −205.4  2.6  44.1    GMSSA  327.4  −204.4  4.3  44.3  2007 October  GEDM  216.3  −425.1  1.0  43.2    GMSSA  196.0  −409.4  0.0  40.9  View Large As another tool to analyse the effectiveness of GMSSA, we use the method suggest by Chen et al. (2006). The idea is based on the fact that GRACE measurement errors are approximately at the same level over both land and ocean, but the surface mass variability is stronger in the continents than that in the oceans. The optimal filter should maximize the quotient of the latitude weighted rms of the continental and oceanic signals which represent the signal-to-noise ratio. Additionally, in order to reduce the leakage of the continental signal only the ocean points farther than 400 km from the coast are included. As a result, the rms ratio from GEDM is 1.45 which is lower than that from GMSSA, that is, 1.72. So GMSSA is able to obtain greater rms ratio and retain more interesting geophysical signals which can improve the precision and reliability of inversed water storage. It is inevitable that some useful geophysical signals will be filtered out together with noises in the filtering procedure. We will evaluate to which extent the filter might remove true geophysical signal. Here, we address this problem based on the study of rms of the linear trend of surface mass changes from GRACE and noise-free synthetic data. The trend is determined by least-squares fitting the GRACE and synthetic monthly Stokes coefficients during 2003 January to 2013 December with an offset, a linear term that represents the trend, an annual term and a half-yearly term as Rangelova et al. (2012). In addition, the synthetic data consist of global grids with ocean bottom pressure over the oceans and with GLDAS over land for the same period as GRACE, and a surface mass change rate corresponding to the gravity variation of a GIA model from Paulson et al. (2007) have been added back to all grids. As GLDAS simulated values over ice covered regions (e.g. Greenland and Antarctica) are not realistic, these regions are avoided in the comparison with GRACE. The synthetic data are converted to spherical harmonic coefficients up to degree 60, and the degree-1 coefficients are set to zero. Due to the influence of the noises in GRACE data, the rms of trend is higher in the unfiltered data but lower after using a filtering technique. However, an excessive reduction of rms would mean an undesirable attenuation of the geophysical signals. The optimal filter should minimize the rms of the global linear trend and maximize it in the noise-free synthetic data. Besides, the filters attenuate the signal in both data, but noises reduce only in real GRACE data. As the GLDAS data used here does not include groundwater and separate surface water components (such as rivers and lakes) which make the linear trend rms from synthetic data lower than those from GRACE data. The optimal filter also should show the difference rms (Drms) of linear trend between GRACE and synthetic data minimum (Belda et al. 2015). We present in Fig. 6 the trend result computed by no-filtering, GEDM and GMSSA, respectively, The GEDM and GMSSA both use Gaussian radii 400 km. As the result in Fig. 6(a) is computed without filtering, it shows obvious north–south stripes noises and a little more geophysical signals over large basins. Comparing to the GEDM results that present in Fig. 6(b), the GMASS results (Fig. 6c) are better coincide with coastlines which can be easily found in the Greenland, Antarctica, Scandinavia, etc. which indicates that GMSSA can more effectively reduce the signal leakage errors than GEDM. Over Antarctic, there are some differences between the results of GMSSA and GEDM, for example, we can see the signals of the apparent highs and lows at the Amundsen Sea Embayment in the unfiltering and GMSSA results, but it does not appear in the GEDM results. Table 2 lists rms of the trend of global mass changes except Greenland and Antarctic applying GEDM and GMSSA filtering, respectively. We see that different approaches do not lead to significantly different rms for the GRACE data. When applying the GEDM and GMSSA to the noise-free synthetic data, the rms value of later is larger than that of former, which means that the GMSSA can preserve more signals. Furthermore, the Drms value between the result of GRACE and the synthetic result from GMSSA is 6.46, which is slightly smaller than that from GEDM (6.79). The results of Fig. 6 and Table 2 verify the reliability of the new method. Figure 6. View largeDownload slide GRACE global EWH rate during 2003–2013 determined from (a) unfilter, (b) GEDM and (c) GMSSA. Figure 6. View largeDownload slide GRACE global EWH rate during 2003–2013 determined from (a) unfilter, (b) GEDM and (c) GMSSA. Table 2. Rms of the trend of global mass changes except Greenland and Antarctic (unit in mm). Data  Unfilter  GEDM  GMSSA  GRACE  15.22  8.84  8.83  Synthetic  4.15  2.05  2.37  Data  Unfilter  GEDM  GMSSA  GRACE  15.22  8.84  8.83  Synthetic  4.15  2.05  2.37  View Large 4.2 Analysis on TWS changes in Xinjiang The mean water storage change series over Xinjiang from GRACE, GLDAS and CPC are shown in Fig. 7. We can find that these three series have the obvious seasonal variations with the maximum in summer (June–August) and the minimum in winter (December–February). GRACE-derived water storage changes have the great anomalies in 2008 and 2009 which are significantly smaller than other years, and the anomalies can also be found in GLDAS and CPC results. This phenomenon is mainly caused by the continuous high temperature and less precipitation from 2008 April to 2009 December in Xinjiang (Cao et al. 2015). After a continuous drought event taking place from summer to winter in 2008, Xinjiang again experienced severe drought in 2009, which led to the minimum annual precipitation since 1961. However, from 2009 December to 2010 February, the strongest snowstorm in the last 60 yr caused the whole 2010 year's annual precipitation exceeding the historical extremes (Chen et al. 2011). Those results in the water storage rapidly increase in the study area that can be obviously found in the result of GRACE, GLDAS and CPC in Figs 7 and 8. This indicates that GRACE has a good prospect in detecting drought and flood events. We also find that there was more obviously amplitude anomaly in 2008 November determined from GRACE than that from GLDAS and CPC, which can further indicate that GRACE has the ability to monitor the regional water storage changes. Figure 7. View largeDownload slide TWS changes in Xinjiang from GRACE, GLDAS and CPC. Figure 7. View largeDownload slide TWS changes in Xinjiang from GRACE, GLDAS and CPC. Figure 8. View largeDownload slide Interannual variations of water storage after removing the annual and semi-annual signals. Figure 8. View largeDownload slide Interannual variations of water storage after removing the annual and semi-annual signals. Over the course of the research time, the water storage change amplitude in Xinjiang, estimated by GRACE data, drifted from −43.5 to 31.2 mm. The maximum was present in 2011 May and the minimum took place in 2008 November. In particular, the amplitude derived from GRACE is obviously greater than those from GLDAS and CPC. These differences among the estimates can either be attributed to model deficiencies, such as missing surface and groundwater components in GLDAS and inadequate snow water equivalent in CPC, or due to uncertainties related to the GRACE data, like aliasing, data processing, leakage error, etc. Moreover, comparing GRACE results with hydrological models shows a phase lag of around 1–2 month which is mainly caused by the soil water interception and the regulation and storage capacity of lake, glacier and groundwater. The correlation coefficient between GRACE and GLDAS is 0.86 which is greater than that between GRACE and CPC of 0.72 within a 95 per cent confidence interval. This may be attributed to variances in the model structure, the forcing data choices and the components of the TWS changes. Although the hydrological models are not a perfect reproduction of TWS changes, it is sufficiently enough to capture the magnitude and variability of terrestrial hydrology in the study area. The strong correlations between GRACE and these two models can verify the reliability of GRACE-inversed results. The annual (first and second components) and semi-annual (eighth and ninth components) signals are removed in the GRACE MSSA filtering process in order to study the interannual variation of water storage in Xinjiang (Fig. 8). From Figs 7 and 8, we can find that the water storage in Xinjiang has the relatively large fluctuation form 2003 to 2013. To better illustrate the TWS changes in Xinjiang, we divided the research time span into four short intervals according to our visual detection of the variations in the TWS. The first interval spans from 2003 January to 2005 December during which the water storage presented the rising trend. The second interval covers the period from 2006 January to 2008 December during which the water storage dropped quickly. The third interval lasts from 2009 January to 2010 December during which the water storage rose obviously. The last interval starts from 2011 January to 2013 December during which the water storage first decreased and then swelled. GRACE-derived interannual change agreed relatively well with those from hydrological models, although there are small deviations in some time quantum, such as the interval spans from 2003 January to 2005 December. However, The GRACE-derived results are much notable, for example, the drought event in 2008 and the high precipitation in 2010 are all obviously detected by GRACE. The water storage changes determined from GRACE excluding the annual and semi-annual signals are fitted with the least-squares method to estimate the water storage rates in Xinjiang from 2003 to 2013 whose spatial distribution is shown in Fig. 9. The results show that the water storage overall appears the rising trend in the south Xinjiang with the maximum rising rate of 9.5 mm yr−1 and the descending trend in the north Xinjiang with the minimum rate of −4.4 mm yr−1. The water storage shows the gradual increase from the north to the south in the south Xinjiang. Figure 9. View largeDownload slide TWS change rates in Xinjiang from GRACE. Figure 9. View largeDownload slide TWS change rates in Xinjiang from GRACE. Contrary to the south Xinjiang, the north Xinjiang wholly has the descending trend of water storage. The descending rates of water storage are bigger close to Tianshan. This may be certainly relative to the great groundwater exploitation for the agriculture and living water demand over the piedmont region. But the water storage has the rising trend on the north Altai area. The water storage wholly appears the descending trend over the Tianshan except slightly rising in the east area. The central Tianshan has the minimum rate of −4.4 mm yr−1. Tianshan breeds the abundant water resources including glaciers and permanent snows which seriously affect the regional water storage change. So the descending water storage in Tianshan can indirectly respond to the accelerated melting of glaciers and snows (Li et al. 2011). 4.3 Comparison of GRACE-derived water storage and actual precipitation The regional water mass balance can be depicted by   $$\frac{{dW(t)}}{{dt}} = {\rm{precipitation}}\left( t \right){\rm{ }} - {\rm{ evaporation}}\left( t \right){\rm{ }} - {\rm{ runoff}}\left( t \right)$$ (9)where dW(t)/dt is the change in water content. According to eq. (10), the precipitation is an important factor causing the TWS change. However, there is no one-to-one relation between water storage and precipitation, as the evaporation and runoff can also lead to changes. In the arid land of Xinjiang, surface and groundwater resources are scarce, and precipitation is the primary water source. Thus analysing precipitation and water storage information on their overall relationship is important to interpret the GRACE results. Fig. 10 shows the actual precipitation and the GRACE-derived water storage change in Xinjiang. We can find that the water storage and the precipitation all have the obvious seasonal variations with the identical change procedure. The precipitations mainly concentrated from May to September during which the GRACE-derived water storage was high, and there are less precipitations from October to February of the following year during which low water storage was detected by GRACE. In 2008 and 2009, the precipitations were significantly less than the same period. Correspondingly, the GRACE-derived water storages were also very less at the same period. In particular, we can see that precipitation is larger in 2007 than that in 2006. But from 2006, the GRACE-derived water storage is descending. As the evaporation and runoff both play an important role in regional water allocations. Snow/ice melt water were also a crucial factor to supply the local water storage besides precipitation in Xinjiang. The GRACE results reflect the TWS changes that caused by all factors. Despite less precipitation, the joint effects of all factors result in a higher TWS in 2006 than that of 2007. Over the course of the research time, the GRACE TWS changes showed a relatively good agreement with precipitation, which could be explained that the precipitation is a main factor to control the water storage in Xinjiang. Figure 10. View largeDownload slide TWS and precipitation changes in Xinjiang from 2003 to 2013. Figure 10. View largeDownload slide TWS and precipitation changes in Xinjiang from 2003 to 2013. 5 CONCLUSIONS In this paper, we combine MSSA with Gaussian filter to process the GRACE products to extract the water storage signals over Xinjiang, a northwest arid area of China, from 2003 to 2013. GMSSA can efficiently denoise the high-degree Stokes coefficients and decorrelate GRACE data which reduce the north–south striping in the GRACE EWH fields. In addition, this method can also reduce the leakage errors thanks to the enhanced signal-to-noise ratio by comparing with GEDM, so it can retain more useful geophysical signals from GRACE products. However, the method can also inevitably remove some interesting signals. How to reasonably determine the window M and the number of singular values in MSSA and optimally select the filtering window in the Gaussian filter should be further studied. TWS change has obvious seasonal and periodic variations with the maximum in summer and the minimum in winter. The GRACE-derived water storage is basically consistent with those from these two hydrological models and actual precipitation, but the amplitudes of water storage change from GRACE are greater than those from GLDAS and CPC, which can indicate that the glacier, snow, groundwater, river and lake can have a great effect on TWS change in Xinjiang. The interannual variation of water storage in Xinjiang has large fluctuations. The water storage rapidly dropped in 2006–2008 and quickly ascended in 2009–2010, which is identical to the extreme weathers in the period. Therefore, we conclude that GRACE can effectively monitor drought and flood events. The spatial pattern of water storage change in Xinjiang is also determined from GRACE. The water storage mostly shows a descending trend in the northern Xinjiang and a rising trend in the southern Xinjiang. Tianshan Mountain has the largest decreasing rate of water storage with the minimum rate of −4.4 mm yr−1 in the middle Tianshan. This study can help to scientifically manage water resources which can promote the ecological environment protection, and the social and economic development in Xinjiang. GRACE can monitor the global and local water storage change from space observation. With the extension of observation time and subsequent implementation of higher precision gravimetry satellites, the spatial and temporal resolution of water storage inversion will be further improved. This can provide more scientific data to optimize the allocation of water resources and promote the sustainable development of water resources in the arid area. Acknowledgements We are very grateful to the anonymous reviewers for their insightful comments and suggestions. We thank ISDC for providing GRACE data. This study is supported by the National Basic Research Program of China (grant no. 2013CB733302), the National Natural Science Foundation of China (grant nos 41374009 and 41774001), the Special Project of Basic Science and Technology of China (grant no. 2015FY310200) and the SDUST Research Fund (grant no. 2014TDJH101). REFERENCES Belda S., García-García D., Ferrándiz J.M., 2015. On the decorrelation filtering of RL05 GRACE data for global applications, Geophys. J. Int. , 200( 1), 173– 184. https://doi.org/10.1093/gji/ggu386 Google Scholar CrossRef Search ADS   Bettadpur S., 2012. UTCSR Level-2 Processing Standards Document for Level-2 Product Release 05, Rev.4.0 , Center for Space Research, University of Texas at Austin. Cao Y., Nan Z., Cheng G., 2015. GRACE gravity satellite observations of terrestrial water storage changes for drought characterization in the arid land of Northwestern China, Remote Sens. , 7( 1), 1021– 1047. https://doi.org/10.3390/rs70101021 Google Scholar CrossRef Search ADS   Chao N.F., Wang Z.T., Sun J., 2015. The inversion of terrestrial water changes by non-isotropic combination filter, Acta Geod. Cartograph. Sin. , 44( 2), 174– 182. Chen J.L., Wilson C.R., Seo K.W., 2006. Optimized smoothing of Gravity Recovery and Climate Experiment (GRACE) time-variable gravity observations, J. geophys. Res. , 111( B6), 115- 139. https://doi.org/10.1029/2005JB004064 Google Scholar CrossRef Search ADS   Chen J.L., Wilson C.R., Blankenship D., Tapley B.D., 2009. Accelerated Antarctic ice loss from satellite gravity measurements, Nature Geosci. , 2( 12), 859– 862. https://doi.org/10.1038/ngeo694 Google Scholar CrossRef Search ADS   Chen Y., Jiang Y.A., Mao W.Y., Wang S., 2011. 2009/2010 winter snow disaster in Northern Xinjiang under global warming, Adv. Clim. Chang. Res. , 7( 2), 104– 109. Cheng M., Tapley B.D., 2004. Variations in the Earth's oblateness during the past 28 years, J. geophys. Res. , 109( B9), B09402, doi:10.1029/2004JB003028. https://doi.org/10.1029/2004JB003028 Deng H., Chen Y., 2017. Influences of recent climate change and human activities on water storage variations in Central Asia, J. Hydrol. , 544, 46– 57. https://doi.org/10.1016/j.jhydrol.2016.11.006 Google Scholar CrossRef Search ADS   Duan X., Guo J., Shum C.K., van der Wal W., 2009. On the postprocessing removal of correlated errors in GRACE temporal gravity field solutions, J. Geod. , 83( 11), 1095– 1106. https://doi.org/10.1007/s00190-009-0327-0 Google Scholar CrossRef Search ADS   Fan Y., van den Dool H., 2004. Climate Prediction Center global monthly soil moisture data set at 0.5° resolution for 1948 to present, J. geophys. Res. , 109( D10), D10102, doi:10.1029/2003JD004345. https://doi.org/10.1029/2003JD004345 Google Scholar CrossRef Search ADS   Gao C.C., Lu Y., Zhang Z.Z., Shi H.L., Zhu C.D., 2015. Ice sheet mass balance in Antarctica measured by GRACE and its uncertainty, Chin. J. Geophys. , 58( 3), 780– 792. Golyandina N., Nekrutkin V., Zhigljavsky A.A., 2001. Analysis of Time Series Structure: SSA and Related Techniques , CRC Press. Google Scholar CrossRef Search ADS   Guo J., Mu D., Liu X., Yan H., Dai H., 2014. Equivalent water height extracted from GRACE gravity field model with robust independent component analysis, Acta Geophys. , 62( 4), 953– 972. https://doi.org/10.2478/s11600-014-0210-0 Google Scholar CrossRef Search ADS   Guo J., Mu D., Liu X., Yan H., Sun Z., Guo B., 2016. Water storage changes over the Tibetan Plateau revealed by GRACE mission, Acta Geophys. , 64( 2), 463– 476. https://doi.org/10.1515/acgeo-2016-0003 Google Scholar CrossRef Search ADS   Hassan A., Jin S., 2016. Water storage changes and balances in Africa observed by GRACE and hydrologic models, Geod. Geodyn. , 7( 1), 39– 49. https://doi.org/10.1016/j.geog.2016.03.002 Google Scholar CrossRef Search ADS   Hassani H., 2007. Singular spectrum analysis: methodology and comparison, J. Data Sci. , 5( 2), 239– 257. Houborg R., Rodell M., Li B., Reichle R., Zaitchik B.F., 2012. Drought indicators based on model-assimilated Gravity Recovery and Climate Experiment (GRACE) terrestrial water storage observations, Water Resour. Res. , 48( 7), 2515– 2521. https://doi.org/10.1029/2011WR011291 Google Scholar CrossRef Search ADS   Jekeli C., 1981. Alternative Methods to Smooth the Earth's Gravity Field , Department of Geodetic Science and Surveying, Ohio State University. Ju X.L., Shen Y.Z., Zhang Z.Z., 2013. Antarctic ice mass change analysis based on GRACE RL05 data, Chin. J. Geophys. , 56( 9), 2918– 2927. Ju X., Shen Y., Zhang Z., 2014. GRACE RL05-based ice mass changes in the typical regions of Antarctica from 2004 to 2012, Geod. Geodyn. , 5( 4), 57– 67. https://doi.org/10.3724/SP.J.1246.2014.04057 Google Scholar CrossRef Search ADS   Klees R., Revtova E.A., Gunter B.C., Ditmar P., Oudman E., Winsemius H.C., Savenije H.H.G., 2008. The design of an optimal filter for monthly GRACE gravity models, Geophys. J. Int. , 175( 2), 417– 432. https://doi.org/10.1111/j.1365-246X.2008.03922.x Google Scholar CrossRef Search ADS   Kusche J., 2007. Approximate decorrelation and non-isotropic smoothing of time-variable GRACE-type gravity field models, J. Geod. , 81( 11), 733– 749. https://doi.org/10.1007/s00190-007-0143-3 Google Scholar CrossRef Search ADS   Li K.M., Li Z.Q., Gao W.Y., Wang L., 2011. Recent glacial retreat and its effect on water resources in eastern Xinjiang, Chin. Sci. Bull. , 56( 33), 3596– 3604. https://doi.org/10.1007/s11434-011-4720-8 Google Scholar CrossRef Search ADS   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, Chin. J. Geophys. , 56( 6), 1843– 1849. Lu F., You W., Fan D.M., Huang Q., 2015. Chinese continental water storage and ocean water mass variations analysis in recent ten years based on GEACE RL05 data, Acta Geod. Cartograph. Sin. , 44( 2), 160– 167. Luo Z.C., Li Q., Zhong B., 2012. Water storage variations in Heihe River Basin recovered from GRACE temporal gravity field, Acta Geod. Cartograph. Sin. , 41( 5), 676– 681. Mu D.P., Guo J.Y., Sun Z.C., Kong Q.L., 2014. Equivalent water height from GRACE gravity model based on principal component analysis, Prog. Geophys. , 29( 4), 1512– 1517. Mu D., Yan H., Feng W., Peng P., 2017. GRACE leakage error correction with regularization technique: case studies in Greenland and Antarctica, Geophys. J. Int. , 208( 3), 1775– 1786. Oropeza V., Sacchi M., 2011. Simultaneous seismic data denoising and reconstruction via multichannel singular spectrum analysis, Geophysics , 76( 3), V25– V32. https://doi.org/10.1190/1.3552706 Google Scholar CrossRef Search ADS   Paulson A., Zhong S., Wahr J., 2007. Inference of mantle viscosity from GRACE and relative sea level data, Geophys. J. R. astr. Soc. , 171( 2), 497– 508. https://doi.org/10.1111/j.1365-246X.2007.03556.x Google Scholar CrossRef Search ADS   Peltier W.R., 2004. Global glacial isostasy and the surface of the ice-age Earth: the ICE-5G (VM2) model and GRACE, Annu. Rev. Earth Planet. Sci. , 32( 1), 111– 149. https://doi.org/10.1146/annurev.earth.32.082503.144359 Google Scholar CrossRef Search ADS   Plaut G., Vautard R., 1994. Spells of low-frequency oscillations and weather regimes in the Northern Hemisphere, J. Atmos. Sci. , 51( 2), 210– 236. https://doi.org/10.1175/1520-0469(1994) Google Scholar CrossRef Search ADS   Rangelova E., van der Wal W., Braun A., Sideris M.G., Wu P., 2007. Analysis of Gravity Recovery and Climate Experiment time-variable mass redistribution signals over North America by means of principal component analysis, J. geophys. Res. , 112( F3), F03002, doi:10.1029/2006JF000615. https://doi.org/10.1029/2006JF000615 Google Scholar CrossRef Search ADS   Rangelova E., van der Wal W., Sideris M.G., Wu P., 2010. Spatiotemporal analysis of the GRACE-derived mass variations in North America by means of multi-channel singular spectrum analysis, in Gravity, Geoid and Earth Observation. International Association of Geodesy Symposia , Vol. 135, no. 7, pp. 539– 546, ed. Mertikas S., Springer, Berlin, Heidelberg. https://doi.org/10.1007/978-3-642-10634-7_72 Rangelova E., Sideris M.G., Kim J.W., 2012. On the capabilities of the multi-channel singular spectrum method for extracting the main periodic and non-periodic variability from weekly GRACE data, J. Geodyn. , 54( 2), 64– 78. https://doi.org/10.1016/j.jog.2011.10.006 Google Scholar CrossRef Search ADS   Rodell M. et al.  , 2004. The global land data assimilation system, Bull. Am. Meteor. Soc. , 85( 3), 381– 394. https://doi.org/10.1175/BAMS-85-3-381 Google Scholar CrossRef Search ADS   Seoane L., Ramillien G., Frappart F., Leblanc M., 2013. Regional GRACE-based estimates of water mass variations over Australia: validation and interpretation, Hydrol. Earth Syst. Sci. , 17( 12), 4925– 4939. https://doi.org/10.5194/hess-17-4925-2013 Google Scholar CrossRef Search ADS   Shen Y., Guo J., Liu X., Kong Q., Guo L., Li W., 2018. Long-term prediction of polar motion using a combined SSA and ARMA model, J. Geod. , 2018, doi:10.1007/s00190-017-1065-3. Steffen H., Denker H., Müller J., 2008. Glacial isostatic adjustment in Fennoscandia from GRACE data and comparison with geodynamical models, J. Geodyn. , 46( 3–5), 155– 164. https://doi.org/10.1016/j.jog.2008.03.002 Google Scholar CrossRef Search ADS   Swenson S., Wahr J., 2006. Post-processing removal of correlated errors in GRACE data, Geophys. Res. Lett. , 33( 8), L08402, doi:10.1029/2005GL025285. https://doi.org/10.1029/2005GL025285 Google Scholar CrossRef Search ADS   Tamisiea M.E., 2011. Ongoing glacial isostatic contributions to observations of sea level change, Geophys. J. Int. , 186( 3), 1036– 1044. https://doi.org/10.1111/j.1365-246X.2011.05116.x Google Scholar CrossRef Search ADS   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, J. geophys. Res. , 103( B12), 30205– 30229. https://doi.org/10.1029/98JB02844 Google Scholar CrossRef Search ADS   Wang H., Xiang L., Jia L., Wu P., Steffen H., Jiang L., Shen Q., 2015. Water storage changes in North America retrieved from GRACE gravity and GPS data, Geod. Geodyn. , 6( 4), 267– 273. https://doi.org/10.1016/j.geog.2015.07.002 Google Scholar CrossRef Search ADS   Wang X., de Linage C., Famiglietti J., Zender C.S., 2011. Gravity Recovery and Climate Experiment (GRACE) detection of water storage changes in the Three Gorges Reservoir of China and comparison with in situ measurements, Water Resour. Res. , 47( 12), W12502, doi:10.1029/2011WR010534. https://doi.org/10.1029/2011WR010534 Google Scholar CrossRef Search ADS   Zhang J.B., 2001. Strategic countermeasures for sustainable use of water resources of Xinjiang, Arid Land Geogr. , 24( 3), 217– 222. Zotov L.V., Shum C.K., 2010. Multichannel singular spectrum analysis of the gravity field data from GRACE satellite, in Proceedings of the 4th Gamow International Conference on Astrophysics and Cosmology after Gamow and the 9th Gamow Summer School “Astronomy and Beyond: Astrophysics, Cosmology, Radio Astronomy, High Energy Physics and Astrobiology”, Ukraine , Vol. 1206, no. 1, pp, 473– 479, eds Chakrabarti S.K., Bisnovatyi-Kogan G.S., Zhuk A.I., AIP Publishing. Zotov L.V., 2012. Application of multichannel singular spectrum analysis to geophysical fields and astronomical images, Adv. Astron. Space Phys. , 1( 2), 82– 84. © The Author(s) 2018. Published by Oxford University Press on behalf of The Royal Astronomical Society.

### Journal

Geophysical Journal InternationalOxford University Press

Published: Apr 1, 2018

## You’re reading a free preview. Subscribe to read the entire article.

### DeepDyve is your personal research library

It’s your single place to instantly
that matters to you.

over 18 million articles from more than
15,000 peer-reviewed journals.

All for just $49/month ### Explore the DeepDyve Library ### Search Query the DeepDyve database, plus search all of PubMed and Google Scholar seamlessly ### Organize Save any article or search result from DeepDyve, PubMed, and Google Scholar... all in one place. ### Access Get unlimited, online access to over 18 million full-text articles from more than 15,000 scientific journals. ### Your journals are on DeepDyve Read from thousands of the leading scholarly journals from SpringerNature, Elsevier, Wiley-Blackwell, Oxford University Press and more. All the latest content is available, no embargo periods. DeepDyve ### Freelancer DeepDyve ### Pro Price FREE$49/month
\$360/year

Save searches from
PubMed

Create lists to

Export lists, citations