Change of precipitation characteristics in the water-wind erosion crisscross region on the Loess Plateau, China, from 1958 to 2015

Precipitation plays an important and crucial role in processes in the water-wind erosion crisscross region of the Loess Plateau than in other parts of the region. We analyzed precipitation data and standardized precipitation index (SPI) at 14 representative synoptic stations from 1958 to 2015 used trend-free prewhitening, linear trend estimation, Spearman’s rho test, the Mann-Kendall trend test, the Mann-Kendall abrupt change test and rescaled range analysis. The following conclusions were drawn. First, the analysis of monthly precipitation at all stations suggested that precipitation during the rainy season (July, August, September), especially rain in July and August, exhibited a general decreasing trend, while both increasing and decreasing trends were observed in other months. Moreover, the annual precipitation of all stations continued to exhibit decreasing trends except Wuzhai. Erosive rainfall frequency in the rainy season and the annual scale was weakly reduced but erosive force of single rainfall has been enhanced. Second, the SPI exhibited different increasing degrees in winter, while decreasing trends were observed in other seasons. Additionally, the annual-scale SPI at most stations exhibited a stable and sustained downward trend. Therefore, this region is currently associated with a drought trend, and the drought degree will likely continue to increase.


Results
Analysis of precipitation characteristics. Time series of annual precipitation collected at 14 synoptic stations in the study area from 1958 to 2015 (Fig. 1) and annual precipitation characteristics from each station ( Table 1) show that precipitation exhibited a skewed distribution during the study period. The spatial differences in precipitation were obvious in the study area. Notably, the average annual precipitation at the 14 synoptic stations ranged from 234.59 (Jingyuan) to 493.28 mm (Xingxian), and the highest and lowest amounts of annual precipitation were 844.4 mm (Xingxian station in 1964) and 114.9 mm (Youyu station in 1962), respectively. The interannual variations in annual precipitation at 14 synoptic stations were also significant. The coefficients of variation (CVs) of annual precipitation were higher than 20% at the synoptic stations, and the highest value occurred at Hequ (31.00%), while the lowest value was 21.57% at Youyu.
The mean monthly precipitation and mean monthly air temperature data from 14 synoptic stations in the study area are shown in Fig. 2. Figure 2 shows that precipitation is mainly concentrated in July-September in the study area, with rain and high temperatures in the same period. The total amount of precipitation in these three months can account for more than 60% of the total annual precipitation; thus, this period is considered the rainy season in the study area.

Precipitation trend analysis
Serial correlation coefficient and trend-free prewhitening. In this study, formula (8) is used to calculate the autocorrelation coefficient of precipitation data at monthly, seasonal and annual scales. The results are shown in Table 2. According to formula (10), we conclude that lag-1 ∈ (−0.2748, 0.2398) (n = 58). Thus, the null hypothesis that the time series data show no autocorrelation was accepted. The data in Table 2 show that very Scientific REPORTS | 7: 8048 | DOI: 10.1038/s41598-017-08600-y few stations exhibit autocorrelation at the monthly, seasonal and annual time scales, and the TFPW method is used to analyze the precipitation time series with autocorrelation (dark area in the table). The values fall between −0.2748 and 0.2398; therefore, autocorrelation is eliminated from the original precipitation data. The next step is performing linear trend estimation, Spearman's rho test and the M-K tests 29 . Table 3 presents the trend test for monthly precipitation at the 14 synoptic stations in the study area from 1958 to 2015. The information in the table shows that the trend of precipitation at each station differs in different months. Overall, the number of synoptic stations that exhibited increasing trends or decreasing trends in January, May and December was relatively constant. Additionally, the majority of stations exhibited increasing trends in February and June. Only Hequ and Xingxian stations exhibited decreasing trends in February, and only Suide exhibited a decreasing trend in June. The stations exhibited downward trends in the remaining months. Notably, precipitation declined in July, August and September (rainy season) at most stations and in July and August at all stations. Among the stations, precipitation at Jingyuan, Yuzhong and Xiji stations in August decreased significantly by 0.978 mm/year, 0.646 mm/year and 0.862 mm/year, respectively. However, only a few months of precipitation trends were significant at the 95% level. All of the significant changes in the rainy season were decreasing trends, and the trends in other months were significant but weak increases. However, Xingxian exhibited a significant decreasing trend in October. The significant monthly precipitation trend at each synoptic station is shown in Fig. 3.

