Spatial–Temporal Evolution of Vegetation NDVI in Association with Climatic, Environmental and Anthropogenic Factors in the Loess Plateau, China during 2000–2015: Quantitative Analysis Based on Geographical Detector Model

: In the Loess Plateau (LP) of China, the vegetation degradation and soil erosion problems have been shown to be curbed after the implementation of the Grain for Green program. In this study, the LP is divided into the northwestern semi-arid area and the southeastern semi-humid area using the 400 mm isohyet. The spatial–temporal evolution of the vegetation NDVI during 2000–2015 are analyzed, and the driving forces (including factors of climate, environment, and human activities) of the evolution are quantitatively identiﬁed using the geographical detector model (GDM). The results showed that the annual mean NDVI in the entire LP was 0.529, and it decreased from the semi-humid area (0.619) to the semi-arid area (0.346). The mean value of the coefﬁcient of variation of the NDVI was 0.1406, and it increased from the semi-humid area (0.1165) to the semi-arid area (0.1926). The annual NDVI growth rate in the entire LP was 0.0079, with the NDVI growing faster in the semi-humid area (0.0093) than in the semi-arid area (0.0049). The largest increments of the NDVI were from grassland, farmland, and woodland. The GDM results revealed that changes in the spatial distribution of the NDVI could be primarily explained by the climatic and environmental factors in the semi-arid area, such as precipitation, soil type, and vegetation type, while the changes were mainly explained by the anthropogenic factors in the semi-humid area, such as the GDP density, land-use type, and population density. The interactive analysis showed that interactions between factors strengthened the impacts on the vegetation change compared with an individual factor. Furthermore, the ranges/types of factors suitable for vegetation growth were determined. The conclusions of this study have important implications for the formulation and implementation of ecological conservation and restoration strategies in different regions of the LP.


