Analysis of changes in drought and terrestrial water storage in the Tarim River Basin based on principal component analysis

Climate change may be accelerating the hydrological cycle in many regions of the world (Alan et al. 2003; Zhang et al. 2013a), and increases in extreme precipitation, persistent dry conditions, floods, and tropical cyclones have significantly increased the risk of extreme weather and climatic events since the 1950s (Mirza 2002; Zhang et al. 2013b; Apurv et al. 2015). In this respect, hydro-meteorological hazards (such as floods and droughts), which result in death and the loss of property and are usually driven by climatic extremes (Li et al. 2015), are attracting considerable research attention (e.g., Beniston & Stephenson 2004; She et al. 2015).

Drought is also a widespread environmental hazard that is currently the focus of a considerable amount of research (Li et al. 2016b). The increased demand for water and the oscillation of hydro-meteorological fluxes in relation to climate change have intensified the impacts of drought (Mishra & Singh 2010). Many previous studies have been conducted on the multiple aspects of drought events, such as drought forecasting (Mishra & Desai 2005), drought frequency analysis (Mishra & Singh 2009), and drought spatiotemporal analysis (Mishra & Singh 2010). In this respect, meteorological drought, agricultural drought, and hydrological drought can be characterized using standardized indices that represent a scaled continuum of extreme wetness to extreme dryness (McKee et al. 1993), whereas agricultural drought can also be expressed using the Palmer drought severity index (PDSI), which describes the level of soil moisture deficit (Palmer 1968).

The main sources of water in the Tarim River Basin (TRB) are glaciers and snow, and the water system within the basin is thus extremely sensitive to climate change (Chen et al. 2007). With an increase in the area of upstream arable land, there has been an associated increase in the amount of water consumed; this has caused a corresponding decrease in runoff and drying of the river downstream (Chen et al. 2014). As a result, agricultural areas are competing for water with the needs of the environment (Chen et al. 2016). Additionally, an extensive amount of desert vegetation has been lost due to expanding artificial oases, which has created a lifeless desert (Chen et al. 2016) that is characterized by reduced biodiversity and encroaching desertification. Although the central government of China spent approximately ¥10.7 billion (US$1.3 billion) on river reclamation and Taitemahu Lake in 2000 (Tao et al. 2008), the program was not significantly effective. However, the government has also proposed a ‘one belt, one road’ policy to strengthen the development of western China and to increase its economic connection with Europe and Africa. Hence, detection of the spatiotemporal changes associated with drought and their potential factors in the TRB (i.e., teleconnection indices) are significant to determine for its core status within the new policy.

Although drought in the TRB has been previously studied (e.g., Zou et al. 2005; Zhang et al. 2012), this study aimed to detect spatiotemporal changes in drought and their potential associated factors within the TRB by examining an adapted SPI, standardized soil moisture index, and terrestrial water storage (TWS) data obtained from the Gravity Recovery and Climate Experiment (GRACE).

The TRB is located in the southern part of Xinjiang Province, China, and is bordered by the Tien Shan Mountains, Pamir Plateau, and Kunlun Mountains (Li et al. 2016b). The Tarim River (TR) is the longest and largest endorheic river in northwest China, and it is the major water source of the TRB (Li et al. 2016b) (Figure 1). Light precipitation and strong potential evaporation occur in the basin, and an extremely arid climate is widespread throughout the TRB (Sun et al. 2012). With the exception of alpine rainfall, seasonal snowpack and glacier melt dominate the runoff in the TR. The monthly temperature ranges from 20 °C to 30 °C in July and from −20 °C to 10 °C in January (Chen et al. 2009).

Approximately four to six billion cubic meters of annual runoff is generated from glacier/snow melt in the TRB, which accounts for 48.2% of the total water volume in the basin (Chen et al. 2007). Approximately ten million people live in the basin, and eight million of these live around oases or in the middle and lower river reaches (Li et al. 2016a). The population of the three largest headstreams (i.e., Aksu River, Hotan River, and Yerkand River) accounts for 73.2%, 23.2%, and 3.6% of the total, respectively (Chen et al. 2007). Oases in the upper and lower reaches of the TR were enlarged steadily and significantly for agricultural activities during the 1950s and 1990s (Li et al. 2016b).