Trend analysis of precipitation at different time scales.
The results of the trend analysis of seasonal and annual precipitation at each synoptic station are shown in Table 4. According to the information in the table, the variations in precipitation in different seasons differ at different synoptic stations. Notably, two trends can be observed in the spring, autumn and winter at each synoptic station, and the only summer trend is a decreasing trend. In spring, the number of synoptic stations that exhibited each type of trend was relatively equal, while more stations exhibited decreasing trends in autumn, and more stations exhibited increasing trends in winter. All significant trends occurred in winter (Haiyuan station, Guyuan station and Xiji station) and were increasing trends with rates of 0.082 mm/year (Haiyuan), 0.071 mm/ year (Guyuan) and 0.058 mm/year (Xiji). The annual precipitation trends indicate that annual precipitation at almost all stations exhibited a decreasing trend (only Wuzhai station exhibited a non-significant increase trend).  Although there no significant increasing or decreasing trend was observed, the results suggest that annual precipitation is associated with a decreasing trend in the water-wind erosion crisscross region on the Loess Plateau.
Trend analysis of erosive rainfall. Daily rainfall greater than 12 mm has been used as the standard for describing erosive rainfall on the Loess Plateau [24][25][26] . In the study area, the erosive rainfall occurred in the rainy season of each synoptic station was more frequent than non-rainy season (Fig. 4). The change trend of frequency of erosive rainfall (rainy season, non-rainy season and annual) was analyzed by linear trend estimation, and the results showed that there were differences in the frequency of occurrence of erosive rainfall at different timescales, for specific performance: it presented a decrease trend in rainy season and a increase trend in non-rainy season, however, most of synoptic stations also showed a decrease trend but the trend was weak. In terms of precipitation depth, the relationship between average annual erosive rainfall (x) and annual precipitation (y) in the study area were analyzed, and the result showed that there was a not quite obvious linear relationship between these two parameters from 1958 to 2015 (y = 0.568x-60.01, R = 0.684, P < 0.001). However, the average annual precipitation in the whole study area showed a decreasing trend but the erosive rainfall was gradually increased (Fig. 5).
The frequency of erosive rainfall in the rainy season and the annual scale was reduced but the erosive precipitation depth in the corresponding period was increased, this indicated the erosive force of single rainfall in the study area has been enhanced that during the study period. SPI analysis. Over a timescale of 1-6 months, the SPI can reflect agricultural drought conditions, and over a longer timescale (6 months or longer), it can be used to research hydrological drought in groundwater, rivers, and reservoirs 3, 6, 7 . In the SPI trend analyses at seasonal and annual scales, the SPI-3 (3-month scale) values in May, August, November and February represent spring, summer, autumn and winter, respectively, and SP-12 (12-month scale) values in December represent the whole year. The trend analyses of seasonal and annual SPI (Table 5) showed that the seasonal SPI differed at each station. Specifically, in the spring, summer and autumn, the SPI exhibited a decreasing trend, while the SPI generally increased in winter. Thus, spring, summer, autumn are becoming more arid, and the winter is becoming more humid. Annual SPI exhibited a decreasing trend at all stations except Wuzhai and Datong. Although the trends are not significant, the SPI values in the study area are trending toward drought conditions, and the occurrence of drought will likely continue to increase. The M-K abrupt analysis of annual SPI at each station ( Fig. 6) showed that the SPI at Wuzhai station decreased from 1960 to 2012. The SPI then abruptly increased and continued to increase until 2015. The annual SPI in Datong fluctuated, including increases from 1990-2000 and after 2013. The UF k curves of other stations are always below the y = 0 line, except those of Wuzhai and Datong. Thus, the annual SPI trend is continuously and stably decreasing at most stations, but the annual SPI trends at Wuzhai and Datong are unstably increasing.
The M-K trend test and Spearman's rho test were mainly used to analyze the current trends of time series, whereas the R/S analysis was mainly performed to evaluate the long-term correlation of time series 30,31 . The combination of the three methods can not only clarify the current trends of time series but also predict future trends to a certain extent. As illustrated by the Hurst index values of the annual SPI at each station (Fig. 7), only Xingxian (H = 0.485), Youyu (H = 0.398) and Datong (H = 0.380) exhibit trends that are predicted to change, and other stations will likely maintain current trends. The annual SPI at Xingxian and Youyu may increase, which is different from the current trend, while the annual SPI at Datong is likely to decrease in the future. Based on Figs 6 and 7, although the annual SPIs of the above three stations may exhibit different trends in the future, most of the stations still exhibit stable and sustained downward trends. Therefore, the degree of drought will likely continue to increase in the future.   Table 2. Lag-1 serial autocorrelation coefficient of precipitation at different time scales. Bold indicates that autocorrelations in the precipitation time series were corrected.
SPI is calculated based on monthly precipitation. In this study, a cluster analysis based on the Ward method was performed using monthly precipitation data collected at 14 synoptic stations from 1958-2015. According to the results of the cluster analysis ( Fig. 8), three stations from different cluster types but with similar values of annual precipitation were chosen to analyze SPI-1 (1-month scale), SPI-3 and SPI-12 trends: Wuqi (470.25 mm), Wuzhai (474.15 mm) and Guyuan (455.54 mm). As shown in Fig. 9, both monthly scale SPI-1 and seasonal scale SPI-3 are sensitive to short-term precipitation changes. Additionally, these values fluctuate considerably and fully reflected the short-term characteristics of the study area with respect to impermanent droughts and frequent variations. Moreover, the annual SPI-12 shows that the periodicity of drought and flood fluctuation is more obvious and stable.
The value of SPI-12 in December can be used to determine the dry and wet conditions in the subsequent year. The frequencies of dry, wet and normal years at each station are shown in Fig. 10. The figure shows that the