Introduction
Vegetation plays an indispensable role in regional terrestrial ecosystems, and constitutes an essential link between soil, water, and atmosphere. Moreover, vegetation is the Qin River, Fen River, etc. The study area includes the provinces (or autonomous regions) of Qinghai, Gansu, Ningxia, Inner Mongolia, Shaanxi, Shanxi, and Henan. Throughout history, the overall vegetation coverage has been low and soil erosion has been very serious on the LP due to the special natural and geographical environment, complex vegetation types, and severe impacts of human activities. Since the 1980s, a series of measures for soil and water conservation have been taken to control the soil erosion. Particularly, after 2000, when the GGP was implemented, the regional ecological environment has significantly improved.
Remote Sens. 2021, 13,4380 4 of 30 the soil erosion, floods, landslides, debris flows and other disasters. The total annual radiation is 50.2 × 10 4 J/cm 2 to 67.0 × 10 4 J/cm 2 , with a long illumination time and high radiation. About 200 rivers have their headwaters on the LP, including the Tao River, Zuli River, Qingshui River, Huangfuchuan River, Kuye River, Wuding River, Beiluo River Wei River, Qin River, Fen River, etc. The study area includes the provinces (or autonomous regions) of Qinghai, Gansu, Ningxia, Inner Mongolia, Shaanxi, Shanxi, and Henan. Throughout history, the overall vegetation coverage has been low and soil erosion has been very serious on the LP due to the special natural and geographical environment, complex vegetation types, and severe impacts of human activities. Since the 1980s, a series of measures for soil and water conservation have been taken to control the soil erosion. Particularly, after 2000, when the GGP was implemented, the regional ecological environment has significantly improved. The study area was divided into semi-arid and semi-humid areas using the 400 mm isohyet (calculated using the annual mean precipitation data for the LP during [2000][2001][2002][2003][2004][2005][2006][2007][2008][2009][2010][2011][2012][2013][2014][2015] to carry out the comparative analysis. It should be noted that the 400 mm isohyet is of great significance to the physical geography and socio-economic regionalization in China, and it coincides with the Huhuanyong Line. In addition to dividing the climatic zones, this line divides the forest vegetation from the grassland vegetation and the densely populated area from the sparsely populated area [31,32].

Geographical Detector Model (GDM)
The GDM is proposed based on spatial differentiation theory and geographic information system (GIS) spatial analysis technique [27]. It is usually employed to study the factors affecting spatial hierarchical heterogeneities and the underlying mechanisms. The model assumes similar spatial pattern between the independent and dependent variables if these variables are spatially correlated [33]. The GDM consists of four detectors, including the factor, interaction, risk and ecological detectors. The study area was divided into semi-arid and semi-humid areas using the 400 mm isohyet (calculated using the annual mean precipitation data for the LP during [2000][2001][2002][2003][2004][2005][2006][2007][2008][2009][2010][2011][2012][2013][2014][2015] to carry out the comparative analysis. It should be noted that the 400 mm isohyet is of great significance to the physical geography and socio-economic regionalization in China, and it coincides with the Huhuanyong Line. In addition to dividing the climatic zones, this line divides the forest vegetation from the grassland vegetation and the densely populated area from the sparsely populated area [31,32].

Geographical Detector Model (GDM)
The GDM is proposed based on spatial differentiation theory and geographic information system (GIS) spatial analysis technique [27]. It is usually employed to study the factors affecting spatial hierarchical heterogeneities and the underlying mechanisms. The model assumes similar spatial pattern between the independent and dependent variables if these variables are spatially correlated [33]. The GDM consists of four detectors, including the factor, interaction, risk and ecological detectors.
(1) The factor detector can be utilized to analyze the spatial heterogeneity and to quantify the explanatory power (measured by the q value) of different impact factors X to dependent variable Y. where q ∈ [0, 1], and with the increase in the q value, the explanatory power is expected to be stronger; h represents the stratum of variable Y (or factor X), h ∈ [1, L]; N h is the unit number of stratum h; N is the unit number in all the strata; σ 2 h represents the variance of variable Y of stratum h; σ 2 is the variance of variable Y of all the strata.
(2) The interaction detector can be used to determine the explanatory power of interaction between two factors (say X 1 and X 2 ) on the spatial heterogeneity of variable Y, and to judge whether the interactive effect on variable Y would be enhanced or weakened. The steps are as follows. (1) Compute the respective q value of factors X 1 and X 2 to variable Y (q X1 and q X2 ). (2) Compute the q value of the interaction between X 1 and X 2 to Y (q X1&X2 ).
(3) Compare q X1 , q X2 , and q X1&X2 . If max(q X1 ,q X2 ) < q X1&X2 < q X1 + q X2 is true, the factors X 1 and X 2 enhance each other (bi-enhance). If q X1&X2 > q X1 + q X2 is true, the nonlinearity of two factors is enhanced (nonlinear enhancement). The interaction criteria are presented in Table 1. Table 1. Definition of interaction detector.

Description
Interaction min(q X1 ,q X2 ) < q X1&X2 < max(q X1 ,q X2 ) Weaken, uniq X1&X2 < min(q X1 ,q X2 ) Weaken, nonlinear q X1&X2 > max(q X1 ,q X2 ) Enhance, biq X1&X2 > q X 1 + q X2 Enhance, nonlinear q X1&X2 = q X1 + q X2 Independent (3) The risk detector is utilized to judge the difference of significance between the attribute mean values of two strata with the t statistic: where Y h is the attribute mean value within stratum h; n h is the number of samples within stratum h; Var is the variance.
(4) The ecological detector is developed to compare the explanatory powers of two factors, X 1 and X 2 , on variable Y with the F statistic: where N X1 and N X2 are the sample numbers of factors X 1 and X 2 , respectively; SIV X1 and SIV X2 are the sum of the internal variances of the strata from factors X 1 and X 2 , respectively; and L 1 and L 2 are the strata numbers of factors X 1 and X 2 , respectively.

Data Description
In this study, multi-source data (Table 2) were collected, mainly including the following.
(1) Vector data: The shape file data for the LP were downloaded from the Resource and Environment Science and Data Center (RESDC), Chinese Academy of Sciences (available from http://www.resdc.cn, accessed on 27 October 2021).
(2) Vegetation NDVI data: The SPOT NDVI data were used in this study (Figure 2), which were obtained from the RESDC. The time span is 2000-2015, the time step is 1 year, and the spatial resolution is 1 km × 1 km. It should be noted that the annual maximum synthesis method was used to obtain the annual value of NDVI, and therefore, there was no area with NDVI < 0.   According to [34], the NDVI of the LP was divided into 5 intervals using the equal interval method, i.e., 0 to 0.2, 0.2 to 0.4, 0.4 to 0.6, 0.6 to 0.8, and 0.8 to 1.0, which represent 5 types of vegetation, i.e., low, relatively low, medium, relatively high, and high vegetation coverages, respectively.
(3) Meteorological data: The monthly mean precipitation and temperature from 101 meteorological stations on the LP and in its surrounding areas during 2000-2015 were collected for use, which were downloaded from the China Meteorological Data Service Center (available from http://data.cma.cn/site/index.html, accessed on 27 October 2021). The annual value of the data from each meteorological station was calculated. Moreover, ANUSPLIN v4.3 was used to perform spatial interpolation of the meteorological station data, which were divided into 6 intervals using the equal interval method. Therefore, annual precipitation and temperature grid data (Figures 3 and 4) with the same projection mode, spatial resolution, and time series as the NDVI data were obtained. It should be noted that ANUSPLIN is a widely used software application for spatial interpolation, which was developed on the basis of the thin plate smoothing splines theory. The main feature of ANUSPLIN is that the topographic information can be used in the spatial interpolation of the meteorological data [35].
(4) DEM data: The DEM data were downloaded from the Shuttle Radar Topography Mission (available from http://srtm.csi.cgiar.org/, accessed on 27 October 2021). According to [28], the natural breakpoint method was used to divide the elevation data into 6 categories, including 90 to 790 m, 790 to 1228 m, 1228 to 1611 m, 1611 to 2136 m, 2136 to 2963 m, and 2963 to 4914 m (Figure 5a).
Based on the DEM data, the slope and slope aspect data for the study area were extracted using GIS spatial analysis tools. According to the policy of the Grain for Green [36], the slope was divided into 6 categories, including 0 • to 5 • , 5 • to 10 • , 10 • to 15 • , 15 • to 20 • , 20 • to 25 • , >25 • (Figure 5b).
Based on the DEM data, the slope and slope aspect data for the study area were extracted using GIS spatial analysis tools. According to the policy of the Grain for Green [36], the slope was divided into 6 categories, including 0° to 5°, 5° to 10°, 10° to 15°, 15° to 20°, 20° to 25°, >25° ( Figure 5b).
(5) Geomorphic type data: These data were obtained from the RESDC. Referring to the Geomorphological Atlas of the People's Republic of China (1:1 million) (2009), the geomorphology was divided into 6 types, including plain, platform, hill, small undulating mountain, medium undulating mountain, and large undulating mountain ( Figure 5d).
(6) Soil type data: These data were downloaded from the RESDC. Referring to the 1:1 Million Soil Map of the People's Republic of China (1995), the soil was divided into 11 types, including alpine soil, anthropogenic soil, saline alkali soil, hydrogenetic soil, semihydrogenetic soil, primary soil, desert soil, arid soil, calcareous soil, semi-eluvial soil, and eluvial soil (Figure 5e).
(7) Vegetation type data: These data were obtained from the RESDC. Referring to the 1:1 Million Vegetation Atlas of China (2001), the vegetation was divided into 8 types, including cultivated vegetation, meadow, grass, grassland, desert, shrub, broad-leaved forest, and coniferous forest (Figure 5f). (5) Geomorphic type data: These data were obtained from the RESDC. Referring to the Geomorphological Atlas of the People's Republic of China (1:1 million) (2009), the geomorphology was divided into 6 types, including plain, platform, hill, small undulating mountain, medium undulating mountain, and large undulating mountain ( Figure 5d).
(6) Soil type data: These data were downloaded from the RESDC. Referring to the 1:1 Million Soil Map of the People's Republic of China (1995), the soil was divided into 11 types, including alpine soil, anthropogenic soil, saline alkali soil, hydrogenetic soil, semi-hydrogenetic soil, primary soil, desert soil, arid soil, calcareous soil, semi-eluvial soil, and eluvial soil ( Figure 5e).
(7) Vegetation type data: These data were obtained from the RESDC. Referring to the 1:1 Million Vegetation Atlas of China (2001), the vegetation was divided into 8 types, including cultivated vegetation, meadow, grass, grassland, desert, shrub, broad-leaved forest, and coniferous forest (Figure 5f).
(8) Socio-economic data: These data were obtained from the RESDC. The GDP density and population density data for China in 2000, 2005, 2010, and 2015 were obtained. The GDP density and population density were divided into 6 categories, respectively ( Figures 6 and 7).
(9) Land-use type data: These data were downloaded from the RESDC. The spatial resolution is 1 km × 1 km. There are 6 land-use types, including unused land, construction land, water bodies, grassland, woodland, and farmland ( Figure 8), with a total of 36 transfer combinations.
The influences of the various factors on the spatial-temporal evolution of the NDVI on the LP were calculated and analyzed with the GDM. The various factors were classified into climatic, environmental, and anthropogenic aspects. For the climatic factors, precipitation, and temperature were selected as proxy variables. For the environmental factors, altitude, slope, slope aspect, geomorphic type, soil type, and vegetation type were selected as proxy variables. For the anthropogenic factors, the proxy variables were the GDP density, population density, and land-use type.        (9) Land-use type data: These data were downloaded from the RESDC. The spatial resolution is 1 km × 1 km. There are 6 land-use types, including unused land, construction land, water bodies, grassland, woodland, and farmland ( Figure 8), with a total of 36 transfer combinations. The influences of the various factors on the spatial-temporal evolution of the NDVI on the LP were calculated and analyzed with the GDM. The various factors were classified into climatic, environmental, and anthropogenic aspects. For the climatic factors, precipitation, and temperature were selected as proxy variables. For the environmental factors, altitude, slope, slope aspect, geomorphic type, soil type, and vegetation type were selected as proxy variables. For the anthropogenic factors, the proxy variables were the GDP density, population density, and land-use type.
The following points should be noted: (1) Although environmental factors do not change in the short-term, these factors play a key role in vegetation growth and changes. Thus, they were also considered in this study according to [37]. (2) The data for each factor were discretized since the input of the GDM needs to be discrete data. (3) The data for the factors were converted into grid data with the 1 km × 1 km spatial resolution (same with the NDVI data) in order to facilitate calculation and analysis. (4) Due to the large area of the LP, the GDM's calculation capacity would be exceeded if all of the grid data were used. Therefore, the sampling method was used and 32,000 sampling points were evenly selected in space to carry out the calculations and analysis (Figure 1), and the consistency of the results of each calculation could be ensured.  The following points should be noted: (1) Although environmental factors do not change in the short-term, these factors play a key role in vegetation growth and changes. Thus, they were also considered in this study according to [37]. (2) The data for each factor were discretized since the input of the GDM needs to be discrete data. (3) The data for the factors were converted into grid data with the 1 km × 1 km spatial resolution (same with the NDVI data) in order to facilitate calculation and analysis. (4) Due to the large area of the LP, the GDM's calculation capacity would be exceeded if all of the grid data were used. Therefore, the sampling method was used and 32,000 sampling points were evenly selected in space to carry out the calculations and analysis (Figure 1), and the consistency of the results of each calculation could be ensured. (5) Two types of inputs were considered in the analysis of the explanatory power of different factors for the spatial-temporal evolution of the NDVI on the LP using the GDM. First, the typical annual values of 2000, 2005, 2010, and 2015 were used to drive the GDM to analyze the influence of various factors on the spatial distribution state of the NDVI. Second, the differences between the annual mean values from 2000 to 2005 and those from 2010 to 2015 were calculated, which were used to drive the GDM to analyze the influence of various factors on the spatial distribution change of the NDVI.

Spatial-Temporal Evolution in NDVI
(1) Grid scale area and area transfer of the NDVI and the driving forces: The vegetation types on the LP are shown in Figure 5f. In the semi-arid area, the vegetation types were mainly grassland, cultivated vegetation, and desert, accounting for 15.9, 7.2 and 5.0% of the total area, respectively. While in the semi-humid area, cultivated vegetation, grassland, shrub, and broad-leaved forest were widely distributed, accounting for 37.8, 9.2, 6.7 and 5.4% of the total area, respectively.
The grid scale area transfer matrix of the NDVI on the LP during 2000-2015 is shown in Table 3. In 2000, over the entire LP, the coverage areas of the low, relatively low, medium, relatively high, and high vegetation accounted for, respectively, 10.0, 31.5, 29.3, 28.1 and 1.3% of the total area; in the semi-arid area, these types of vegetation accounted for 9.4, 16.0, 5.1, 1.2 and 0.01%, respectively; and in the semi-humid area, they accounted for 0.5, 15.5, 24.2, 26.9 and 1.3%, respectively. In 2015, over the entire LP, the coverage areas of the low, relatively low, medium, relatively high, and high vegetation accounted for 8.7, 20.3, 22.9, 31.5 and 16.7%, respectively; in the semi-arid area, they accounted for 8.5, 14.5, 5.0, 3.3 and 0.4%, respectively; and in the semi-humid area, they accounted for 0.3, 5.8, 17.9, 28.1 and 16.3%, respectively. The proportion of the coverage areas with an NDVI of greater than 0.6 increased significantly, and the proportions of the relatively low vegetation converted to medium vegetation, medium vegetation converted to relatively high vegetation, and relatively high vegetation converted to high vegetation accounted for the largest proportions. The areas of these transitions accounted for, respectively, 10.8, 15.8 and 14.1% of the total area. The land-use type on the LP during 2000-2015 is shown in Figure 8, and the grid scale area transfer matrix is shown in Table 4. In 2000, over the entire LP, the proportions of the farmland, woodland, grassland, water bodies, construction land, and unused land were 33.3, 14.7, 41.5, 1.4, 2.5 and 6.7%, respectively; in the semi-arid area, they were 7.1, 1.2, 16.6, 0.7, 0.8 and 5.3%, respectively; and in the semi-humid area, they were 26.2, 13.5, 24.9, 0.8%, 1.7 and 1.4%, respectively. In 2015, over the entire LP, the proportions of farmland, woodland, grassland, water bodies, construction land, and unused land were 32.4, 15.1, 41.2, 1.5, 3.3 and 6.6%, respectively; in the semi-arid area, they were 6.9, 1.3, 16.4, 0.7, 1.1 and 5.2%, respectively; and in the semi-humid area, they were 25.4, 13.8, 24.8, 0.8, 2.1 and 1.4%, respectively. Overall, the area of the change in land-use types on the LP was not significant, and it was mainly in the 0.94% decrease in farmland, the 0.44% increase in woodland, and the 0.82% increase in construction land. To a certain extent, these changes reflect the impacts of human activities in the study area, such as the implementation of the GGP and urbanization-related development.
The difference in the NDVI of the different land-use type transfers on the LP during 2000-2015 is shown in Figure 9. It should be noted that Figure 9 was obtained by subtracting the NDVI data for 2000 from that for 2015 on the grid scale, and the differences in the NDVI of the different land-use type transfers were then counted. Figure 9a shows the mean difference in the NDVI of all grids of different land-use type transfers. Most of the transfers had a positive impact on the NDVI, except for the transfers of farmland and water bodies to construction land, which decreased the NDVI slightly (reduction rates of −0.0173 and −0.0016, respectively). Moreover, the largest increment of the NDVI was due to farmland transferred into woodland (increment of 0.1411), woodland transferred into farmland (increment of 0.1359), and unused land transferred into farmland (increment of 0.1338). Figure 9b shows the total difference in the NDVI of all grids of different land-use type transfers. It can be seen that the largest increments of the NDVI were due to grassland, farmland, and woodland, and these land-use types did not change. The increments were 28,539.55, 25,545.42, and 12,450.86, respectively. The probable reasons are as follows. On the LP, the precipitation is limited and there exists a gap between agricultural water demand and supply. In the past, there used to be large areas of farmland that could not be irrigated adequately and were greatly affected by the precipitation. With the implementation of the GGP, some infertile areas of farmland that could not be irrigated adequately were returned to the woodland or grassland, and what remained was adequately irrigated or fertile areas of farmland that were less affected by the precipitation. Moreover, with the development of agricultural technologies and the optimization and adjustment of planting structures, the crop yields on the LP were continuously increased. (2) Spatial changes in the NDVI on the grid scale: The spatial changes in the NDVI on the grid scale on the LP during 2000-2015 are shown in Figure 10. Figure 10a indicates that the annual mean value of the NDVI on the LP was 0.529, with an uneven spatial distribution, i.e., decreasing from the southeastern to the northwestern areas. The high and low values were mainly distributed in the semi-humid and semi-arid areas, with annual mean values of 0.619 and 0.346, respectively. Figure 10b shows that the mean value of the coefficient of variation of the NDVI was 0.1406 on the LP, with an uneven spatial distribution, i.e., increasing from the southeastern semi-humid area (mean value of 0.1165) to the northwestern semi-arid area (mean value of 0.1926). There was 7.0% of the total area with 0 < CV < 0.05, mainly distributed in the eastern and southern regions. (2) Spatial changes in the NDVI on the grid scale: The spatial changes in the NDVI on the grid scale on the LP during 2000-2015 are shown in Figure 10. Figure 10a indicates that the annual mean value of the NDVI on the LP was 0.529, with an uneven spatial distribution, i.e., decreasing from the southeastern to the northwestern areas. The high and low values were mainly distributed in the semihumid and semi-arid areas, with annual mean values of 0.619 and 0.346, respectively. Figure 10b shows that the mean value of the coefficient of variation of the NDVI was 0.1406 on the LP, with an uneven spatial distribution, i.e., increasing from the southeastern semi-humid area (mean value of 0.1165) to the northwestern semi-arid area (mean value of 0.1926). There was 7.0% of the total area with 0 < CV < 0.05, mainly distributed in the eastern and southern regions. The area with 0.05 ≤ CV < 0.1 accounted for 25.9% of the  Figure 10c,d show that of the total area, the positive area accounted for 91.5% while the negative area accounted for 8.5%, which was obtained by subtracting the NDVI in 2000 from that in 2015. The NDVI in most areas of the LP (92.8% of the total area) increased from 2000 to 2015, indicating that the ecological environment in the region significantly improved. In addition, 70.0% of the area passed the significance test (p < 0.05); and 68.6% of the area increased significantly, while 1.4% of the area decreased significantly. The decrease in the NDVI was mainly concentrated in the semi-arid area as well as the urban area with rapid economic development and a large population concentration. the negative area accounted for 8.5%, which was obtained by subtracting the NDVI in 2000 from that in 2015. The NDVI in most areas of the LP (92.8% of the total area) increased from 2000 to 2015, indicating that the ecological environment in the region significantly improved. In addition, 70.0% of the area passed the significance test (p < 0.05); and 68.6% of the area increased significantly, while 1.4% of the area decreased significantly. The decrease in the NDVI was mainly concentrated in the semi-arid area as well as the urban area with rapid economic development and a large population concentration. (3) Temporal changes in the NDVI on the regional scale: The temporal changes in the NDVI on the regional scale on the LP during 2000-2015 are shown in Figure 11. The NDVI on the LP exhibited a fluctuating upward trend. The annual NDVI growth rates of the entire LP, the semi-arid area, and the semi-humid area were 0.0079, 0.0049, and 0.0093, respectively, indicating that the growth rate of the NDVI in the semi-humid area was higher than in the semi-arid area after the GGP was implemented. (3) Temporal changes in the NDVI on the regional scale: The temporal changes in the NDVI on the regional scale on the LP during 2000-2015 are shown in Figure 11. The NDVI on the LP exhibited a fluctuating upward trend. The annual NDVI growth rates of the entire LP, the semi-arid area, and the semi-humid area were 0.0079, 0.0049, and 0.0093, respectively, indicating that the growth rate of the NDVI in the semi-humid area was higher than in the semi-arid area after the GGP was implemented.

Individual Effects of Factors on NDVI
(1) Results of typical years: Using the typical annual values of 2000, 2005, 2010, and 2015 to drive the factor detector, the individual effects of the factors (represented by the q value) could be obtained, as shown in Figure 12.
(2) Results of annual mean differences: Using the differences between the annual mean values from 2000 to 2005 and those from 2010 to 2015 to drive the factor detector, the individual effects of the factors (represented by the q value) could be obtained, as shown in Figure 13.

Individual Effects of Factors on NDVI
(1) Results of typical years: Using the typical annual values of 2000, 2005, 2010, and 2015 to drive the factor detector, the individual effects of the factors (represented by the q value) could be obtained, as shown in Figure 12.
(a) (2) Results of annual mean differences: Using the differences between the annual mean values from 2000 to 2005 and those from 2010 to 2015 to drive the factor detector, the individual effects of the factors (represented by the q value) could be obtained, as shown in Figure 13.
In summary, the explanatory powers of the factors on the spatial distribution state of the NDVI in typical years and the spatial distribution change of the NDVI during 2000-2015 were different. In addition, the explanatory powers were different for the entire LP, and for the semi-arid/humid areas. With respect to the spatial distribution state of the NDVI in typical years, the decisive climatic factor was precipitation, the decisive environmental factors were geomorphic type, soil type, and vegetation type, and the decisive anthropogenic factor was land-use type. With respect to the spatial distribution change of the NDVI during 2000-2015, in the semi-arid area, the climatic and environmental factors were the decisive factors, including precipitation, soil type, and vegetation type; the impacts of anthropogenic factors, such as the GDP density, land-use type, and population density, were more significant in the semi-humid area.

Interactive Effects between Factors on NDVI
(1) Results of typical years: Using the typical annual values of 2000, 2005, 2010, and 2015 to drive the interaction detector, the interactive effects of the factors (represented by the q value) could be obtained, as shown in Figure 14. For the entire LP and the semi-arid/humid areas, the interactive effects between factors were greater than their individual effects, indicating that none of the factors acted independently, but they had a certain enhancement effect, including nonlinear enhancement and bi-enhancement.
For the entire LP, 26.4% of the interactive factor combinations exhibited nonlinear enhancement and 73.6% exhibited bi-enhancement. The interactive effect between precipitation and vegetation type was the strongest (mean annual?q = 0.7034), followed by precipitation and soil type (mean annual?q = 0.6987).
For the semi-arid area, 32.3% of the interactive factor combinations exhibited nonlinear enhancement and 67.7% exhibited bi-enhancement. The interactive effect between precipitation and soil type was the strongest (mean annual?q = 0.6375), followed by precipitation and vegetation type (mean annual?q = 0.6235).
For the semi-humid area, 44.1% of the interactive factor combinations exhibited nonlinear enhancement and 55.9% exhibited bi-enhancement. The interactive effect between precipitation and geomorphic type was the strongest (mean annual?q = 0.3961), followed by precipitation and vegetation type (mean annual?q = 0.3701).
(2) Results of annual mean differences: Using the differences between the annual mean values from 2000 to 2005 and those from 2010 to 2015 to drive the interaction detector, the interactive effects of the factors (represented by the q value) could be obtained, as shown in Figure 15. (2) Results of annual mean differences: Using the differences between the annual mean values from 2000 to 2005 and those from 2010 to 2015 to drive the interaction detector, the interactive effects of the factors (represented by the q value) could be obtained, as shown in Figure 15.
For the entire LP, 90.9% of the interactive factor combinations exhibited nonlinear enhancement and 9.1% exhibited bi-enhancement. The interactive effect between soil type and temperature was the strongest (q = 0.2194), followed by soil type and vegetation type (q = 0.2164).
For the semi-arid area, 83.6% of the interactive factor combinations exhibited nonlinear enhancement and 16.4% exhibited bi-enhancement. The interactive effect between precipitation and soil type was the strongest (q = 0.2438), followed by soil type and land-use type (q = 0.2379).
For the semi-humid area, 61.8% of the interactive factor combinations exhibited nonlinear enhancement and 38.2% exhibited bi-enhancement. The interactive effect between GDP density and geomorphic type was the strongest (q = 0.2470), followed by GDP density and soil type (q = 0.2300). For the entire LP, 90.9% of the interactive factor combinations exhibited nonlinear enhancement and 9.1% exhibited bi-enhancement. The interactive effect between soil type and temperature was the strongest (q = 0.2194), followed by soil type and vegetation type (q = 0.2164).
For the semi-arid area, 83.6% of the interactive factor combinations exhibited nonlinear enhancement and 16.4% exhibited bi-enhancement. The interactive effect between precipitation and soil type was the strongest (q = 0.2438), followed by soil type and land-use type (q = 0.2379).
For the semi-humid area, 61.8% of the interactive factor combinations exhibited nonlinear enhancement and 38.2% exhibited bi-enhancement. The interactive effect between GDP density and geomorphic type was the strongest (q = 0.2470), followed by GDP density and soil type (q = 0.2300). Sens. 2021, 13

Ranges or Types of Factors for NDVI
Using the typical annual values for 2000, 2005, 2010, and 2015 to drive the risk detector, the ranges or the types of factors for the NDVI could be obtained, as shown in Figure  16. The ranges or the types of factors for the NDVI were different in different years in the entire LP, and in the semi-arid/humid areas.
For the entire LP, the suitable precipitation range for the NDVI was greater than 650 mm (as the precipitation increased, the NDVI increased), the temperature range was smaller than 0 °C (as the temperature increased, the NDVI decreased first and then increased), the altitude ranged from 90 to 790 m or from 2963 to 4914 m (as the altitude increased, the NDVI decreased first and then increased), the slope ranged from 15° to 20°, the slope aspect was no slope aspect, the geomorphic type was large undulating mountain, the soil type was eluvial soil, the vegetation type was broad-leaved forest, the GDP density ranged from 500 to 1000 CNY/km 2 , the population density ranged from 200 to 500 people/km 2 , and the land-use type was woodland.
For the semi-arid area, the suitable precipitation for the NDVI ranged from 450 to 550 mm, the temperature range was smaller than 0 °C, the altitude ranged from 2963 to 4914 m, the slope ranged from 20° to 25°, the slope aspect was no slope aspect, the geomorphic type was large undulating mountain, the soil type was anthropogenic soil, the vegetation type was cultivated vegetation, the GDP density ranged from 500 to 1000 CNY/km 2 , the population density ranged from 200 to 500 people/km 2 , and the land-use type was farmland.
For the semi-humid area, the suitable precipitation range for the NDVI was greater than 650 mm, the temperature range was smaller than 0 °C, the altitude ranged from 2963 to 4914 m, the slope ranged from 15° to 20°, the slope aspect was no slope aspect, the geomorphic type was large undulating mountain, the soil type was eluvial soil, the vegetation type was broad-leaved forest, the GDP density ranged from 3000 to 5000 CNY/km 2 ,

Ranges or Types of Factors for NDVI
Using the typical annual values for 2000, 2005, 2010, and 2015 to drive the risk detector, the ranges or the types of factors for the NDVI could be obtained, as shown in Figure 16. The ranges or the types of factors for the NDVI were different in different years in the entire LP, and in the semi-arid/humid areas.
For the entire LP, the suitable precipitation range for the NDVI was greater than 650 mm (as the precipitation increased, the NDVI increased), the temperature range was smaller than 0 • C (as the temperature increased, the NDVI decreased first and then increased), the altitude ranged from 90 to 790 m or from 2963 to 4914 m (as the altitude increased, the NDVI decreased first and then increased), the slope ranged from 15 • to 20 • , the slope aspect was no slope aspect, the geomorphic type was large undulating mountain, the soil type was eluvial soil, the vegetation type was broad-leaved forest, the GDP density ranged from 500 to 1000 CNY/km 2 , the population density ranged from 200 to 500 people/km 2 , and the land-use type was woodland.
For the semi-arid area, the suitable precipitation for the NDVI ranged from 450 to 550 mm, the temperature range was smaller than 0 • C, the altitude ranged from 2963 to 4914 m, the slope ranged from 20 • to 25 • , the slope aspect was no slope aspect, the geomorphic type was large undulating mountain, the soil type was anthropogenic soil, the vegetation type was cultivated vegetation, the GDP density ranged from 500 to 1000 CNY/km 2 , the population density ranged from 200 to 500 people/km 2 , and the land-use type was farmland.
For the semi-humid area, the suitable precipitation range for the NDVI was greater than 650 mm, the temperature range was smaller than 0 • C, the altitude ranged from 2963 to 4914 m, the slope ranged from 15 • to 20 • , the slope aspect was no slope aspect, the geomorphic type was large undulating mountain, the soil type was eluvial soil, the vegetation type was broad-leaved forest, the GDP density ranged from 3000 to 5000 CNY/km 2 , the population density ranged from 200 to 500 people/km 2 , and the land-use type was woodland.

Differences of Significance between Factors on NDVI
Using the typical annual values for 2000, 2005, 2010, and 2015 to drive the ecological detector, the differences of significance between factors for the NDVI could be obtained, as shown in Table 5. For the entire LP, the differences were not significant between the precipitation and other factors in the different years (i.e., relatively certain). The significant differences between the soil type and the other factors did not vary in the different years (i.e., relatively certain). The soil type showed no significant differences with precipitation, GDP density, population density, and land-use type, but it exhibited significant differences with other factors. The differences of significance varied between the other pairs of factors in the different years (i.e., relatively uncertain). Table 5. Statistical significance of factors derived from the ecological detector with inputs of typical annual values.
indicates that there is a significant difference in the effect of two factors on NDVI (confidence is 95%); N indicates no significant difference; N/Y indicates that there is (or is no) significant difference across different years.

Differences of Significance between Factors on NDVI
Using the typical annual values for 2000, 2005, 2010, and 2015 to drive the ecological detector, the differences of significance between factors for the NDVI could be obtained, as shown in Table 5. For the entire LP, the differences were not significant between the precipitation and other factors in the different years (i.e., relatively certain). The significant differences between the soil type and the other factors did not vary in the different years (i.e., relatively certain). The soil type showed no significant differences with precipitation, GDP density, population density, and land-use type, but it exhibited significant differences with other factors. The differences of significance varied between the other pairs of factors in the different years (i.e., relatively uncertain). Table 5. Statistical significance of factors derived from the ecological detector with inputs of typical annual values.
Note: Y indicates that there is a significant difference in the effect of two factors on NDVI (confidence is 95%); N indicates no significant difference; N/Y indicates that there is (or is no) significant difference across different years.

Comparison between Semi-Arid and Semi-Humid Areas
Climatic factors influence the environmental conditions of the LP; while environmental factors determine the range and intensity of the human activities. The three types of factors are not independent of each other, and they interact with each other to influence the spatial distribution of the NDVI and its evolution over time on the LP. Based on the above results, there were significant spatial differences in the NDVI and its driving forces on both sides of the 400 mm isohyet on the LP.
For the semi-arid area to the northwest of the 400 mm isohyet, a number of observations were made. (1) The analysis of the spatial-temporal evolution revealed that the precipitation is relatively small and the temperature is relatively low, with annual mean values of 294 mm and 7.6 • C, respectively. The geomorphic types are mainly hills and plains (accounting for 12.3 and 11.0%, respectively). The soil types are mainly primary soil and arid soil (accounting for 11.1 and 8.0%, respectively). The main vegetation types are grassland, cultivated vegetation, and desert (accounting for 15.9, 7.2 and 5.0%, respectively). The grassland and farmland are the main land-use types (accounting for 16.4% and 6.9% in 2015, respectively). Moreover, the GDP density is relatively small (772 CNY/km 2 in 2015). The population density is relatively small (96 people/km 2 in 2015). The vegetation coverage is relatively low (annual mean NDVI of 0.346), and the coefficient of variation is relatively large (mean value of 0.1926). The growth rate of the NDVI was 0.0049/a after 2000 when the GGP was implemented. (2) The analysis based on the GDM revealed the following: The individual factors that determine the changes in the spatial distribution of the NDVI are the climatic and environmental factors, such as precipitation (q = 0.1476), soil type (q = 0.1373), vegetation type (q = 0.1353), temperature (q = 0.0500), and geomorphic type (q = 0.0417). The main interactive factors are precipitation and soil type (q = 0.2438), and soil type and land-use type (q = 0.2379). The suitable precipitation for the NDVI ranged from 450 to 550 mm, the temperature range was smaller than 0 • C, the altitude ranged from 2963 to 4914 m, the slope ranged from 20 • to 25 • , the slope aspect was no slope aspect, the geomorphic type was large undulating mountain, the soil type was anthropogenic soil, the vegetation type was cultivated vegetation, the GDP density ranged from 500 to 1000 CNY/km 2 , the population density ranged from 200 to 500 people/km 2 , and the land-use type was farmland.
Observations were also made regarding the semi-humid area to the southeast of the 400 mm isohyet. (1) The analysis of the spatial-temporal evolution revealed that both the precipitation and the temperature are relatively high, with annual mean values of 531 mm and 8.3 • C, respectively. The geomorphic types are mainly hills and medium undulating mountains (accounting for 19.4 and 13.3%, respectively). The soil types are mainly primary soil and semi-eluvial soil (accounting for 30.9 and 14.9%, respectively). The main vegetation types are cultivated vegetation and grassland (accounting for 37.8 and 9.2%, respectively). The farmland and grassland are the main land-use types (accounting for 25.4 and 24.8% in 2015, respectively). Moreover, the GDP density is relatively large (834 CNY/km 2 in 2015). The population density is relatively large (212 people/km 2 in 2015). The vegetation coverage is relatively high (annual mean NDVI of 0.619), and the coefficient of variation is relatively small (mean value of 0.1165). The growth rate of the NDVI was 0.0093/a after the implementation of the GGP in 2000. (2) The analysis based on the GDM revealed the following. The individual factors that determine the changes in the spatial distribution of the NDVI are the anthropogenic factors, such as GDP density (q = 0.1583), land-use type (q = 0.0761), and population density (q = 0.0520). The main interactive factors are GDP density and geomorphic type (q = 0.2470), and GDP density and soil type (q = 0.2300). The suitable precipitation range for the NDVI was greater than 650 mm, the temperature range was smaller than 0 • C, the altitude ranged from 2963 to 4914 m, the slope ranged from 15 • to 20 • , the slope aspect was no slope aspect, the geomorphic type was large undulating mountain, the soil type was eluvial soil, the vegetation type was broad-leaved forest, the GDP density ranged from 3000 to 5000 CNY/km 2 , the population density ranged from 200 to 500 people/km 2 , and the land-use type was woodland.
The conclusions above have important implications for ecological conservation and restoration in the different regions. The impacts of the climatic, environmental, and anthropogenic factors need to be comprehensively considered when formulating and implementing policies. The NDVI is sensitive to the climatic and environmental factors related to water changes in the semi-arid area where there is relatively little water. Therefore, appropriate vegetation types should be selected, regional vegetation structures should be optimized, and high water consumption vegetation should be gradually replaced by low water consumption vegetation such as grass and drought-tolerant plants. In the semihumid area, the regional water conditions basically meet the needs of vegetation growth because of the suitable climatic and natural conditions. Moreover, human activities have a relatively large disturbance effect on the NDVI. Thus, more positive human activities should be encouraged, including Grain for Green, afforestation, urban greening, and agricultural modernization projects. Furthermore, negative human activities should be controlled, including disordered urban expansion, population expansion, and overgrazing.

Connections and Distinctions with Other Studies
(1) Vegetation NDVI: In this study, it was found that the NDVI on the entire LP has increased significantly after 2000 when the GGP was implemented, consistent with the findings of [26,37,38].
(2) Climatic factors: In this study, it was found that among all factors, precipitation has the greatest influence on vegetation growth in the semi-arid area of the LP, consistent with the findings of [39][40][41]. In the semi-humid area, precipitation is not the decisive factor controlling vegetation growth because of the relatively good water conditions. Several studies have suggested that temperature is also a key factor affecting the vegetation growth. In general, when the temperature is low, the physiological activity of the vegetation is low. An increase in temperature promotes photosynthesis and vegetation growth, but as the temperature increases further, the growth of vegetation is inhibited by the accelerated evaporation and soil drying [42,43]. It was found in this study that conditions are suitable for vegetation growth when the temperature smaller than 0 • C on the LP. This is the case because the vegetation growth is influenced by many factors other than temperature. It can be seen from Figures 2-5 and 8 that the NDVI value is large in the southwest corner of semi-arid area of the LP, where the precipitation is large, the land-use types are mainly woodland and grassland, and the altitude ranges from 2963 m to 4914 m (the high altitude causes the low temperature). The finding is consistent with the results of [29]. It was also found that compared with other factors, temperature has relatively small explanatory power in terms of the spatial distribution and its change in the NDVI on the LP, which agrees with the findings and results of [39].
(3) Environmental factors: In this study, it was found that the environmental factors such as topographic types have an impact on the vegetation growth and restoration, consistent with the findings and results of [44]. The study revealed that the explanatory power of interactions between the factors would be increased for the spatial distribution and its change in the NDVI, which is consistent with the practical situation. For instance, the soil moisture is lower for steeper slopes, inhibiting vegetation growth. In this study, it was found that the NDVI value is large when the slope is between 15 • and 20 • , while it decreases when the slope becomes greater.
(4) Human factors: Certain studies have suggested that the regional land-use types have changed since the implementation of the GGP [45,46], which is not consistent with the results of this study. The reason for this may be that the spatial resolutions of the remote sensing data used in these studies are different. In this study, it was found that the conservation and restoration of the original woodland and grassland have been strengthened since the implementation of the GGP on the LP, resulting in a significant increase in the NDVI throughout the entire region. Additionally, many studies have used the NDVI to predict the grain yield because of the significant correlation between the grain yield and the NDVI [47,48]. The LP is a dry farming area, which has been affected by the scarcity of precipitation for a long time, resulting in a low grain yield [49]. With the agricultural technology developing recently, the optimization and adjustment of the crop planting structure has resulted in a continuous improvement in the grain yield on the LP. Therefore, the NDVI has increased significantly [19].

Possible Future Work
There is further research that can be carried out, including the following. (1) The factors affecting the spatial-temporal evolution of the NDVI are extremely complex. Due to data limitations, certain factors are not currently considered, such as agricultural fertilization, irrigation area, and CO 2 concentration [50,51], but can be further improved with the acquisition of new data in the future. (2) The interaction detector of the GDM only considers the interactions between two factors. There are interactive effects involving more than two factors. The combined impacts of multiple factors would provide insights on the evolution of vegetation, which can be considered in the future by introducing other methods, such as the analytical hierarchy process (AHP) and principal component analysis (PCA) [52]. (3) It is necessary to perform further collection of data with a relatively higher resolution through field investigations and monitoring to determine the practical situation of the surface undulations in the study area. (4) The cumulative and time-lag effects of various factors and indicators (e.g., standardized precipitation and evapotranspiration index, SPEI) should also be evaluated in the future [53]. (5) It is possible for the GDM to identify the impacts of driving forces on the evolution of vegetation from the perspective of spatial analysis and mathematical statistics. In order to further explore the evolutionary mechanism of vegetation, some methods from the field of landscape studies can be introduced in the future [54].

Conclusions
The main conclusions of this paper include the following.
(1) The spatial-temporal evolution characteristics of the NDVI on the entire LP, and in the semi-arid/humid areas during 2000-2015 were analyzed via the linear regression, coefficient of variation, and transfer matrix models.

•
The proportion of the total area with an NDVI of greater than 0.6 increased significantly, and the proportions of relatively low vegetation (NDVI of 0.2 to 0.4) converted to medium vegetation (NDVI of 0.4 to 0.6), medium vegetation converted to relatively high vegetation (NDVI of 0.6 to 0.8), and relatively high vegetation converted to high vegetation (NDVI of 0.8 to 1.0) accounted for the largest proportions.

•
The annual mean value of the NDVI on the LP was 0.529, decreasing from the southeastern semi-humid area (0.619) to the northwestern semi-arid area (0.346). The mean value of the NDVI coefficient of variation was 0.1406 on the LP, increasing from the southeastern semi-humid area (0.1165) to the northwestern semi-arid area (0.1926).

•
The NDVI on the LP exhibited an upward trend. The annual growth rate of the NDVI in the entire LP was 0.0079, and the growth rate in the semi-humid area (0.0093) was higher than in the semi-arid (0.0049) area after the GGP was implemented.

•
The area of the change in land-use types on the LP was not significant. Overall, a positive impact on the NDVI was found by the changes in the land-use type. The largest increments of the NDVI were due to grassland, farmland, and woodland, and these land-use types did not change.
(2) The GDM was adopted to quantitatively identify the impact of multiple factors on the spatial-temporal evolution of the NDVI on the entire LP, and in the semi-arid/humid areas.

•
Using the factor detector, it was found that in the semi-arid area, the climatic and environmental factors were the decisive factors influencing the spatial distribution changes of the NDVI during 2000-2015, including precipitation, soil type, and vegeta-tion type. The impacts of anthropogenic factors, such as the GDP density, land-use type, and population density, were more significant in the semi-humid area. • Using the interaction detector, it was found that the explanatory power of interactions between factors were greater than their individual effects, exhibiting two types of nonlinear enhancement and bi-enhancement. For the semi-arid area, 83.6% of the interactive factor combinations exhibited nonlinear enhancement and 16.4% exhibited bi-enhancement. The interactive effect between precipitation and soil type was the strongest. For the semi-humid area, 61.8% of the interactive factor combinations exhibited nonlinear enhancement and 38.2% exhibited bi-enhancement. The interactive effect between GDP density and geomorphic type was the strongest. • Using the risk and ecological detectors, the ranges or types of various factors that are suitable for vegetation growth and the differences of significance between factors for the NDVI on the LP were determined.
The conclusions of this study have important implications for policy makers and administrative managers in terms of the formulation and implementation of ecological conservation and restoration strategies in the different regions. In the semi-arid area, appropriate vegetation types should be selected, regional vegetation structures should be optimized, and high water consumption vegetation should be gradually replaced by low water consumption vegetation such as grass and drought-tolerant plants. In the semihumid area, more positive human activities should be encouraged, including Grain for Green, afforestation, urban greening, and agricultural modernization projects. Furthermore, negative human activities should be controlled, including disordered urban expansion, population expansion, and overgrazing.