Spatial and Temporal Characteristics of Precipitation and Potential Influencing Factors in the Loess Plateau before and after the Implementation of the Grain for Green Project

Abstract

Since the implementation of the Grain for Green Project (GFGP) in the 1990s, the warming and wetting trend in the Loess Plateau is becoming statistically significant in the context of climate change. However, the correlation between precipitation increase and the regional vegetation restoration is still controversial. To explore the main factors influencing the regional precipitation change, this study selected five potential influencing factors including potential evapotranspiration (PET), normalized difference vegetation index (NDVI), precipitable water (PW), surface temperature (ST), and water vapor transport (WVT). We used the statistical methods to analyze the spatial-temporal distribution of precipitation before and after the GFGP and to quantify the relative influence degree of different factors to precipitation change. The results show that: (1) The precipitation increased significantly (95% confidence level) after the GFGP, with an increase rate of 4.96 mm a⁻¹; (2) from the perspective of spatial-temporal distribution, the precipitation in the southern part of the Loess plateau was significantly increasing with an increase rate of 20–50 mm in the period of 2000–2014; (3) the relative influence degree of NDVI to precipitation increased after the GFGP, and the annual precipitation (PRE_A) and summer precipitation (PRE_S) was more influenced by NDVI (relative influence degree of 30.18% and 31.37%, respectively) compared with winter precipitation. In winter, the PW and the PET are the main influencing factors for the precipitation change with relative influence degrees of 30.13% and 27.64%, respectively. Based on this study, we speculate that the warming and wetting trend of the Loess Plateau in recent years is not only closely related to global climate change, but also significantly affected by local climate change brought by vegetation restoration. The above conclusions are important for future ecological restoration and water resources management in the water-scarce Loess Plateau.

1. Introduction

The Loess Plateau is located in a semi-arid and semi-humid region. The main land use types are farmland and grassland [1]. At the same time, this area is the largest loess-covered area in the world. With the impact of climate change (extreme precipitation and drought) and human activities (overgrazing and deforestation) [2,3], a series of ecological and environmental problems, such as soil erosion, drought, and water scarcity, and siltation of rivers and lakes on the Loess Plateau have been aggravated [3]. Therefore, the Chinese government implemented the “Grain for Green Project” (GFGP) in the 1990s in the Loess Plateau. The implementation of the GFGP has led to dramatic changes in land use/land cover, with the vegetation cover increasing by 25% from 1999 to 2013 [4], and the change rate of normalized difference vegetation index (NDVI) in the new century is 10 times higher than that before 1999 [1]. With the evapotranspiration increase and soil drying induced by vegetation recovery, observed stream flow has also declined by approximately 50% from the 1950s to 2015 [5]. The improved vegetation cover, in combination with engineering measures (e.g., check dam, terrace, level furrow, and fish-scale pits), has improved soil properties, reduced soil erosion, and saved incoming water to prevent sediment load from entering the river [6]. As a result, sediment discharge into the Yellow River has declined by approximately 90% over the past 60 years [5]. In general, the vegetation restoration effectively reduced the surface runoff and river sediment load [7], but it also brought out some problems of water resources.

Direct local-scale observations suggest a considerable reduction in surface water yield over the region after revegetation, including soil desiccation, stunted or “small old” trees, and even plant mortalities. Extra planted vegetation consumes water, and when evapotranspiration exceeds precipitation, a soil water imbalance occurs, which plays a crucial role in the formation of dried soil layers. The dried soil layer is, in turn, deemed a serious obstacle to sustainable vegetation restoration, suggesting a negative feedback of revegetation on local water yield and subsequently on vegetation growth. For a long time, scholars believed that the Loess Plateau tends to be warm and dry, but in recent years, especially since the 1990s, there has been a significant warming and wetting trend in this region [2,8,9,10,11]. The debate on the ecological sustainability of the uneven distribution of water resources on the Loess Plateau continues, but the sustainability of afforestation in the Loess Plateau region needs to be judged according to the local hydrology and meteorology. The knowledge of the relationship between local precipitation and afforestation characteristics would indeed be highly valuable for a more comprehensive evaluation of the impacts of revegetation on surface water yield, thus guiding the ecological restoration implementation.

Precipitation plays a critical role in the human life and hydrological cycle [12,13]. The change in precipitation leads to an aggravation of natural hazards such as floods, drought, and soil erosion [14]. At present, there are many studies on the influence of atmospheric circulation [15], meteorological factors [16], and aerosol factors [17,18] on precipitation. However, there are few studies paying attention to the impact of vegetation restoration and the local microclimate on precipitation. Gao et al. [2] analyzed the changes in actual evapotranspiration (AET), normalized difference vegetation index (NDVI), and precipitation in the Loess Plateau during 1990–2014, and they found that the precipitation increase may cause an increase in AET, and in turn, the AET increase may also influence the precipitation. Li et al. [9] analyzed the impacts of the Chinese Grain for Green program and climate change on vegetation in the Loess Plateau during 1982–2015, and found that precipitation has mainly contributed to the increase in NDVI in the southeastern region; meanwhile, they found that with the increase in vegetation coverage, the change in underlying surface also had a feedback effect on regional precipitation. Therefore, to further analyze the precipitation change and to identify its potential driving factors in the Loess Plateau, it is urgent to comprehensively consider both climate factors and vegetation restoration factors, and quantitatively analyze the main contribution factors of precipitation under the condition of global climate change and fast land use/cover change. At present, the analysis methods of precipitation driving factors mainly include regression analysis [19], machine learning [15], and principal component analysis [20]. Regression analysis cannot solve the problem of multicollinearity among predictors. Principal component analysis is mainly used to extract important information to explain data structure [21], while principal component regression analysis can be used for driver factor analysis and is also widely used in predictive analysis. Chen et al. [21] analyzed the spatial variation in driving factors of irrigation water consumption based on principal component regression, and pointed out that the water consumption structure, irrigation technique, and planting structure were major influential factors in most provinces of China. Tang et al. [22] integrated principal component analysis with statistically based models for analysis of causal factors and landslide susceptibility mapping; the results showed that rainfall and land use were essential in predicting both loess landslide and rockfall occurrences. Zeng [23] analyzed the influence of NDVI, EVI, precipitation, and temperature on surface water fluctuation based on regression analysis, and pointed out that vegetation cover is an important factor in controlling permanent water changes. Therefore, we can use principal component regression analysis to eliminate multicollinearity in independent variables, and the weights coefficient can be used to measure the influence of the factors to precipitation.