Discussion and Conclusions
Precipitation is closely related to hydraulic erosion and drought, but drought often causes severe wind erosion. Compared with other regions of the Loess Plateau, precipitation plays a more important role in the water-wind erosion crisscross region. Additionally, precipitation is the main form of soil water recharge on the Loess Plateau, and soil water is the main limiting factor of vegetation restoration. Fully understanding the characteristics of precipitation and dry/wet conditions in the water-wind erosion crisscross region can improve the management of regional water resources, the water use efficiency, the efficiency of vegetation restoration, and the ecological environment in the region. To understand the characteristics of precipitation and dry/wet conditions, this study analyzed these characteristics and variations in the SPI at different time scales based on observations from 14 synoptic stations in the study area. The analysis included trend-free prewhitening, linear trend estimation, Spearman's rho test, the M-K trend test, the M-K abrupt change test and rescaled range analysis. The variations in precipitation differed at different time scales. Notably, monthly precipitation in different months varied at different stations. These differences may be related to the fractured topography of the study area 32 , but rainy season precipitation, especially that in July and August, generally exhibited a decreasing trend. Winter precipitation at most stations in the study area displayed an increasing trend and most stations exhibited decreasing trends in the spring, summer autumn. Specifically, the annual precipitation at Wuzhai station began to shift to a non-significant increasing trend, and the rest of the stations exhibited decreasing trends. The water-wind erosion crisscross region on the Loess Plateau is located in the monsoon region of China. Previous studies showed that the monsoon climate in China is influenced by the East Asian summer monsoon 33,34 , and the East Asian summer monsoon is strongest in July and then begins to weaken until it ends 35 . Recent studies have shown that the East Asian summer monsoon is weakening 36,37 , which may be one of the main reasons for the decreases in precipitation in July and August and in summer and autumn. Additionally, the winter monsoon index is negatively correlated with precipitation in Northwest China. The weakening of the winter monsoon 38 , which reduces the transport of dry and cold air to the study area, may increase precipitation to a certain extent.
The SPI is calculated from monthly precipitation. Therefore, the results of the trend analysis of SPI series were similar to the results of the precipitation series analysis at different scales. The frequencies of wet and dry years vary from station to station, but the frequency of normal years is the highest at each station. Seasonal SPI exhibited a decreasing trend in spring, summer and autumn at most stations and an increasing trend in winter. The annual SPI at most stations in the study area displayed a stable and sustained downward trend, which suggests that the water-wind erosion crisscross region of the Loess Plateau is currently associated with a drought trend, and the drought degree will likely continue to increase.
Analyses of precipitation data and SPI series can improve the efficiency of water resource use for hydroelectric and agricultural production. Although the study area was always in the trend of drought but the drier-and-hotter trend has been slightly restrained to varying levels since the implementation of the "Grain for Green" policy 39,40 . The trend of drought may cause more wind erosion and the increase in erosive precipitation depth enhanced the hydraulic erosion force, but the erosion status of the study area has been effectively inhibited by increasing the vegetation coverage since the "Grain for Green" policy 41,42 , it can be seen that vegetation restoration is an effective ecological means for improving the environmental conditions of study area to a certain extent, at the same time, vegetation cover can also reduce soil erosion 43,44 . Better precipitation and temperature conditions are also more suitable for the progress of vegetation restoration. The research areas studied as part of the "China National Scientific and Technical Innovation Research Project for 13th Five Year Plan" include the water-wind erosion crisscross region of the Loess Plateau. Vegetation restoration technology and the relationship between the Normalized Difference Vegetation Index 45 and climate change should be studied in this area to further improve the regional environment.

