Spatial–Temporal Variation Characteristics of Multiple Meteorological Variables and Vegetation over the Loess Plateau Region

: Ecological restoration and climate change in the Loess Plateau region have become research hotspots. Climate change and anthropogenic activities have led to spatial–temporal pattern variations in vegetation and extreme climatic indices and meteorological factors. Therefore, obtaining a better understanding is necessary of the internal relations between vegetation and meteorological factors. In this paper, the interplay between vegetation index and various factors, including extreme climatic indices and meteorological factors, during a long ‐ term time series is investigated using Mann–Kendall trend analysis, and Pearson correlation coefficient analysis. The mechanisms of interaction between vegetation growth and various factors in the Loess Plateau are then analyzed. Results reveal that (i) the rapid growth of vegetation during 2000–2015 has made a major contribution to the growth trend of the Loess Plateau in the past 33 years (1982–2015). During 2000–2015, the increase of vegetation may inhibit the increase of extreme warm index and the decrease of extreme cold index; (ii) a warm and dry climate developed with decreasing relative humidity and increasing temperature; (iii) the normalized vegetation index (NDVI) is strongly correlated with extreme climatic indices and meteorological factors, especially precipitable water vapor (PWV), with a correlation coefficient of 0.94; and (iv) the daily temperature range, diurnal temperature range and sunshine duration (SSD) exerted different time ‐ delay effects on vegetation growth in the Loess Plateau. The above findings provide an essential theoretical basis for ecological policy formulation in the Loess Plateau.


Introduction
The fifth assessment of the Intergovernmental Panel on Climate Change (IPCC) reported that the global average surface temperature increased by 0.78 °C over the past 10 years (2003-2012). Anthropogenic activities are the main cause of global warming [1]. Climate change has also been closely linked to extreme weather events, such as extremely high or low temperatures, heavy precipitation, hurricanes, and persistent droughts [2]. In recent years, these extreme events have occurred with increasing frequencies, ranges, and intensities and thus exert a profound impact on the ecological environment and social production of various countries.
The continued increase in the frequency and intensity of extreme temperatures and precipitation events have become a research focus. For example, the summer heat wave in 2003 caused severe drought in the European continent [3]. Omondi et al. [4] found that the number of cold nights decreased in the African triangle region from 1961 to 2010. Pepler et al. [5] revealed severe flooding after heavy rains in the eastern coast of Australia during the summer of 2012. Yilmaz et al. [6] predicted that the intensity of precipitation would increase by 23% relative to the current condition from 2080 to 2099. Extreme disastrous events, such as extreme precipitation, high temperature and floods, are mainly caused by changes in the global climate [7][8][9]. Such conditions have caused severe damages to the global ecosystem and are a threat to human survival.
The Normalized Difference Vegetation Index (NDVI) can accurately reflect surface vegetation coverage. It can be monitored using satellite-derived NDVI to track the impact of global climate change on the ecological environment [10][11]. Zhou and Zhang [12] found that different climate conditions correspond to different distributions of vegetation types. Some studies have also investigated the relationship between extreme weather and vegetation using NDVI. For example, Rao et al. [9] and Gornall et al. [13] found that the frequent occurrence of extreme events exerts a negative effect on vegetation growth and leads to a decrease in vegetation coverage over the coastal areas of China. Bokhorst et al. [14] revealed that frequent extreme temperature events during winter seriously affect vegetation growth over the arctic subregion. The impact of extreme rainfall on vegetation coverage over the Qinghai-Tibet Plateau is less severe than the influence of average rainfall [15]. Over central Asia, the annual maximum and minimum air temperatures significantly increased between 1958 and 2012, and increased winter temperatures contributed the most to the annual warming, thereby increasing the growing season length (GSL) and productivity [16]. However, only a few extreme climate factors are considered during the analysis of their correlations with the NDVI [17] and the relationships between NDVI and extreme climatic indices are rarely discussed. Therefore, carefully investigating the mechanisms of extreme weather events in response to climate change is necessary. However, whether a time lag effect on vegetation growth caused by meteorological factors (temperature, rainfall, etc.) exists remains unclear. In addition, the influence of time delay on vegetation growth caused by extreme climatic indices is not considered when analyzing the correlations between NDVI and these meteorological factors [17].
The Loess Plateau, which is located in Northwestern China, is a typical ecologically vulnerable region. [18]. Only a few studies have addressed the relationship between extreme weather and vegetation conditions. Sun et al. [19] studied the interannual change trends of extreme climatic factors in the Loess Plateau. Zhao et al. [17] analyzed the variation of 12 extreme climatic indices and NDVI in the Loess Plateau and found that extreme rainfall and temperature occurrences are significantly correlated with different vegetation cover types. Most previous studies use the average NDVI of the entire Loess Plateau to reflect the overall change [4,6,17,20], discounting spatial variations within the region, or temporal variation over time. Moreover, climate change is not only a simple vegetation state change caused by temperature change. As an important component of the lower atmosphere, precipitable water vapor (PWV) changes dramatically with variations in temperature. PWV is among the conditions necessary for precipitation weather evolution and an important indicator influencing the occurrence of extreme weather [21][22][23][24]. Several studies have shown that the true rate of nowcasting precipitation forecasting with PWV can reach 85% [24], which reveals the strong correlation between precipitation and PWV. According to the Clausius-Clapeyron equation, a high temperature corresponds to a strong atmospheric storage capacity [25], which reveals the close relationship between temperature and atmospheric water vapor. Manandhar et al. [26] also proved that the relative humidity (RHU) and SSD are two primary factors influencing precipitation.
In this work, we investigated the potential relationship of vegetation growth with extreme weather events and meteorological parameters (PWV, T, RHU, and SSD). Variations in the relationship at different periods were considered to better understand the time-varying characteristics of the multiple meteorological factors and extreme climatic indices and vegetation over Loess Plateau. In addition, we evaluated the time-lag effect of extreme weather changes and meteorological factors on vegetation growth in detail. We aim to further understand the extreme weather change impact on the ecological environment in the Loess Plateau.