In the past decades, jointly influenced by climate change and ecological restoration, the spatial and temporal pattern of precipitation in the Loess Plateau was undergoing a significant change. However, the potential influencing factors and their influence on local precipitation are still unclear and not well revealed. Therefore, this study aims to analyze the main influencing factors of precipitation in the condition of vegetation recovery and climate change. Based on previous studies and the local conditions of the Loess Plateau, we selected five potential influencing factors, including potential evapotranspiration (PET), normalized difference vegetation index (NDVI), precipitable water (PW), surface temperature (ST), and water vapor transport (WVT). The paper is structured as follows: (1) To collect the precipitation data and analyze the spatial and temporal characteristics of precipitation before and after the GFGP in the Loess Plateau; (2) to analyze the factors and reveal the spatial and temporal characteristics of each factor; (3) to use the statistical methods to quantitatively analyze the main influencing factors of precipitation before and after the GFGP, and discuss the impact of vegetation restoration on precipitation. The above findings in this research may provide the basic reference for the precipitation prediction in the future. Meanwhile, it may be a great significance for the ecological restoration and water resource management of the Loess Plateau.

2. Study Area and Data Collection

2.1. Study Area

The Loess Plateau is situated in the upper and middle reaches of the Yellow River (as shown in Figure 1), with a latitude between 33°43′ and 41°16′ N and longitude between 100°54′ and 114°33′ E. The plateau is covered with a deep layer of loess with an average thickness of 50–200 m [24]. The overall terrain is high in the northwest and low in the southeast, ranging from 94 to 5000 m and with an area of 648,700 km², approximately.

The Loess Plateau belongs to two temperate zones (warm temperate zone in the south and middle temperate zone in the north). It has the typical continental monsoon climate characteristics. The annual average temperature is from 4.3 °C to 14.3 °C, and the annual precipitation decreases gradually from the southeast to northwest, ranging from 200 to 750 mm. The total precipitation in the rainy season (June–September) accounts for about 65% of the annual precipitation [25].

2.2. Data

This study mainly used the following data:

(1) Meteorological data. Meteorological data were collected from the China Meteorological Data Network (http://data.cma.cn/), and the temporal resolution of the data is 1 day. Meteorological data mainly include temperature, air pressure, wind speed, relative humidity, sunshine duration, and ground temperature data of 52 meteorological stations in the Loess Plateau, as shown in

Table 1. The names and coordinates of 52 meteorological stations.

Station	Name	Lon	Lat	Station	Name	Lon	Lat
52,765	Menyuan	101.62	37.38	53,740	Hengshan	109.23	37.93
52,787	Wushaoling	102.87	37.20	53,754	Suide	110.22	37.50
52,866	Xining	101.77	36.62	53,764	Lishi	111.10	37.50
52,876	Minhe	102.85	36.32	53,772	Taiyuan	112.55	37.78
52,895	Jingyuan	104.68	36.57	53,806	Haiyuan	105.65	36.57
52,983	Yuzhong	104.15	35.87	53,810	Tongxin	105.90	36.98
52,984	Linxia	103.18	35.58	53,817	Guyuan	106.27	36.00
52,986	Lintao	103.87	35.37	53,821	Huanxian	107.30	36.58
52,996	Huajialing	105.00	35.38	53,853	Xixian	110.95	36.70
53,336	Wulatehouqi	108.52	41.57	53,863	Jiexiu	111.92	37.03
53,446	Baotou	109.85	40.67	53,868	Linfen	111.50	36.07
53,463	Huhehaote	111.68	40.82	53,903	Xiji	105.72	35.97
53,478	Youyu	112.45	40.00	53,915	Pingliang	106.67	35.55
53,513	Linhe	107.42	40.75	53,923	Xifengzhen	107.63	35.73
53,519	Huinong	106.77	39.22	53,929	Changwu	107.80	35.20
53,529	Etuokeqi	107.98	39.10	53,942	Luochuan	109.50	35.82
53,543	Hedong	109.98	39.83	53,959	Yuncheng	111.02	35.03
53,564	Hequ	111.15	39.38	53,975	Yangcheng	112.40	35.48
53,614	Yinchuan	106.22	38.48	56,080	Hezuo	102.90	35.00
53,615	Taole	106.70	38.80	56,093	Minxian	104.02	34.43
53,646	Yulin	109.70	38.23	57,034	Wugong	108.22	34.25
53,663	Wuzhai	111.82	38.92	57,046	Huashan	110.08	34.48
53,664	Xingxian	111.13	38.47	57,051	Sanmenxia	111.20	34.80
53,705	Zhongning	105.67	37.48	57,067	Lushi	111.03	34.05
53,723	Yanchi	107.40	37.78	57,071	Mengjin	112.43	34.83
53,738	Wuqi	108.18	36.83	57,077	Luanchuan	111.60	33.78

and Figure 1. The precipitation data were collected from 368 rainfall gauge stations in and around the Loess Plateau region, and the spatial distribution is shown in Figure 1. The temporal resolution of the data is 1 day.

(2) Reanalysis Dataset. NCEP/NCAR reanalysis dataset was created through the cooperative efforts of the National Centers for Environmental Prediction (NCEP) and National Center for Atmospheric Research (NCAR), available at (https://psl.noaa.gov/data/gridded/data.ncep.reanalysis.html). The spatial resolution of the dataset is 2.5° × 2.5°, and the temporal resolution is 1 day. The dataset used in this study included u/v wind, and specific humidity for 1000, 925, 850, 700, 600, 500, 400, and 300 hPa; precipitation water (PW).

(3) Remote sensing data. Before 1999, the NDVI data were from the National Earth System Science Data Center, National Science & Technology Infrastructure of China (http://www.geodata.cn). After 1999, the NDVI data were from MODIS datasets (https://ladsweb.modaps.eosdis.nasa.gov/search/).

3. Methodologies

3.1. Technical Framework

This paper aims to study the change in precipitation before and after the implementation of the Grain for Green Project (GFGP) in the Loess Plateau and to analyze the main influencing forces behind this. First, the Mann–Kendall (M-K) test and Hurst analysis method are used to study the trend and persistence of precipitation before and after the GFGP in the Loess Plateau at the spatial-temporal scale. Then, five main factors were selected, including PET, NDVI, PW, ST, and WVT, and their changes were analyzed at spatial-temporal scales. Finally, the statistical method was implemented to quantitatively evaluate the main influencing factors of precipitation before and after the GFGP. The research framework of this article is as shown in Figure 2.

3.2. Trend Analysis

The Mann–Kendall (M-K) test and Hurst analysis method were introduced to analysis the changing trends of precipitation and the other five factors on the Loess Plateau. Based on the M-K test, the trends of the five factors and precipitation before and after the GFGP are analyzed, and the Hurst analysis method is used to analyze the continuity of precipitation after the GFGP.

3.2.1. Mann–Kendall Trend Test with Trend-Free Pre-Whitening

Mann–Kendall is a nonparametric analysis used to estimate the significance of the trend of hydrological and meteorological [26]. The evaluating variable Z is often used to estimate whether a series has a significant trend or not, which was calculated as follows:

$$Z=\left\{\begin{array}{c}\frac{S-1}{\sqrt{Var\left(S\right)}} when S>0\\ 0 when S=0\\ \frac{S+1}{\sqrt{Var\left(S\right)}} when S<0\end{array}\right.$$

(1)

$$S={\displaystyle \sum }_{i=1}^{n-1}{\displaystyle \sum }_{j=i+1}^{n}sign\left(X_j-X_i\right)$$

(2)

$$sign\left(X_j-X_i\right)=\left\{\begin{array}{c}+1 when \left(X_j-X_i\right)>0\\ 0 when \left(X_j-X_i\right)=0 \\ -1 when \left(X_j-X_i\right)<0\end{array}\right.$$

(3)

where $X_i$ and $X_j$ are the sequential values, and $n$ is the length of the time series.

The variance of S is given by the following equation:

$$Var\left(S\right)=\frac{n\left(n-1\right)\left(2n+5\right)}{18}$$

(4)

where $n$ is the length of the data series.

The evaluating variable $Z$ indicates the trend of the time series; when $Z>0$, the trend is increasing; otherwise, it is decreasing. When |$Z$| ≥ 1.64, the trend is significant at 90% confidence; when |$Z$| ≥ 1.96, the trend is significant at 95% confidence [27].

3.2.2. Hurst Exponent and Rescaled Range (R/S) Analysis

The Hurst analysis is a useful method to understand the properties of a time series [28,29] and is widely used in hydrology, climatology, and ecology. Rescaled range (R/S) analysis is the method widely used to calculate the H exponent. The basic principles of the R/S analysis method are as follows:

Suppose there is a time series ξ (t) = {ξ |t = 0, 1, 2, ⋅⋅, N}, and the time series is cut into m (m = N/n) continuous subsequences with a length of n (n < N); then, each time subsequence can be presented as:

$${\xi }_n={\displaystyle \sum }_{i=1}^n{\xi }_i$$

(5)

The average value of a subsequence can be presented as:

$$\overline{{\xi }_n}=\frac{1}{n}{\xi }_n$$

(6)

Then, the ith cumulated dispersion is:

$$X\left(i,n\right)={\displaystyle \sum }_{t=1}^i\left({\xi }_i-\overline{i{\xi }_n}\right)$$

(7)

The range is:

$$R\left(n\right)=\begin{array}{c}max\\ 0\le i\le n\end{array}X\left(i,n\right)-\begin{array}{c}min\\ 0\le i\le n\end{array}X\left(i,n\right)$$

(8)

where 1 ≤ i ≤ n.

The standard deviation is:

$$S\left(n\right)={\left[\frac{1}{n}{\displaystyle \sum }_{i=1}^n{\left({\xi }_i-\overline{{\xi }_n}\right)}^2\right]}^{1/2}$$

(9)

When the characteristics of the time series ξ (t) = {ξ⋅|t = 0, 1, 2, ⋅⋅, N} is a constant scale, the rescaled range of R/S for any length n can be presented as:

$$\frac{R\left(n\right)}{S\left(n\right)}={(\alpha n)}^H$$

(10)

where α is a constant, and H is the Hurst exponent and estimated by fitting the above formula. When H = 0.5, the time series is random, indicating that the future trend is independent of the past. When 0.5 < H < 1, the future change trend is consistent with the past. The closer the value is to 1, the stronger the continuity; when 0 < H < 0.5, the future change trend is opposite to the past. The closer the value is to 0, the stronger the anti-continuity.

3.3. Calculation of the Five Main Factors

In this study, we aim to analyze the main driving factors of precipitation in the condition of vegetation recovery and climate change, and comprehensively considered five factors. PET is selected as the climate factor, NDVI and ST as the underlying surface change factor, and PW and WVT as the water vapor factor of atmospheric circulation. Moreover, CO₂ is also a potential influencing factor worthy of consideration. However, the previous study showed that the elevated CO₂ was not the most primary factor for ET change in the Three-North Region including the Loess Plateau in China [30]. Meanwhile, there are no gridded data of CO₂ concentration from 1985 to 2014 in the world. Therefore, we decided not to consider the effect of CO₂ in this study. Finally, the five main factors are calculated as follows.

3.3.1. Potential Evapotranspiration (PET)

Potential evapotranspiration (PET) is a typical indicator that reflects the comprehensive effect of the climate [31]. Therefore, we selected PET as the climate factor. PET is estimated through the Penman–Monteith (P–M) method that is based on the meteorological data and parameters, which is widely used and recommended by numerous international organizations. The formula is as follows:

$$\mathrm{PET}=\frac{0.408\Delta \left(R_n-G\right)+\gamma ·\frac{900}{273+T}·u_{2·}\left(e_a-e_d\right)}{\Delta +\gamma \left(1+0.34·u_2\right)}$$

(11)

where PET is the potential evapotranspiration, in mm/day; $R_n$ is the net solar radiation, in MJ/m²; G is the soil heat flux, in MJ/m²;$\Delta$ is the slope of the saturation vapor pressure curve, in kPa/°C;$\gamma$ is the psychometric constant, in kPa/°C; T is the daily temperature, in °C; $e_a$ is the saturation vapor pressure, in kPa; $e_d$ is the actual water vapor pressure, in kPa; $u_2$ is the daily average wind speed at 2 m height, in m/s. The data of temperature, air pressure, wind speed, relative humidity, and sunshine duration were collected from the 52 meteorological stations on the Loess Plateau. $R_n$ refers to the surface ingoing energy minus the surface outgoing energy, which can be expressed as:

$$R_n=\left(1-\alpha \right)K_{in}+\left(L_{in}-L_{out}\right)-\left(1-e_0\right)L_{in}$$

(12)

where $K_{in}$ is the incoming short-wave radiation (W/m²); $L_{in}$ and $L_{out}$ are the incoming long-wave radiation and outcoming long-wave radiation (W/m²), respectively; $\alpha$ is the surface short-wave albedo; and $e_0$ is the land surface emissivity.

3.3.2. Normalized Difference Vegetation Index (NDVI)

Vegetation coverage in the Loess Plateau has increased significantly after the GFGP. The impact of vegetation changes on precipitation should be considered in this study. Therefore, we selected Normalized Difference Vegetation Index (NDVI) as one of the five factors. NDVI data before 1999 were collected from the National Earth System Science Data Center, National Science & Technology Infrastructure of China (http://www.geodata.cn), which has a spatial resolution of 0.05° and a temporal resolution of 8 days; NDVI data after 1999 were collected from the MODIS datasets (https://ladsweb.modaps.eosdis.nasa.gov/search/), with a spatial resolution of 250 m and a temporal resolution of 16 days.

In order to keep the spatial consistency of the NDVI, we used the maximum value composite method [32] to synthesize the 8 days and 16 days NDVIs into monthly NDVIs, and resampled them to 0.05° × 0.05° resolution based on WGS-1984. Annual NDVI was defined as the average monthly composite NDVI of 12 months, summer NDVI was defined as the average monthly composite NDVI from April to September, and winter NDVI was defined as the average monthly composite NDVI from October to March.

3.3.3. Precipitable Water (PW)

Precipitable water (PW) plays a decisive role in latent heat transportation, cloud formation, and the occurrence of precipitation [33,34,35], which is defined as the depth of liquid water produced after all the water vapor in the unit air column has condensed and landed. In this study, we selected PW as one of the main factors. Zhao et al. [36] found that NCEP reanalysis data were highly correlated with the observational data from the perspective of climatological statistics. Therefore, the data of PW were collected from the NCEP/NCAR reanalysis dataset with a spatial resolution of 2.5° × 2.5° and a temporal resolution of one day.

3.3.4. Surface Temperature (ST)

Temperature is the main factor that affects precipitation [37,38]. After the GFGP, the vegetation coverage of the Loess Plateau changed greatly, and vegetation variation can affect the surface temperature (ST) by altering surface albedo and evapotranspiration [39]. Therefore, we chose ST as the temperature factor of the underlying surface change. The surface temperature comes from the China Meteorological Data Network (http://data.cma.cn/), including surface temperature data at 52 sites in the Loess Plateau, with a temporal resolution of one day.

3.3.5. Water Vapor Transport

Precipitation originates from local evaporation moisture and external water vapor transport (WVT) [40]. Water vapor transport (WVT) represents the amount of water vapor transfer through a unit area in a unit interval. In this study, we selected WVT as the external water vapor factor. WVT is calculated as follows [41,42]:

$$Q_x=-\frac{1}{g}{{\displaystyle \int }}_{p_s}^{p}uqdp$$

(13)

$$Q_y=-\frac{1}{g}{{\displaystyle \int }}_{p_s}^{p}vqdp$$

(14)

where $Q_x$ and $Q_y$ indicate the zonal and meridional water vapor transportation, respectively, $g$ is the gravitational acceleration, $p_s$ is the surface pressure, $p$ is the pressure at 300 hPa, $q$ is the specific humidity, and u/v is the zonal and meridional wind.

According to the spatial resolution of the NCEP/NCAR reanalysis dataset, this study divided the study area into 22 grids of 2.5° × 2.5°, as shown in Figure 1. The east boundary of grids 21, 17, 11, and 5 is defined as the east boundary. The south boundary of grids 7, 8, 3, 4, and 5 is the south boundary of the region. The west boundary of grids 19, 15, 7, and 3 is the west boundary. The north boundary of grids 7, 8, 19, 20, and 21 is the north boundary. The water vapor transport of the east, south, west, and north boundary is the sum of water vapor transport volume at the corresponding grid boundary. The zonal water vapor net input is equal to the results of water vapor flux of the western boundary subtracted by those of the eastern boundary, and the meridional component is equal to the southern boundary minus the northern boundary [41,43]. Then, the net water vapor input (WVT) of the Loess Plateau is calculated.

3.4. Principal Component Regression Analysis

In order to eliminate the multi-collinearity among various factors, principal component regression analysis was used to quantitatively analyze the impact of each factor on precipitation. Principal component analysis can transform multiple indicators into fewer well-represented comprehensive indicators that reduce dimensionality and simplify the data structure [44,45,46].

Suppose that there are p indicators involved with a certain object, denoted by $X_1$,$X_2$,...,$X_p$. These p indicators compose a p-dimensional vector: X = ($X_1$,$X_2$,...,$X_p$). After dimension reduction, new variables $F_1$, $F_2$, …, $F_m$ (m ≤ p) are obtained:

$$\left\{\begin{array}{c}{F}_1=l_{11}{x}_1+l_{12}{x}_2+\dots +l_{1p}{x}_p\\ F_2=l_{21}{x}_1+l_{22}{x}_2+\dots +l_{2p}{x}_p\\ \dots \mathrm{..}\\ F_m=l_{m1}{x}_1+l_{m2}{x}_2+\dots +l_{mp}{x}_p\end{array}\right.$$

(15)

In this study, the principal component analysis method is used to convert five factors into several independent principal components. Then, combined with multiple linear regression, the principal component with a contribution rate more than 90% is selected for regression analysis [21,22], and the relative influence degree of each factor is calculated by the weight of the regression coefficient, which is calculated as follows:

$$\mathrm{Y}=a_1{X}_1+a_2{X}_2+\dots +a_{\mathrm{P}}{X}_{\mathrm{P}}$$

(16)

$${\eta }_i=\frac{\left|a_i\right|}{\left|a_1\right|+\left|a_2\right|+\left|a_3\right|+\dots \left|a_p\right|}$$

(17)

where Y is the dependent variable, $X_{\mathrm{i}}\mathrm{and}{a}_{\mathrm{i}}$ are the independent variables and the corresponding regression coefficients, and ${\eta }_i$ is the relative influence degree of the change in $X_{\mathrm{i}}$ to the change in Y. We eliminated the influence of dimensionality through data standardization [47], and the statistical calculations in this study were through SPSS 20.

4. Results

4.1. The Spatial-Temporal Characteristic of Precipitation Changes in the Loess Plateau

4.1.1. The Temporal Variability of Precipitation

In this study, we defined the period before the GFGP as base Period I (1985–1999), and the period after the GFGP as change Period II (2000–2014). The precipitation in the Loess Plateau is concentrated in summer and autumn, and the seasonal precipitation fluctuates greatly. Therefore, we defined the summer scale as April to September and the winter scale as October to March. Then, the temporal variability of precipitation is shown in Figure 3.

As shown in Figure 3, during Period I, the annual precipitation (PRE_A), summer precipitation (PRE_S), and winter precipitation (PRE_W) all showed a decreasing trend. The PRE_A decreased at a rate of 2.22 mm a⁻¹, which is most significant. The PRE_S and PRE_W decreased at a rate of 1.56 and 1.19 mm a⁻¹, respectively. However, based on the M-K trend test, the PRE_A, PRE_S, and PRE_W decreasing trend is not statistically significant. In Period II, the PRE_A and PRE_S showed a significant increasing trend. The PRE_A increased at a rate of 4.96 mm a⁻¹ and the PRE_S increased at a rate of 5.63 mm a⁻¹. Moreover, based on the M-K test, the Z values of both PRE_A and PRE_S were 2.43 with a 95% confidence level. The PRE_W showed a slight decreasing trend at a rate of 0.47 mm a⁻¹. In general, precipitation in the Loess Plateau was mainly affected by the PRE_S, and which had increased significantly after the GFGP.

4.1.2. Changes and the Continuity of Station-Scale Precipitation in the Loess Plateau before and after GFGP

Precipitation in the Loess Plateau has increased significantly since 1999. In order to further analyze the trend and its significance of precipitation at each site scale before and after the GFGP, this study adopted the M-K test to calculate the precipitation trends at 368 meteorological stations in Period I and Period II. Furthermore, the continuity of the trend of precipitation after the GFGP was tested based on the Hurst analysis method.

As shown in Figure 4a,b, based on the Z value from the M-K test, there were 267 stations showing decreasing trends (Z < 0) in Period I, accounting for 72.6% of all stations. Except for a few stations in the northeast area, almost all the stations showed decreasing trends, in which 17 stations were with a 90% confidence level. In Period II, there were 311 stations showing an increasing trend (Z > 0), accounting for 84.5% of all stations, in which 37 stations were with a 90% confidence level. The central and northern parts of stations showed a significant increasing trend. At the same time, based on Hurst analysis, the H value of precipitation in the central, northern, and southwestern Loess Plateau after 1999 was all greater than 0.55 (Figure 4c), which means that the precipitation had a strong continuity under the current conditions. Generally, precipitation in the Loess Plateau was dramatically different before and after the GFGP. After 1999, precipitation in the central and the northeastern region had increased significantly.

4.1.3. The Spatial Variability of Precipitation

The spatial distribution and variation in precipitation on the Loess Plateau in Period I and Period II, and the difference in precipitation between the two periods are shown in Figure 5.

The precipitation on the Loess Plateau decreased from the southeast to northwest, and the spatial distribution and the difference distribution of the PRE_S were highly consistent with the PRE_A. In Period II, the PRE_A was 144.73–713.18 mm, decreasing from the southeast to northwest, and the PRE_S decreased from the southeast to northwest with 132.09–604.81 mm. From the difference in the precipitation in Figure 5g–i, the PRE_A and the PREs in the southeast increased significantly, with the rates mainly between 20 and 50 mm, and mostly all the regions showing increasing trends. The variation in the PRE_W was mainly between −8 and 2.5 mm. To sum up, after implementing the GFGP, the PRE_A and PRE_S in the southeast and central regions increased, while the spatial variation in the PRE_W changed little.

4.2. The Spatial-Temporal Characteristic of the Main Influencing Factors

4.2.1. The Annual and Seasonal Variation in the Main Factors

The annual and seasonal variation in the five factors on the Loess Plateau from 1985 to 2014 are shown in Figure 6, and the significance of the changing trend of each factor is tested based on M-K analysis.

As shown in Figure 6, the PET_A and PET_S changed from the increase trend in Period I to the decrease trend in Period II, and the PET_W increased all the time. The NDVI_A and NDVI_W in Period I showed a slight upward trend; moreover, in Period II, NDVI_A, NDVI_S, NDVI_W all showed an increasing trend with a 95% confidence level. It should be noted that the NDVI showed a gradually increasing trend from 2000 to 2003, but increased afforestation would increase transpiration and lead to water shortages; vegetation growth is restricted, which may be one of the reasons that NDVI showed a slightly decreasing between 2004 and 2006. After 2009, the NDVI gradually increased, which was mainly due to the surviving vegetation beginning to flourish after adapting to the local environment. Therefore, there may be a lag period for vegetation restoration. The PW_A and PW_W both showed a decrease trend before and after the GFGP, but PW_S kept increasing. The ST showed an increasing trend, but the increase rate of ST in Period II was lower than that in Period I. The WVT_A, WVT_S, and WVT_W showed a statistically significant downward trend with a 95% confidence level in Period I, but increasing in Period II. Together, from the temporal trend of each factor, the PET, NDVI, and WVT changed before and after the GFGP.

4.2.2. The Spatial Variability of the Main Factors

To further analyze the spatial distribution and variation of each factor in Period I and Period II, the spatial distribution of each factor in the annual and seasonal scales and spatial difference between the two periods were mapped, respectively, in Figure 7a–c. It is worth noting that the WVT is the water vapor flux value of the entire Loess Plateau, which cannot present spatial distribution.

As shown in Figure 7a–c, from the perspective of spatial distribution, the spatial distribution, NDVI, and PW of the annual, summer, and winter scales showed a decreasing trend from the southeast to northwest, similar to the spatial distribution of precipitation. The higher values of PET are distributed in the northwest, central, and southeast. The lower values of ST are distributed in the southwest. From the perspective of spatial variation, the NDVI in the southeast increased significantly in the annual, summer, and winter scales, similar to those of precipitation. At the same time, PET increased in the southwest and northeast, ST increased in the whole area, while PW decreased in the whole area. To sum up, NDVI was similar to precipitation in spatial distribution and variation. Meanwhile, PET, NDVI, and ST all increased after the GFGP.

4.3. Identification of the Main Influencing Factors of Precipitation Change

4.3.1. Correlation Analysis

This study drew scatter plots to analyze the correlation of precipitation and five factors before and after the GFGP, as shown in Figure 8. From the perspective of the annual and summer scales, the precipitation was positively correlated with PW and WVT, and negatively correlated with PET and ST. PW and WVT are water vapor factors, and the increase in water vapor had a promoting effect on the occurrence of precipitation. Moreover, NDVI was positively correlated with precipitation after the GFGP, and the correlation coefficient r was 0.82 and 0.85. This may be due to the increase in vegetation coverage having enhanced the regional water cycle and then increased the precipitation. In terms of winter, the precipitation was positively correlated with PW and WVT, and negatively correlated with other factors. The moisture sources of PRE_W mainly come from atmospheric circulation. Therefore, with the decrease in PW and WVT, the winter precipitation will also decrease. In summary, we found that the five factors had an obvious linear correlation with precipitation through correlation analysis. Therefore, in the next part, the principal component regression analysis will be used to further quantify the relative influence degree of each factor to precipitation before and after the GFGP.

4.3.2. The Relative Influence Degree of the Five Factors to Precipitation Change

Based on the correlation of precipitation and main factors, we further calculated the influence degree of each factor to precipitation through principal component regression analysis. The extraction of principal components and evaluation parameters of the principal component regression analysis are shown in

Table 2. The extraction of principal components and evaluation parameters.

Scale		Cumulative Contribution Rate/%	Number of Principal Components	Adjusted-R2	Sig.
Annual	Period I	98	4	0.50	0.026
	Period II	98	4	0.45	0.038
Summer	Period I	98	4	0.67	0.004
	Period II	96	4	0.67	0.003
Winter	Period I	96	4	0.58	0.011
	Period II	96	4	0.39	0.048

. We converted five original factors into four linearly uncorrelated variables, and the cumulative contribution rate was greater than 96%. That is, the selected principal component could retain almost all the information of the original factors. At the same time, the significance test parameter (sig.) showed that all the regression models had passed the significance test (α = 0.05), and the adjusted R² was more than 0.39. In general, the principal component regression analysis can be used to measure the impact of factors on precipitation.

The influence degree of factors to precipitation before and after the GFGP is shown in Figure 9. In terms of annual scale, PW had the greatest influence on precipitation in Period I, with a relative influence degree of 41.33%; in Period II, NDVI and PW had the greatest impact on precipitation, with relative influence degrees of 30.18 and 26.63%, respectively. In terms of the summer season, PET and PW influenced precipitation most in Period I, with relative influence degrees of 32.17% and 22.63%, respectively; while in Period II, NDVI and PW had the greatest impact on precipitation, with relative influence degrees of 31.37 and 24.26%, respectively. In terms of the winter season, WVT and PW had the greatest influence on precipitation in Period I, with relative influence degrees of 41.21% and 25.50%, respectively; while in Period II, PW and PET had the greatest impact on precipitation, with relative influence degrees of 30.13% and 27.64%, respectively.

Based on the results of statistical analysis, we compared the main factors of precipitation before and after the GFGP and found that, in terms of annual scale, the relative influence degree of NDVI to precipitation changed significantly from 14.80% before the GFGP to 30.18% after the GFGP. In terms of summer scale, the relative influence degree of NDVI to precipitation increased significantly from 9.78% before the GFGP to 31.37% after the GFGP. Therefore, we speculate that vegetation restoration had a certain promoting effect on precipitation. The substantial increase in vegetation coverage in summer had increased the water consumption through transpiration and evaporation; then, the regional water consumption increased significantly along with the accelerated local water cycle, which might lead to the significant increase in the precipitation. PW and PET were the main factors affecting PRE_W after the GFGP. According to the scatter plot of Figure 8, PRE_W was positively correlated with PW and negatively correlated with PET; PW decreased and PET increased during Period II, which jointly leads to the decrease in PRE_W. In summary, with the increase in precipitation in recent years, the Loess Plateau presented a trend of warming and wetting, and we speculated that this might be due to global climate change on the one hand and the changes in the local climate brought by vegetation restoration on the other hand.

5. Discussion

5.1. Identification of Main Factors in Sub-Regions

After the implementation of the GFGP in the Loess Plateau in 1999, the vegetation coverage in the central and western part of the loess Plateau increased significantly [9], while the NDVI in the northwest part obviously did not increase [48]. The vegetation restoration speed was uneven in spatial distribution. Therefore, to further analyze the contribution rate of each factor to precipitation in different regions of the Loess Plateau, the study area was divided into 22 sub-regions according to the 2.5° × 2.5 grid, and the most influence factors of each sub-region were identified based on the principal component regression method. The results are shown in Figure 10.

The influence degree of factors to PRE_A are shown in Figure 10a,b. PW was the main influencing factor in 63.6% of the regions, with an average relative influence degree of 36.2% in Period I. During Period II, NDVI was the main influencing factor in 59.1% of the regions, with an average relative influence degree of 35.5%. In the central area, the relative influence degree of NDVI to precipitation was higher than that before the GFGP. The main influencing factors for the PRE_S are shown in Figure 10c,d; during Period I, PET was the main influencing factor in 68.2% of the regions, with an average relative influence degree of 35.1%. In Period II, the main influencing factor in 72.7% of the regions was NDVI, with an average relative influence degree of 36.7%, and which had increased significantly in the central and eastern regions. The influence degree of factors of the PRE_W are shown in Figure 10e,f. During Period I, PW was the main influence factor in 68.2% of the regions, with an average influence degree of 34.3%. During Period II, the main influencing factor was PW, with an average relative influence degree of 38.1%.

From the results, this study found that NDVI was the main influencing factor to the annual and summer precipitation after the GFGP. The mechanism may be that, with the increase in vegetation coverage, the water consumption of vegetation and evapotranspiration of soil water increase. The regional water cycle accelerates, and the local microclimate has improved, which has a positive feedback effect on precipitation [49,50]. As precipitation increases, the water demand for vegetation growth is satisfied. Therefore, the vegetation coverage is improved, which means NDVI and precipitation have a close reciprocal relationship [50,51]. The winter precipitation was mainly affected by PW and PET. In winter, the source of water vapor was affected by atmospheric circulation, and the occurrence of precipitation was also greatly affected by climate.

5.2. The Comparison between Our Findings and Those of the Previous Studies

To analyze the change in precipitation before and after the implementation of the Grain for Green Project (GFGP) in the Loess Plateau, we found that precipitation had increased significantly after the GFGP. This conclusion is consistent with previous research with respect to climatic variables. For example, Zhao et al. [10] analyzed the spatiotemporal evolution of drought in the Loess Plateau from 1998 to 2014 based on TRMM (The Tropical Rainfall Measuring Mission) multi-satellite precipitation data. The results showed that precipitation on the Loess Plateau increased significantly from 1998 to 2014, and that drought severity in the Loess Plateau was a declining trend; Gao et al. [2] pointed out that the precipitation increased significantly from 1990 to 2014 in the Loess Plateau; Cheng et al. [14] analyzed the variability in average precipitation from 77 meteorological stations during 1960–2011, and found that the precipitation in the Loess Plateau had a decreasing trend before 2000, but an increasing trend after 2000. However, Tang et al. [52] made an analysis of the spatial-temporal changes of precipitation on the Loess Plateau based on 170 meteorological data. The result pointed out that the precipitation showed a decreasing trend during 1965–2014. The conflict between this result and our study may be due to the different data sources of precipitation, and more importantly, the different research period selected; Tang ignored the increasing trend of precipitation after the implementation of the GFGP.

In this study, we found that vegetation restoration had a great influence on precipitation after the GFGP, and speculated that vegetation restoration improves the local microclimate and promotes precipitation; the previous related studies have similar conclusions to ours. For example, Gao et al. [2] found that during 1990–2014, NDVI, precipitation, and AET have increasing trends; moreover, AET has a similar distribution to precipitation. Regional precipitation likely contributes to an increase in AET, and the AET increase may also influence the precipitation. Li et al. [9] evaluated the changes in NDVI in the Loess Plateau around the past 30 years and concluded that the study area became warmer, more humid, and in which drought decreased after the implementation of the GFGP. NDVI in the southeast of the Loess Plateau was mainly affected by precipitation, and the change in NDVI also had a feedback effect on regional precipitation. Wang et al. [48] studied the spatiotemporal changes of soil moisture content (SMC) and vegetation in the Loess Plateau, as well as the impact of precipitation on soil moisture content (SMC) and vegetation coverage. The results showed a significant correlation between annual precipitation and vegetation coverage during 2000 to 2015 in the northeast and west of the Loess Plateau. In summary, this study quantitatively analyzed the factors affecting precipitation on the Loess Plateau, and the results are basically consistent with previous findings.

6. Conclusions

Since the implementation of the Grain for Green Project (GFGP) in 1999, precipitation has increased significantly and the region features a warmer and wetter climate. In order to quantitatively analyze the impact of vegetation restoration on precipitation, this study mainly analyzed the temporal and spatial variation of precipitation and the potential influencing factors before and after the GFGP in the Loess Plateau. This study compared the influence degree of NDVI to precipitation before and after the GFGP, found that vegetation restoration had a great influence on precipitation after the GFGP, and speculated that vegetation recovery may have a positive feedback effect on regional precipitation. The main findings are concluded as follows:

(1) Precipitation in the Loess Plateau increased significantly after the GFGP. Before the GFGP, precipitation at 72.6% of all stations was in a decreasing trend, while 84.5% of stations showed an increasing trend after the GFGP and in the central and north stations at a 90% confidence level.

(2) With the implementation of the GFGP and climate change, the annual precipitation (PRE_A) increased at a rate of 4.96 mm/a, and the summer precipitation (PRE_S) increased at a rate of 5.63 mm/a. The PRE_A and PRE_S both showed increasing trends with a 95% confidence level. The spatial distribution of precipitation decreased from the southeast to northwest. After the GFGP, the PRE_A and PRE_S in the southeast increased significantly by around 20–50 mm, while the change in the winter precipitation (PRE_W) was not significant.

(3) The influence degree of NDVI on regional precipitation increased after the GFGP. Specifically, the annual precipitation (PRE_A) and the summer precipitation (PRE_S) were mainly influenced by NDVI, with relative influence degrees of 30.00% and 30.62%, respectively, and the PW was the second influencing factor. The winter precipitation (PRE_W) was more influenced by PW and PET, with relative influence degrees of 33.03% and 27.80%, respectively. According to the above results, we speculate that the trend of warming and wetting of the Loess Plateau is not only closely related to global climate change on the one hand, but also significantly affected by the changes in the local climate brought by vegetation restoration on the other.

(4) In this study, the implementation of the Grain for Green Project (GFGP) was used as the time node to analyze the contribution rates of five main factors to precipitation changes. However, statistical methods were used to analyze the feedback of vegetation restoration to precipitation, which makes the interaction mechanism between them difficult to explain. Therefore, it should be necessary to further study the ecological and hydrological processes in this region.