Material and methods
Study area and data. The water-wind erosion crisscross region on the Loess Plateau (Fig. 11) is located in the transition region between the hill and gully loess region and the Maowusu Desert (35°25′N-40°38′N,  103°00′E-13°53′E). Climate change is intense in this area, and vegetation is sparse. The northern part of the water-wind erosion crisscross region on the Loess Plateau is the wind erosion area, and the southern portion is the water erosion area.  The meteorological data used in this study were obtained from the China Meteorological Data Network (http://data.cma.cn/). To ensure synchronization between meteorological data and the long-term series of observations collected at each synoptic station, stations with missing data were excluded. A total of 14 representative stations were selected for this study (Table 6). Additionally, precipitation data from 1958 to 2015 were used.

Study Methods
Standardized precipitation index. The SPI is a drought index recommended by the WMO, and it has been used in many studies 4,46 . The calculation steps are as follows: , the probability distribution of precipitation related to the gamma function, is calculated as follows:  In Eq. (5), w s is the probability weighted moment, where s = 0, 1, 2, and l is the ordinal number of the precipitation sequence x in ascending order.  Table 5. Results of the statistical tests of seasonal and annual SPI in the study period. Z: M-K trend test; S: Spearman's rho test; and b (mm/a): slope of linear regression. Bold type represents trends identified by 2 statistical methods and trends that are statistically significant at the 5% level.
In addition, the slope β is calculated as follows.
The lag-1 serial correlation coefficient of sample data x i (designated as R h ) can be expressed as follows. In the two-sided test, R h was computed using the following equation at the 95% significance level 50 : where n is the sample size. For the lag-1 serial correlation coefficient (r 1 ) of the detrended series ′ X t , if r 1 falls within the range calculated by Eq. (10), the sample data are considered to be serially independent. In this case, linear trend estimation,  Spearman's rho test and Mann-Kendall test methods are directly applied to ′ X t . Otherwise, the sample data ′ Y t can be obtained by the following equation.
This prewhitening procedure after detrending the series is referred to as the trend-free prewhitening procedure. After applying the TFPW procedure, the residual series should be an independent series. c) The identified trend (T t ) and the residual ′ Y t are combined as follows.
The blended series (Y t ) includes a trend and noise and is no longer influenced by the serial correlation. Then, the M-K test is applied to the blended series to assess the significance of the trend.
Linear trend estimation. x i denotes a climate variable with a sample size of n, and t i denotes the time corresponding to x i . A linear regression equation can be established between x i and t i as follows.  where n is the length of the time series; R(X i ) is the rank of the time series at the observed value X i ; and positive and negative values of Z D indicate upward and downward trends, respectively. |Z D | > 2.08 rejects the hypothesis of no trend at the 5% significance level.