Study Area and Data Used
The Loess Plateau extends from 33° 41′ N to 41° 16′ N in latitude and from 100° 52′ E to 114° 33′ E in longitude. Given its low vegetation coverage, the ecology of Loess Plateau is very fragile, and the soil erosion is severe [27]. The Loess Plateau features a temperate continental monsoon climate, wherein rainfall is concentrated in July-September. In this paper, the Loess Plateau is divided into 10 subregions (Figure 1) based on differences in soil type and topography (Table 1). This geographical zoning map of the Loess Plateau was acquired from the Loess Plateau Scientific Data Center (http://loess.data.ac.cn). Considering that areas IX and X do not include the selected meteorological stations (right panel of Figure 1), only the eight remaining areas are discussed in this paper.   Table 2). The input meteorological parameters of 83 stations were obtained from the China Meteorological Administration (http://data.cma.cn/). PWV can be derived through many techniques, such as sounding measurement, very long baseline interferometry, microwave radiation, global navigation satellite system, and reanalysis data [28]. In this study, we used PWV data from the European Centre for Medium-Range Weather Forecasts (ECMWF) reanalysis (https://www.ecmwf.int/). This product assimilates multi-source data, such as radiosonde data, satellite radiation data, and satellite altimetry data [29]. The temporal and spatial resolutions of the PWV product used were four times a day and 0.125° × 0.125°, respectively. The PWV values of the meteorological stations were interpolated using the values of the four surrounding grid points derived from ECMWF data.
A long-term NDVI time series for the study region was obtained from the Global Inventory Modeling and Mapping Studies (GIMMS). We downloaded the 15-day composite NDVI3g dataset with a spatial resolution of 8 km over the 35 years of 1981-2015 (https://ecocast.arc.nasa.gov/data/pub/gimms/3g.v1). This dataset has been processed such that radiation and geometric correction were applied to the World Geodetic System 1984(WGS84) coordinate system [30]. The maximum value composite method was adopted to obtain the monthly NDVI time series over Loess Plateau. Then, meteorological parameters could be obtained through interpolation, including PWV, temperature, RHU, and SSD, at the location of the weather stations.

Methodology
To better understand changes in extreme weather and vegetation growth in the Loess Plateau from 1982 to 2015 (period S), we divided the long-time series into two segments using the year 2000 as the cut-off point. This year was selected because China's national policy of "returning farmland to the forest" began in 2000 in the Loess Plateau. Therefore, two sub-periods of 1982-1999 (S1) and 2001-2015 (S2) were separately considered to ensure that specific change patterns in vegetation growth under extreme weather and human intervention in these corresponding periods would not be overridden by long-term trends.
To understand the change trends of the NDVI during the experimental period, the average value, linear trend, and F-test of the NDVI were calculated over the three selected periods. Here, the annual change rate in NDVI was calculated; a positive value was considered to represent an increasing trend of vegetation, whereas a negative value represents a decreasing trend. Three levels of significance were evaluated for the F-test, i.e., P = 0.1, P = 0.05, and P = 0.01. In addition, the percentages of pixels with decreasing and increasing trends for different significance levels were counted.
To analyze the spatiotemporal distribution of PWV, RHU, and SSD in the Loess Plateau over the experimental period, several processes were performed. The 5-year sliding averages of PWV, RHU, and SSD were introduced to directly reflect the change trends of the corresponding time series, and the spatial distribution of the average PWV, RHU, and SSD was calculated over the period of 1982-2015.
To better understand changes in extreme climate indices in the Loess Plateau, the first derivatives of 7 extreme precipitation indices and 13 extreme temperature indices were calculated over the three periods of S, S1, and S2. In addition, the year in which an abrupt change in trend occurred was also calculated using the Mann-Kendall trend and abrupt analysis method. The Mann-Kendall test is a nonparametric method that does not require a certain distribution and is minimally affected by outliers [31].
The Pearson correlation coefficient was introduced to calculate correlations between NDVI and the meteorological factors (i.e., PWV, RHU, and SSD) observed through the weather stations. The formula of calculating Pearson correlation coefficient is give as follows: where r is the correlation coefficient, and r>0 is a positive correlation, r<0 is a negative correlation, the greater the absolute value of r, the stronger the correlation. and are variables, is the number of samples, and are the average values of and , respectively. Considering that the relationship between precipitation and NDVI has been widely investigated, we do not discuss precipitation in this paper. A comprehensive analysis between extreme climatic indices and NDVI was performed to explain the response of vegetation growth to climate change. The F-test was also performed to determine the significance level of vegetation change in the different subregions of the Loess Plateau during the three specified periods.
To better understand the time delay impacts of extreme climatic indices on vegetation growth, the monthly time series of the factors TXn, TNn, TXx, TNx, mean of max temperature (TMAXmean) mean of min temperature (TMINmean), RX1day, RX5day, and diurnal temperature range (DTR) were calculated using the RClimDex model. Their Pearson correlation coefficients were also analyzed along with the corresponding NDVIs for the current month and after 1-3 months during the growing season (April-October) in the Loess Plateau. In this paper, the correlation degree was defined based on previous studies on the relationship between climatic indices and ecological factors [32][33][34]. Table  3 shows the level of correlation.

Vegetation Distribution over the Loess Plateau
The distributions of the average values, linear trends, and F-test results of NDVI over the three selected periods are presented in Figure 2. Figures 2(a)-2(c) show that the vegetation coverage generally decreases from southeast to northwest with a relatively low NDVI of 0-0.4 in the northwestern region. The overall vegetation coverage rate improved during S compared with that during S2, which suggests that the major contributor to the vegetation growth in S is the vegetation growth of S1 (Figure 2(d)-2(f)). The statistical results reveal that the average trends over Loess Plateau during S, S1, and S2 were 0.022/10a, 0.013/10a, and 0.059/10a, respectively. After 2000, the growth rate of vegetation over Loess Plateau was 4.5 times that before 2000, In addition, [35] shows that the vegetation change rate of the whole China during 1982-2012 is 0.002/10a, and the vegetation growth rate of Loess Plateau during S and S2 is 11 and 29.5 times than that of the whole China. From the spatial distribution of NDVI, vegetation decreased in the southern and northern portions of Loess Plateau, especially in area I. Figure 2(g)-2(i) show the F-test distributions vegetation change in the three periods of S, S1, and S2 in the Loess Plateau; Table 4 provides the statistical results of the F-test showing decreasing and increasing trends for different significance levels. Approximately 56.76% of the pixels revealed an increasing trend during S2, while 72.23% of the pixels showed an increasing trend during S (P < 0.01).  According to the spatial distribution of vegetation in the Loess Plateau (Figure 1), in the north and south (area I), the vegetation decreased significantly. The main reason is that there are cities of Zhong, Yinchuan, Shizuishan, Wuhai, Linhai, Baotou, Hohhot and others in the north of area I. These cities have sufficient sunshine, strong evaporation, short summer, and late spring, and are prone to extreme days such as cold waves and sandstorms. In recent years, with the acceleration of urbanization, the area of construction land has increased, and the change of land use has led to the obvious degradation of vegetation in this area in the past 30 years. For the south of zone I, there are 12 cities, such as Shuozhou, Linfen and Yuncheng. The open-pit coal mining in these areas all yearround has caused serious damage to the surface vegetation in this area and is not easy to repair. These reasons directly lead to the significant degradation of vegetation in zone I of the Loess Plateau. Figure 3 presents the average PWV time series with a 5-year sliding average and its spatial distribution from 1982 to 2015. Figure 3(a) shows that the PWV did not undergo evident variations in the Loess Plateau over the past 33 years, and values ranged from 12 mm to 14 mm. However, the spatial distribution of PWV in Figure 3(b) reveals that the atmospheric water decreased from southeast to northwest with values varying from 6 mm to 23 mm, which is partially caused by the topography of Loess Plateau (Figure 4). In addition, an increase in the continental effect toward the northwest increased the dryness of the climate. Figure 4 shows that the variation of elevation is consistent with that of PWV, especially for the eastern portion of Qinghai Province (marked by a red rectangle in Figure 3), which features an elevation of approximately 2250 m above sea level. However, the average PWV was only approximately 6 mm.

RHU and SSD
The time series of the annual average RHU and SSD with the corresponding 5-year sliding averages and spatial distributions are presented in Figures 5 and 6, respectively. Figures 5(a) and 6(a) show that RHU and SSD presented a slight downward trend in the Loess Plateau over the past 34 years. RHU varied from 54% to 64%, whereas SSD fluctuated from 6.3 h to 7.3 h. RHU was relatively high north and west of Loess Plateau, but low in the cities of Wuzhong, Taiyuan, and Changzhi ( Figure 1). SSD was only approximately 4-5 h in the north of Loess Plateau but increased to 7-9 h in the middle of this region. This phenomenon may be associated with the latitudes and topography.

RHU and SSD further showed a significant negative correlation in Figures 5(b) and 6(b) in the Loess
Plateau. For areas with a long SSD, the corresponding RHU was low, whereas areas with a short SSD had a high RHU.

Change in Extreme Climatic Indices over the Loess Plateau Region
The time series of the highest temperature of the daily maximum per month (TXx) is presented in Figure 7, where TXx showed an overall upward trend with increasing global temperature during S. However, while the slope of TXx during the period S1 was positive, the slope in period S2 was negative. This finding suggests that the increase in TXx during S is mainly contributed by that during S1. Compared with the vegetation variation in Figure 2, the vegetation coverage is obviously poor during S1. The statistical results revealed that an increase in NDVI by 0.0015 is associated with a decrease the TXx by 0.0327 °C during S2 period in the Loess Plateau. The trends of extreme climate indices during S, S1, and S2 in the Loess Plateau are presented in Figure 8. All first derivatives of the 7 extreme precipitation indices were higher than 0 and showed an upward trend over the past 34 years (period S). However, the first derivatives of RX1day, RX5day, R95p, and R99p were less than 0 and demonstrated a downward trend over the past 15 years (period S2). The first derivatives of R10, R20, and R25 showed trends completely contrasting those during S1 (less than 0) and S2 (larger than 0). Therefore, the upward trends of these three indices during S2 greatly affected the trends during S, and explain the frequent occurrence of extreme precipitation events during S2.
In terms of extreme temperature cold indices, the first derivatives of TX10P, FD0, and CSDI were consistently less than 0 during the three periods under study, but more so in S1 than in S2. The extreme temperature cold indices show a downward trend during the S, S1 and S2 periods, and the decreased rates are S1 > S > S2, which suggests that the downward trend during S is mainly affected by the discredited trend of period S1. In terms of extreme temperature warm indices, the variation characteristics of TNx were similar to those of TXx. The first derivatives of TX90P, SU25, GSL, TMAXmean, and TMINmean were higher than 0 during S. However, the values derived from S1 were higher than those from S, whereas the values obtained from S2 were less than those from S,which suggests that the upward trend during S is mainly affected by the increased trend of period S1. TXx revealed an abrupt increase at around 1994 (P < 0.01) (Figure 7(b)). The statistical results of Figure 7 reveal that the multi-year average value of TXx over the entire Loess Plateau was 34.36 °C. An extremely high value of 36.37 °C occurred in 2010, while an extremely low value of 32.56 °C was observed in 1989. The change rate of TXx was 0.48°C/10a, and the statistical result of other extreme climate indices are provided in Table 5.

NDVI and Meteorological Factors
Considering that the PWV is highly related to temperature and precipitation, its variation with NDVI was also analyzed. The Pearson correlation coefficients of NDVI with PWV, RHU, and SSD were calculated (Figure 9). Strong correlations existed between NDVI and PWV with a correlation value of 0.94. Many studies have shown that vegetation growth is highly influenced by precipitation [36][37][38][39], and precipitation mainly originates from atmospheric water vapor [24]; thus, a strong correlation exists between NDVI and PWV. Most parts of Loess Plateau are semi-arid regions, and a large part of the moisture in the near-surface air comes from the evapotranspiration of the vegetated land surface. Thus, a moderate correlation of 0.52 between NDVI and RHU was observed. SSD is among the most influential factors contributing to the vegetation photosynthetic active radiation [40], and a good correlation of 0.464 was found between NDVI and SSD.

Correlations between Meteorological Factors and Extreme Climatic Indices
The Pearson correlation coefficients of PWV with extreme temperature indices are provided in Figure 10(a). A strong correlation of above 0.78 was found between PWV and extreme temperature indices, such as TMAXmean, TMINmean, TNn, TXx, TNx, and TXn. The Pearson correlation coefficients of the extreme temperature indices showed a decreasing trend from southwest to northeast, which is possibly associated with the PWV distribution. Xu et al. [25] found that a high temperature corresponds to high atmospheric moisture content. Figure 11 shows the correlation coefficient distribution statistics of extreme temperature index, extreme rainfall index and PWV in the Loess Plateau; five different colors in the figure represent different statistical intervals of correlation coefficient. It can be found that the coefficients of the correlations of PWV with TMAXmean, TMINmean, TXx, TNx, TXn, and TNn were higher than 0.4. Weak correlations existed between PWV and other extreme temperature indices, such as FD0, SU25, WSDI, CSDI, and GSL, accounting for 55.07%, 45.65%, 58.70%, 44.20%, and 63.77% of Loess Plateau, respectively. The increase in extreme cold indices could accelerate the transformation of atmospheric water vapor into the rain, further causing a decrease in atmospheric PWV. Therefore, the extreme cold indices FD0 and CSDI showed negative relationships with PWV, accounting for 42% and 59.4% of Loess Plateau, respectively.  The statistical results of the Pearson correlation coefficients of extreme climate indices with RHU and SSD are also presented in Figures 12 and 13, respectively. Figure 12 shows that the Pearson correlation coefficients of the relationships of RHU with TMAXmean, TMINmean, TXx, TNx, TXn, and TNn are higher than 0.2, accounting for 67.39%, 80.43%, 60.15%, 68.37%, 72.74%, and 82.6% of Loess Plateau, respectively. Extreme temperature cold indices did not show evident correlations with RHU, whereas negative correlations existed between RHU and extreme temperature warm indices (SU25, WSDI, and GSL). Some extreme precipitation indices, such as RX1day and RX5day, showed good correlations with RHU (0.2-1), accounting for 88.40% and 92.75% of Loess Plateau, respectively. Extreme temperature indices, such as TMAXmean, TMINmean, TXx, TNx, TXn, and TNn, also showed positive relationships with SSD with Pearson correlation coefficients ranging from 0.2-1, accounting for 88.5%, 84.05%, 88.41%, 84.05%, 87.68%, and 84.78% of Loess Plateau, respectively.

NDVI and Extreme Climate Indices
Comparisons revealed that extreme temperature indices, such as TX10P, TN10P, TX90P, TN90P, R95p, R99p, R10, R20, R25, CWD, CDD, PRCPTOT, and SDII, were not correlated with NDVI in the Loess Plateau. Figures 14 and 15 respectively present the distribution of Pearson correlation coefficients between NDVI and extreme temperature indices and their percentages. The percentages of the Pearson correlation coefficients of NDVI with cold indices FD0 and ID0 ranged from −1 to −0.2 and accounted for 56.33% and 40.58% of Loess Plateau, respectively, while that of NDVI with warm indices SU25 and TR20 ranged from 0.2 to 1 and accounted for 62.31% and 45.46% of Loess Plateau, respectively. This finding is indirectly supported by analyzing the Pearson correlation coefficients of NDVI with DTR. The percentages of the Pearson correlation coefficients of NDVI with WSDI and GSL ranged from 0.2 to 1 and accounted for 47.58% and 52.17% of the variation in the Loess Plateau, respectively. Statistical results reveal that the Pearson correlation coefficients between NDVI and TXn, TNn, TXx, TNx, TMAXmean, TMINmean, RX1day, and RX5day are very high, with values ranging from 0.4-1 and accounting for 95% of Loess Plateau; by contrast, the percentages of Pearson correlation coefficients between NDVI and DTR ranged from −1 to 0.2 and accounted for 34.78% of Loess Plateau.

Time Lag Effect
The Pearson correlation coefficients between the extreme climatic indices factors (TXn, TNn, TXx, TNx, TMAXmean, TMINmean, RX1day, RX5day, and DTR) and the corresponding NDVI in the current month and after 1−3 months are provided in Figure 16 for the growing season (April−October) in the Loess Plateau. While the effect of the time delay of some factors, such as TXn, TNn, TXx, TNx, TMAXmean, and TMINmean, on vegetation growth was 0 month, that of DTR showed various time lengths in different areas of Loess Plateau ( Table 6). The time lengths of time lag effects of DTR on vegetation growth differed for various vegetation types in the Loess Plateau.  The distribution of time delay impacts on vegetation growth caused by meteorological factors in the Loess Plateau (PWV, RHU, and SSD) is presented in Figure 17. Among other factors, the correlation between PWV and NDVI was the highest when the time delay was 0, accounting for 98.55% of Loess Plateau. The time delay impacts on vegetation growth associated with RHU and SSD differed. Statistical results reveal that only 66.67% of Loess Plateau shows maximum Pearson correlation coefficients between NDVI and RHU when the time delay is 0. The time lag effect of SSD on vegetation growth was distributed from 0 to 3 months, and the response of vegetation to SSD in the northern part of the Loess Plateau was mainly two months, accounting for 30.43% of Loess Plateau. In the middle and south of Loess Plateau, the lag response was mainly 3 months, accounting for 34.78% of the total Loess Plateau area. SSD has no lag effect on 15.22% of the regional vegetation growth in the western region of the Loess Plateau and 19.57% of the vegetation growth in the northeast region.

Discussion
The period 1982-2015 was divided into two shorter sub-periods in this paper to highlight the variation features of various factors within each sub-period. For example, Figure 2 suggests that the major contribution to the vegetation growth in period S comes from period S2 and that implementation of the "returning farmland to the forest" policy effectively promoted vegetation growth and ecological construction in the Loess Plateau. In addition, reductions in vegetation in area I of Loess Plateau are mainly associated with some cities located in the northern portion of Loess Plateau (Figure 1), where urbanization is accelerated, and the construction land has increased. Therefore, changes in land use caused an obvious decrease in vegetation in this region over the past 34 years. In the southern portion of Loess Plateau in area I, a perennial open-pit coal mining has caused serious damage to the surface vegetation. Zhao et al. [17] showed that the average growth rate of vegetation in the whole Loess Plateau during the period of 1982-2013 is 0.025/10a. In this paper, the vegetation change rate of the Loess Plateau from 1982 to 2015 is 0.022/10a, which is similar to Zhao et al. [17]. In addition, the vegetation growth of the Loess Plateau in S2 period is far greater than that of S1 and the vegetation growth of the whole region of China is mainly related to human activities, such as the implementation of the policy of ʺreturning farmland to forestʺ.
Several studies have been carried out to analyze extreme climate indices during the period of 1982-2013 [19,10]. Figure 7 shows that extreme temperature indices (TXx) increased rapidly during S and S1, and the values are larger in S1 period than that in S period, but TXx revealed a downward trend during S2. In addition, the order of vegeration increased rate in the Loess Plateau in three periods are S2>S>S1. Such a phenomenon may suggests that the increase in TXx is may hindered by improved vegetation conditions. Comparison of the correlations between extreme climatic indices and vegetation index during S, S1, and S2 showed that the decrease rate of extreme temperature cold indices and increase rate of extreme temperature warm indices slowed due to the continuous increase in vegetation cover in the Loess Plateau. We found that the slopes of TX10P, FD0, and CSDI are downward during S2 period but with a slowing rate; these values were lower than those from S and S1. By comparing the Figure 2c,f,j, we found that the vegetation cover in the Loess Plateau increased significantly during S2, and the decrease rate of cold index in the Loess Plateau during S2 period may be associated with the increase of vegetation. These findings also revealed that the upward slower trends of the warm indices of extreme temperature (such as TX90P, SU25, GSL, TMAXmean, and TMINmean) and the increase in vegetation coverage during S2 over the Loess Plateau. This finding may suggest that vegetation exerts a negative impact on extreme event occurrences, which is rarely reported.
Additionally, spatial-temporal analysis of PWV is presented in this paper, which has not been previously investigated, in contrast to precipitation, which has been widely analyzed in the Loess Plateau [32,[41][42][43]. PWV revealed good consistency with the elevation distribution of Loess Plateau, agreeing with the results of Gui et al. [28] and Zhao et al. [44]. Our study suggests that the increase in temperature in the Loess Plateau over the past 34 years has caused the climate to become drier, as indicated by the reduced RHU observed.
Zhao et al. [17] showed good correlations between NDVI with extreme climate indices (P < 0.01), such as TXn, TNn, TXx, TNx, TMAXmean, TMINmean, RX1day, and RX5day, a similar conclusion is obtained in this paper. In addition, more extreme climate factors, and time-lag effects are discussed. Among the 9 indices of extreme climate factors, DTR has a lag effect on vegetation growth from 0 to 3 months, and the other 8 extreme climate indexes have an impact on NDVI in the same month. Analysis of the delay effect of meteorological factors and NDVI, it is found that SSD had 0-3 months of lag response to vegetation while the lag effect of PWV and RHU on vegetation is mainly 1 month in the Loess Plateau. The main reasons for this phenomenon may be the characteristics of topography and rapid growth of vegetation in this area; in addition, the level of urbanization may also be another important reason for the change of climate and meteorological factors. Changes in the growing season can affect the exchange and transfer of energy through the biosphere and atmosphere, as well as the global carbon and energy cycles [45]. Therefore, variations in GSL were investigated, and GSL was found to be extended over the past 34 years, which agrees with the conclusions of Liu et al. [4] and Piao et al. [45] that the GSL in Northern China is prolonged. The extreme climatic and meteorological factors change, and their responses to vegetation growth in the recent 33 years may be associated with the vegetation growth rapidly and topography in the Loess Plateau. The surface vegetation cover change and urbanization recent 33 years may be one of the important reasons for climate change in the Loess Plateau during the S2 period [46][47][48]. Additionally, only 53 and 68 meteorological stations in the Loess Plateau were used by Zhao et al. [17] and Sun et al. [19], whereas 83 meteorological stations were selected in this paper. More evenly distributed meteorological stations can better reflect the actual situation of the entire Loess Plateau. Hence, the results obtained in this paper may slightly differ from those of previous studies.

Conclusions
This paper analyzed the relationships of the spatial-temporal variation characteristics of multiple factors with NDVI in the Loess Plateau during S1, S2, and S. Some interesting results were obtained: (i) the growth rate of vegetation coverage increased in the Loess Plateau over the past 15 years due to reforestation, which is approximately 4.5 times higher than that from 1982 to 2000 and 29.5 times than that in the whole of China over the past 34 years; (ii) anthropogenic activities, such as accelerated urbanization and excessive mining easily caused vegetation degradation, especially in area I of Loess Plateau; (iii) The increase of vegetation in the Loess Plateau is associated with the change of extreme weather events and was embodied by a slowing of the decrease rate of extreme temperature cold indices as well as the slowing of increase rate of extreme temperature warm indices; (iv) the regional climate change shifted toward warm and dry over Loess Plateau; and (v) atmospheric water vapor or PWV was strongly correlated with NDVI over Loess Plateau, as confirmed by Pearson correlation coefficients reaching 0.94.
Although some novel findings were discovered in this paper, some results involve qualitative data, and quantitative results could not be obtained because of some limitations. For example, vegetation growth revealed an inhibitory effect on the change in extreme climatic indices, but the specific contribution rate of vegetation growth to the inhibitory effects of extreme climatic indices was not determined. In future research, such problems may be resolved by establishing the corresponding quantitative models.
Author Contributions: Q.Z. and X.M. participated in the design of this study, and they both performed the statistical analysis. Q.Z, L.L., and W.Y. carried out the study and collected important background information. Q.Z. and X.M. processed the corresponding data. All authors have read and agreed to the published version of the manuscript.