The observed meteorological site data in the TRB relate only to a limited time period. Therefore, the Climatic Research Unit (CRU) TS3.24 dataset (http://www.cru.uea.ac.uk/cru/data/hrg/) was applied to explore spatial and temporal variations in precipitation. The spatial resolution for the monthly precipitation product is 0.5° × 0.5°, and SPI was obtained from the reanalysis precipitation product in this study.

Monthly Tropical Precipitation Measuring Mission (TRMM)-V73B43 precipitation rates recorded in 1998–2015 were downloaded from https://trmm.gsfc.nasa.gov/ and used to verify the reanalysis precipitation product in the TRB. These data provide precipitation products that are useful for understanding the distribution and variability in precipitation between 50°N and 50°S. The TRMM3B43 provides monthly precipitation estimates with a high spatial resolution of 0.25° × 0.25° (Kummerow et al. 2000; Huffman et al. 2007).

In this study, monthly Climate Prediction Center (CPC) soil moisture data with a spatial resolution of 0.5° × 0.5° were analyzed to explore water stress based on the SSI. In addition, CPC soil moisture data were assessed to detect climate teleconnection (e.g., El Niño–Southern Oscillation (ENSO)) modes. The data were acquired from the National Oceanic and Atmospheric Administration (NOAA) at http://www.esrl.noaa.gov.

GRACE data beginning in April 2002 were downloaded from http://grace.jpl.nasa.gov/. TWS was computed using Level-2 RL05 data from GRACE spherical harmonics. In this respect, the original C20 coefficients were first replaced with new coefficients from satellite laser ranging (SLR), and the degree-1 coefficients were then identified from the estimation. The modified spherical harmonics coefficients were subsequently moved by Gaussian filtering with a 300-km smooth radius, and the anomaly coefficients were obtained by removing the monthly average coefficients during January 2004 to December 2009 from all the time series. Finally, the anomaly coefficients were converted into TWS based on the method proposed by Wahr et al. (1998). Missing data during the research period were obtained by linear interpolation.

Numerous researchers have found that precipitation and drought conditions in northwest China are associated with climate indices (Wang et al. 2015b; Li et al. 2016a) relating to the North Atlantic Oscillation (NAO), the Arctic Oscillation (AO), ENSO (i.e., El Niño 3.4), and Pacific Decadal Oscillation (PDO). Thus, monthly anomalies of NAO, AO, ENSO, and PDO from the Earth System Research Laboratory (ESRL; http://www.esrl.noaa.gov) recorded in 2002–2015 were selected to analyze the potential factors involved in water storage variations within the TRB using the linear correlation coefficient.

The traditional method of calculating the SPI cannot be applied effectively in areas with zero-precipitation months as this impacts the accumulated probability (Farahmand & AghaKouchak 2015; Vergni et al. 2016). However, a novel method can be used to solve the problem by employing a nonparametric technique that can calculate both the precipitation deficit and the soil moisture deficit (Farahmand & AghaKouchak 2015). To study wet and dry conditions in the TRB, this study used the SPI and the SSI, and a novel and nonparametric method following that of Farahmand & AghaKouchak (2015) was used to obtain the marginal probability of precipitation and soil moisture based on the empirical Gringorten plotting position.

The EOF can be defined by measuring the spatial domain (Westra et al. 2010). If the expression is for the sub-regions only, the variance from the first four patterns can be improved (Ndehedehe et al. 2016; Yang et al. 2017). The positive or negative features of the product between the EOFs and their corresponding principal components (PCs) represent the trends in the grids.

The time series for monthly precipitation and precipitation cycles based on the CRU precipitation product are reported in Figure 2. TRMM precipitation was consistent with CRU precipitation in the TRB from 1998 to 2015 with a determination coefficient of 0.71. Despite the significant increasing trend in spatially averaged precipitation, dry conditions occurred in 1980, 1986, 1994, and 2009 (Figure 2(a)). The annual cycle of precipitation in the TRB peaked in July, whereas the seasonality of Niño 3.4 exhibited the opposite phase to precipitation (Figure 2(b)); these results are similar to those reported in previous studies (Wang et al. 2015a, 2015b).

The variation in seasonal precipitation was also analyzed (Figure 3), and as shown in Figure 3(c) and 3(d), seasonal maximum precipitation occurred during summer (June–July–August) with an average precipitation of approximately 65 mm, whereas less than 50 mm of precipitation occurred during the other seasons. However, the spatial distribution of seasonal precipitation showed a maximum of approximately 170 mm in the western part of the TRB during spring (March–April–May) (Figure 3(a)). Comparatively, maximum amounts of precipitation measuring approximately 120 mm were detected during the summer in the northern and southern parts of the TRB. A significant discrepancy between the edge sections of the TRB and the central part of the TRB was noted (Figure 3).

To analyze the change points (i.e., the turning points or discontinuity points between two adjacent distinct trends) of precipitation in the TRB from 1948 to 2015, the monthly and annual cumulative precipitation values were then evaluated and plotted (Figure 4). A strong precipitation deficit occurred during the 1980s, and a meteorological drought became an agricultural and hydrological drought that affected water storage (surface water, soil moisture, and groundwater). Previous drought studies have reported hydrological drought relating to precipitation deficits at accumulation periods of 6, 12, and 24 months (Hayes et al. 1999; Vicente-Serrano 2006; Santos et al. 2010; Li & Rodell 2015). Therefore, CRU precipitation and CPC soil moisture data in the TRB from 1948 to 2015 were applied to evaluate the drought frequency and severity of drought events at different accumulation periods.

Following the method proposed by Farahmand & AghaKouchak (2015), the CRU precipitation and CPC soil moisture data in the TRB were converted to standardized drought indices (i.e., SPI and SSI) at 3-, 6-, 12-, and 24-month accumulation periods (Figure 5). Results showed that the SPI and SSI series were consistent with each other from 1948 to 2015. However, decreasing trends (p < 0.05) were detected in the SPI and SSI time series from 1948 to 1986, whereas increasing trends (p < 0.05) were noted from 1987 to 2015. Several drought events occurred in 1957–1959, 1961–1963, 1978–1980, 1984–1987, 1995–1997, 2001–2002, and 2008–2010 in all accumulation periods (Figure 5), and a number of drought events in the TRB have also been detected and reported in previous research (Chen et al. 2014; Ye et al. 2014). The levels of some of the drought events produced harmful effects on agrarian systems and food security in the TRB, where agriculture is dominated by oasis agriculture.

Figure 6 shows the temporal and spatial patterns of SPI3. The first four patterns of EOF analysis account for 43.9%, 8.9%, 6%, and 5.6% of the original signal, respectively. Furthermore, EOF1 shows a uniform negative feature (Figure 6(a)) with its center in the TRB, PC1 (Figure 6(b)) exhibits a decreasing trend, and the product between EOF1 and PC1 reports a wetness feature in the region. The second EOF and its PC express annual variations in SPI3 in some sub-regions of the TRB from 1948 to 2015, and although the valleys of PC2 in the positive EOF sub-regions show dry periods, the peaks in the positive EOF regions relate to wet periods. The third EOF (Figure 6(e)) and its PC (Figure 6(f)) show variations at a three-month accumulation period, where a positive center is detected in the Taklimakan Desert and a negative center in the northeastern part of the TRB. Finally, the fourth PCA mode shows insignificant changes in significant seasonal precipitation modes in the TRB.

Figure 7 exhibits the temporal and spatial characteristics of SPI6 in the TRB. Although the first EOF maintains 46.89% of the original signal, spatial discrepancies within the TRB are clearly evident. As the uniformity value in EOF1 is negative (Figure 7(a)), the peaks in the curve of PC1 (Figure 7(b)) (e.g., 1962, 1971, 1981, 1986, 1994, 1998, 2007, and 2009) represent several dry years, and the valleys (i.e., 1963, 1973, 1988, 1993, 1996, 2004, and 2010) denote typical wet years. With respect to the anti-phase distribution of EOF2, the valleys of PC2 in the positive EOF sub-regions represent dry periods, whereas the peaks of PC2 in the positive EOF regions represent wet periods (Figure 7(c) and 7(d)). EOF3 also shows an anti-phase distribution and accounts for 5.92% of SPI6, and the valleys and peaks of PC3 in the positive EOF sub-regions are both dry years. Furthermore, an anti-phase distribution is detected in EOF4, and the valleys and peaks of PC4 in the positive EOF sub-regions are also both dry years.

The temporal and spatial features of SPI12 in the TRB are shown in Figure 8. The first four patterns of EOF analysis represent 46.46%, 9.81%, 6.51%, and 6.09% of the original signal, respectively. The first component (Figure 8(a) and 8(b)) explains the multi-annual variability of SPI, but a negative feature (with its center on the largest negative value) sweeps over the TRB. Thus, when PC1 (Figure 8(b)) exhibits a decreasing trend, it represents a wetness feature in the region. Additionally, for EOF2, a positive center is evident in the eastern part of the TRB (Figure 8(c)) but negative values are centered in the western part. Furthermore, PC2 (Figure 8(d)) exhibits an increasing trend. The valleys of PC2 in the positive EOF sub-regions are dry periods, and the valleys in the negative EOF regions are wet periods. The third EOF (Figure 8(e)) and its corresponding PC (Figure 8(f)) show variation on a monthly scale. A positive center is detected in the Taklimakan Desert, and a negative center can be identified in the northeastern part of the TRB. Finally, EOF4 accounts for 6.09% of SPI12; an anti-phase distribution is detected, and thus the valleys and peaks of PC4 in the positive EOF sub-regions both represent dry years.

A comparison between the EOFs of SPI3, SPI6, and SPI12 reveals analogous features (Figures 6–8). However, the drought frequency shown in Figure 8 is reduced compared with that shown in Figures 6 and 7, where drought frequencies and durations are greater than those observed for the short accumulation period. Therefore, longer drought durations are not linked to the aggregation window. The aggregated SPI at longer accumulation periods (e.g., 12 and 24 months) is usually applied to investigate hydrological droughts (e.g., Santos et al. 2010; Lloyd-Hughes 2012; Ndehedehe et al. 2016). However, the conversion process from meteorological drought to hydrological drought based on surface soil conditions requires the use of data obtained over lengthy periods (Figure 5).

Soil water drought and hydrological drought usually both occur during most drought events (Ndehedehe et al. 2016). Many previous studies have analyzed hydrological and agricultural drought using SPI at multiple accumulation periods (Hayes et al. 1999; Lloyd-Hughes 2012; Li & Rodell 2015). However, the usefulness of quantifying drought through soil moisture in hydrological applications remains unknown (Ndehedehe et al. 2016). In this study, therefore, we converted the CPC soil moisture product to SSI at accumulation periods similar to those of the SPI based on the method proposed by Farahmand & AghaKouchak (2015). The EOF was then applied to extract the PCs of SSI in the grids over the TRB (Figures 9–11).

The PCA results of SSI3 are shown in Figure 9. As a negative EOF1 (Figure 9(a)) is exhibited, peaks in the PC1 curve (Figure 9(b)) represent dry periods (e.g., 1956–1957, 1962, 1971, 1978–1979, 1983–1986, 1994, 2001, and 2008–2009) and valleys in the PC1 curve represent wet periods (e.g., 1950, 1953–1954, 1972–1973, 1982–1983, 1988–1992, 2003–2005, and 2010). The second EOF, which accounts for 10.6% of the original variations, shows an anti-phase pattern in the sub-regions, where the positive center is located in the eastern part of the TRB and the negative kernel remains in the western part. Because of this anti-phase distribution, the peaks of PC2 in the positive EOF sub-regions represent dry periods, and the valleys of the PC2 in the positive EOF regions represent wet periods. For EOF3, which accounts for 6.33% of the original variation, the positive kernel occurs in the southeastern part of the TRB, and the negative center in the northeastern part. For the fourth EOF, which accounts for 5.67% of the variation, the obvious kernel is located in the northern part of the TRB.

Figure 10 shows the results of PCA for SSI6. Compared with the PCs of SSI3, those of SSI6 exhibit a lower frequency and longer duration of drought. As negative EOF1 (Figure 10(a)) is found over the TRB, the significant peaks in the PC1 curve (Figure 10(b)) represent dry periods (e.g., 1957, 1962, 1971, 1978–1979, 1984–1987, 1994, 2002, and 2008–2009) and the valleys in the PC1 curve represent wet periods (e.g., 1950, 1953–1954, 1972, 1982–1983, 1987–1993, and 2002–2006). EOF1 accounts for 51.27% of the original variation. The second EOF shows an anti-phase pattern in the sub-regions of the TRB (Figure 10(c)); in relation to this anti-phase distribution, the peaks of PC2 in the positive and negative regions of EOF2 represent wet and dry periods, respectively (Figure 10(d)). For EOF3, the positive and negative centers are located in the southeastern and northeastern parts of the TRB, respectively, and the peaks of PC3 in the positive and negative regions of EOF3 are wet and dry periods (Figure 10(f)). For the fourth EOF (Figure 10(g)), which accounts for 5.27% of the original signal, the obvious center is located in the northern part of the TRB.

Figure 11 shows the PCA results for SSI12. The first four EOFs of SSI12 account for 50.8%, 11.36%, 7.05%, and 5.52% of the original data. As a negative EOF1 (Figure 11(a)) is detected in the TRB, the significant peaks and valleys in PC1 curve (Figure 11(b)) represent dry and wet periods, respectively. The second EOF shows an anti-phase pattern in the sub-regions of the TRB (Figure 11(c)). Hence, although the peaks of PC2 in the positive regions of EOF2 are wet periods, the peaks of PC2 in the negative regions represent dry periods (Figure 11(d)). For EOF3 (Figure 11(e)), the positive (negative) core is located in the southeastern (northeastern) part of the TRB. In addition, the peaks of PC3 in the positive regions of EOF3 are wet periods and those in the negative regions are dry periods (Figure 11(f)). However, for the fourth EOF, the obvious center is located in the northern part of the TRB (Figure 11(g)).

Similar patterns in EOFs are noted for SSIs at different periods (i.e., SSI3, SSI6, and SSI12). In addition, similar wet or dry periods are found in all monthly periods in the TRB (Figures 9–11). However, several inconsistencies are found in the details between the temporal evolutions of SSI and SPI in the TRB from 1948 to 2015, and the cause needs to be identified in future work. Nevertheless, extreme soil moisture events always follow extreme precipitation events but with a lag of several months, which is consistent with the development of hydrological processes (e.g., Figure 5). Additionally, several limitations are noted in the ability of the soil moisture products modeled by observation data to reflect actual wet or dry conditions, as previously proposed by Ndehedehe et al. (2016) and Dirmeyer et al. (2004). In addition, Dirmeyer et al. (2004) noted that changes in precipitation were not always followed by changes in soil moisture products for the surface characteristics of the soil.

To explore the potential factors of changes in the SPI and SSI, the relationship between SPI/SSI and the teleconnection indices in the TRB from 1948 to 2015 was analyzed (Figure 12 and Table 2). As shown in Table 2, the two largest correlation coefficients (i.e., −0.138 and 0.157) are found between PC1 of SPI in the 6-month period and PDO and between PC1 of SPI in the 12-month period and Niño 3.4. Hence, compared with AO and NAO, the PDO is considered to be the major driver of precipitation for the 6-month period in the TRB from 1948 to 2015, although ENSO/Niño 3.4 had a more significant influence on precipitation in the 12-month period. For soil moisture, correlation coefficients at a 95% confidence level are found between PC1 of the SSIs (i.e., SSI3, SSI6, and SSI12) and Niño 3.4 in the TRB from 1948 to 2015. The three significant correlation coefficients for the different timescales of SSI (i.e., SSI3, SSI6, and SSI12) and Niño 3.4 are 0.239, 0.252, and 0.283, respectively.

Bold numbers represent correlation coefficients within the 95% significance level (p < 0.05).

To investigate these findings further (in addition to precipitation and soil moisture) TWS was also analyzed by EOF (e.g., Figure 13). The first four patterns of accumulative variance, at 71.4%, 17.7%, 5.9%, and 2.9%, were found to account for 97.9% of the total variance (Yang et al. 2017). As EOF1 displayed a unified increase feature, it denoted a more obvious similar feature in the TRB from 2002 to 2015. Hence, the peaks in the curve of PC1 (Figure 13(b)) stand for wet periods (i.e., 2005–2007 and 2010–2011) and the valleys in the curve represent dry periods (i.e., 2002–2004 and 2008–2009). Therefore, the wet and dry periods of TWS are similar to those in the precipitation series and soil moisture series. The second EOF (Figure 13(c)) shows anti-phase changes in TWS in some sub-regions of the TRB. Although the peaks in the curve of PC2 (Figure 13(d)) represent wet periods for the positive sub-regions, the valleys in the curve represent dry periods for the negative sub-regions. The third EOF (Figure 13(e)) and its corresponding PC (Figure 13(f)) show variations on a monthly scale, and both negative and positive centers are found in the study area. The smallest major PCA pattern, which expresses approximately 1.5% of the total TWS variability, denotes an unobvious signal. Overall, regardless of the TWS deficit or surplus, these values are consistent with those of precipitation variations (i.e., deficit or surplus) and soil moisture changes (i.e., deficit or surplus) in the TRB from 2002 to 2015.

Salami et al. (2015) proposed that climate and environmental conditions affect natural water bodies such as rivers and lakes and can thus be used to explore water availability. Hence, the water levels of Bosten Lake in 2002–2012 were analyzed to identify changes in precipitation, soil moisture, and TWS in the TRB. As shown in Figure 14, the water level in Bosten Lake reduced at a rate of −2.8 cm/month from 2002 to 2012. In addition, a significantly lower water level occurred in late 2009 and during the first few months of 2010. This trend in the water level during 2009–2010 was similar to that of soil moisture and GRACE TWS in 2002–2012, even though precipitation showed an increase. However, as noted in previous studies (e.g., Zhang et al. 2015), population growth, shortage of water supplies, and overuse of water resources could account for the changes in TWS and in the water level.

In this study, the first four temporal and spatial modes of SPI, SSI, and TWS were obtained in the TRB using EOF analysis. In addition, teleconnection indices were applied to identify the relationship between SPI/SSI temporal changes and atmospheric circulation. The major conclusions of this study are discussed below.

Although there was a lag time of several months between them, the SPIs and SSIs were relatively consistent in the TRB from 1948 to 2015. Significant similar patterns in the spatial distributions of SPIs and SSIs indicate that the soil moisture in the TRB was significantly affected by precipitation, although actual evapotranspiration was high in the basin. However, there were no apparent identical features between the EOFs of SSIs and those of TWS, which implies that the effect of soil water on TWS was less than that of groundwater. The water level of Bosten Lake was in line with the values of SPIs, SSIs, and TWS from 2008 to 2010. As precipitation, soil moisture, and TWS experienced simultaneously severe water deficits which then recovered, the conditions imply that hydrological fluxes were closely related to each other during extreme hydro-meteorological events in the TRB.

Drought durations were longer on a longer timescale but drought was more frequent on a shorter time scale, according to drought evaluation based on SPI and SSI. Additionally, the climate teleconnection of ENSO had a significant positive relationship with SSI on 3-, 6-, and 12-month scales for PC1 in the TRB from 1950 to 2015, but no significant connection was noted between the other three climatic indices (i.e., AO, NAO, and PDO) and the SPI at all scales for PC1 and PC2 in the TRB in 1950–2015. The results of this study are useful for drought detection and water management related to anthropological activities.

The research is supported by the National Key R&D Program of China (2017YFA0603702). We are also grateful to Professor Xia Jun and Zhan Chesheng for their assistance in helping the second author to finish the main manuscript text, and to Professor Qiao Yunfeng for assisting with figure preparation.