Mann-Kendall trend test.
In the Mann-Kendall (M-K) trend test 51,52 , the null hypothesis H 0 is that the time series data (x 1 , × 2 …, x n ) comprise a sample of n independent and identically distributed random variables.   The alternative hypothesis H 1 is a bilateral test: for all k, j ≤ n and k ≠ j when x k and x j are differently distributed. The statistical variable S is calculated as follows: where Sgn(x) is the sign function. The associated value is determined using the following equation.
S is normally distributed and has a mean of zero, and the variance is denoted as = − + Var n n n ( 1)(2 5)/18. When n > 10, the standard normal statistical variable Z can be calculated as follows.
In the bilateral trend test, if |Z| is ≥ α − Z 1 ( /2) at the given confidence level of α, the null hypothesis H 0 is unacceptable, i.e., the time series exhibits a significant upward or downward trend at the confidence level of α. Specifically, if Z > 0, the series displays an upward or increasing trend, and if Z < 0, the series shows a downward or decreasing trend. If |Z| ≥ 1.28, 1.64, or 2.32, the significant trend test is passed at 90%, 95%, or 99% confidence levels, respectively.  Notably, when x i > x j , r i = +1, and when x i < x j , r i = 0, where j = 1, 2, …, i. The statistical variable is defined based on the assumption of an independent random time series:

Mann-Kendall abrupt change test.
K k k k where UF 1 = 0, and s k and Var(s k ) are the mean and variance of the cumulative value s k , respectively. When × 1 , × 2 …, x n are independent in the same continuous distribution, s k and Var(s k ) can be calculated as follows. UF i , a standard normal distribution, is a statistical series calculated based on the order of time series x. At a given significance level α, |UF i | > U α indicates a significant trend in the series. The above calculation procedure is repeated in reverse order for time series x, and UB k = −UF k (k = n, n−1, …, 1), where UB 1 = 0. Then, UF k > 0 indicates an upward or increasing trend, UF k < 0 indicates a downward or decreasing trend, and exceeding the critical line indicates a significant trend. If the two curves formed by UF k and UB k intersect and if the intersection appears between the boundary lines, the corresponding time of the intersection is the start time of an abrupt change.

SPI Classification
≥2.  Rescaled range analysis. In this study, the rescaled range analysis (R/S) method was used to calculate the Hurst coefficient and predict and quantitatively describe the long-term correlation of the time series 53,54  where Z u is the sequential cumulative deviation in subsequence u i . R u /S u is the rescaled range of each subsequence, and S u is the standard deviation of subsequence u i . The logarithmic processing of the empirical Hurst formula (R N /S N = ωN H ) can be expressed as follows: where H indicates the Hurst index and ω is a constant. H has three value areas with different physical significance: (a) H ∈ (0, 0.5) indicates that the direction of the overall trend in the future is opposite that of the previous trend; (b) H ∈ (0.5, 1) indicates a time series with long-term positive autocorrelation, i.e., a change in the overall direction of the time series will reflect the previous trend; and (c) H = 0.5 indicates a completely uncorrelated series. For time series, the result of each time range is absolutely independent. Therefore, changes in time series are completely random.