Monitoring Vegetation Change and Its Potential Drivers in Inner Mongolia from 2000 to 2019

: Inner Mongolia in China is a typically arid and semi-arid region with vegetation promi-nently affected by global warming and human activities. Therefore, investigating the past and future vegetation change and its impact mechanism is important for assessing the stability of the ecosystem and the ecological policy formulation. Vegetation changes, sustainability characteristics, and the mechanism of natural and anthropogenic effects in Inner Mongolia during 2000–2019 were examined using moderate resolution imaging spectroradiometer normalized difference vegetation index (NDVI) data. Theil–Sen trend analysis, Mann–Kendall method, and the coefﬁcient of variation method were used to analyze the spatiotemporal variability characteristics and sustained stability of the NDVI. Furthermore, a trend estimation method based on a Seasonal Trend Model (STM), and the Hurst index was used to analyze breakpoints and change trends, and predict the likely future direction of vegetation, respectively. Additionally, the mechanisms of the compound inﬂuence of natural and anthropogenic activities on the vegetation dynamics in Inner Mongolia were explored using a Geodetector Model. The results show that the NDVI of Inner Mongolia shows an upward trend with a rate of 0.0028/year ( p < 0.05) from 2000 to 2019. Spatially, the NDVI values showed a decreasing trend from the northeast to the southwest, and the interannual variation ﬂuctuated widely, with coefﬁcients of variation greater than 0.15, for which the high-value areas were in the territory of the Alxa League. The areas with increased, decreased, and stable vegetation patterns were approximately equal in size, in which the improved areas were mainly distributed in the northeastern part of Inner Mongolia, the stable and unchanged areas were mostly in the desert, and the degraded areas were mainly in the central-eastern part of Inner Mongolia, it shows a trend of progressive degradation from east to west. Breakpoints in the vegetation dynamics occurred mainly in the northwestern part of Inner Mongolia and the northeastern part of Hulunbuir, most of which occurred during 2011–2014. The future NDVI trend in Inner Mongolia shows an increasing trend in most areas, with only approximately 10% of the areas showing a decreasing trend. Considering the drivers of the NDVI, we observed annual precipitation, soil type, mean annual temperature, and land use type to be the main driving factors in Inner Mongolia. Annual precipitation was the ﬁrst dominant factor, and when these four dominant factors interacted to inﬂuence vegetation change, they all showed interactive enhancement relationships. The results of this study will assist in understanding the inﬂuence of natural elements and human activities on vegetation changes and their driving mechanisms, while providing a scientiﬁc basis for the rational and effective protection of the ecological environment in Inner Mongolia. breakpoints in the at the using the and to conduct a (4) to predict in in using the and and to quantify the mechanisms of and of and on the using The of this can provide a basis for understanding and and for the of in These of for managing the of factors on in


Introduction
Vegetation is an important component of ecosystems and plays a key role in the carbon cycle [1][2][3]. Distinctive interannual and seasonal variations have been observed while studying vegetation trends [2,4,5]. Vegetation indices can provide substantial information on terrestrial vegetation; they are often used as an important parameter in characterizing surface vegetation for environmental quality assessments. Furthermore, vegetation indices are important for studying the hydrology, ecology, and regional changes [6,7]. Numerous studies have shown that the normalized difference vegetation index (NDVI) correlates well with the biomass and leaf area indices, allowing a good representation of the surface vegetation coverage. This can be used to effectively characterize vegetation activity and productivity and is suitable for representing changes in surface vegetation cover [1,[8][9][10]. With recent advancements in remote sensing technology, moderate resolution imaging spectroradiometer (MODIS) NDVI has become the primary data source for large-scale vegetation cover research because of its advantageous characteristics, such as high quality, wide areal coverage, high temporal resolution, and accessibility [11][12][13].
Currently, researchers mostly focus on the study of the overall trend of NDVI. Sen's slope is often used for quantifying the overall trend of the NDVI time series, as it is a non-parametric method that is less sensitive to outliers and skewed distributions [14]. However, this method can only be used to estimate the magnitude of NDVI changes and cannot determine the significance of the slope, so the results lack statistical significance. Therefore, researchers are increasingly applying the Mann-Kendall test to assess the significance of monotonic trends (linear or nonlinear) using a combination of Mann-Kendall test and Sen's slope algorithm for time series with outliers, to simultaneously assess and quantify monotonic interannual trends. This strategy is more robust than using parametric statistics, such as ordinary least squares. Therefore, the Mann-Kendall method is well suited to examine trends in vegetation [15,16]. However, only studying the interannual variations in NDVI can ignore the changes in the details of NDVI time series, and for a detailed description of NDVI changes, the ephemeral disturbances that occur in a time series need to be assessed [17]. Several methods have been proposed to detect sudden changes in vegetation, including the Landsat-based detection of trends in disturbance and recovery [18], breaks for additive season and trend (BFAST) [19], vegetation change tracker [20], and the detection of breakpoints and estimation of segments in trend methods. Due to the diversity and uncertainty of vegetation change, only high-frequency time series can describe the entire process of vegetation change in a short time interval [21,22]. Therefore, only change detection with dense satellite time series data can meet the requirements for the dynamic monitoring of vegetation cover changes [7]. A trend estimation method based on seasonal trend model (STM), originating from BFAST, considers the seasonality and noise of the NDVI time series. Therefore, it allows changes in the phenological cycle, as well as long-term NDVI trends, to be detected, distinguished, and quantified [23][24][25]. Therefore, it has received wide attention in the analysis of mutations that occur in a time series, and researchers are increasingly using this method for trend and breakpoint analysis of NDVI time series, time series smoothing and interpolation, and surface phenology analysis [25][26][27]. In this study, trends in the NDVI time series were explored using a combination of methods, analyzing both overall and mutational details.
Vegetation trends need to be fully characterized to provide information for decision making. Therefore, the Hurst index estimated by the R/S method was introduced to predict the future trend of vegetation in order to understand whether the vegetation trend has long-range memory. The R/S method provides specific information on correlation and persistence and is an effective index for studying complex processes such as vegetation time series [28,29]. Combining R/S and trend methods for future trend studies can provide insights into the continuity and future direction of vegetation change trends. Multiple researchers have introduced STM, Sen's slope, the Mann-Kendall test, and the Hurst index for studying changes in NDVI [30][31][32]. However, for the Inner Mongolia region, human activities such as national ecological projects, changes in grazing practices, and accelerated urbanization, together with the impact of climate change on vegetation cover, are all nonnegligible factors for the spatial heterogeneity of NDVI [33][34][35]. Therefore, it would be too simple to explore only the spatial and temporal variations in NDVI.
In recent years, the ecological problems in Inner Mongolia have become increasingly severe [36,37]. As an important ecological barrier in the northern region of China, the ecological condition in Inner Mongolia has an important impact in North China, and exploring the driving mechanism of vegetation change is an important prerequisite for solving ecological problems. Climate change and human activities are the main mechanisms studied in vegetation research [38,39]. The Fifth Assessment Report of the Intergovernmental Panel on Climate Change states that many aspects of vegetation phenology, composition, and productivity have been affected by global warming over the past 30 years [40][41][42]. Many studies have attempted to quantify the relationship between vegetation indices and meteorological factors, but most focused on linear relationships and only on two factors, temperature and precipitation, leading to inadequate conclusions and failing to quantify the contribution and interactions of the drivers of vegetation change [43,44]. Additionally, the impact of human activities on vegetation cannot be ignored. Attempts have been made to study this using residual analysis, but this method is crude, allowing the quantification of the overall impact of human activities, but being unable to identify the human activity factor that is specifically responsible. In addition to the methods mentioned above, multiple linear regression methods [43], structural equation models [45], random forest models [46], and other methods have been introduced to study the effects of human activities and meteorological elements on NDVI. However, the influence mechanisms of vegetation are complex and diverse, and it is not sufficient to study complex multivariate nonlinear relationships using simple linear correlations. In order to obtain more objective and accurate evaluation results, we detected NDVI spatial differentiation and related drivers from a multivariate spatial perspective, using a Geodetector approach that can avoid the interference of subjective factors [47,48]. We also explored the potential and interactive effects of anthropogenic and natural factors on vegetation growth.
The Inner Mongolia region is located in the northern frontier of China, with a vast territory, high terrain, and complex and diverse landscapes, and is one of the most sensitive regions in the world to external environmental changes due to its arid and semi-arid climate [49]. Over the past 40 years, temperatures in the region have risen much faster than the global average, and the rapid rise in temperature and dramatic reduction in precipitation have exacerbated the degradation of vegetation [50,51]. In addition, since 2000, the region has entered a phase of rapid socioeconomic development, and climate change and human activities are rapidly restructuring the structure and function of vegetation systems at all levels and simultaneously affecting the ecology and future sustainability of the region [52]. However, current studies on the driving mechanisms of vegetation growth focus on natural factors, and anthropogenic factors are often neglected. Therefore, in this study, we use the Geodetector model combined with geospatial data and climate and human activity factor data, for a comprehensive investigation of the mechanisms that influence NDVI spatial patterns. DEM, slope and land type, and precipitation and temperature data were selected for analyzing the driving factors, after considering the special topographic and climatic conditions of the Inner Mongolia region [53]. In terms of anthropogenic factors, the impact of human activities on the ecological environment is becoming more and more complex as urbanization accelerates, and we chose GDP, land use type, and population data to explore this aspect [54]. In addition, as Inner Mongolia is a major livestock province in China, the impact of grazing intensity on vegetation change should not be neglected. We chose the number of livestock as a representative driver to examine its influence [55].
Overall, the objectives of this study were: (1) To study the interannual variations of the NDVI characteristics during 2000-2019 for each vegetation type in Inner Mongolia using time series analyses and Mann-Kendall tests; (2) to assess the NDVI stability over the last 20 years in Inner Mongolia by calculating the coefficient of variation; (3) to detect

Data Sources
The data used in this research (Table 1) were pre-processed to meet the accuracy requirements of the research and for consistency between different data sources. MODIS reprojection tool (MRT) was used to pre-process the NDVI data. Analysis included reprojection, stitching, and format conversion, while the original HDF format was transformed into the ArcGIS editable TIFF file, and its coordinate system and projection were defined.

Data Sources
The data used in this research (Table 1) were pre-processed to meet the accuracy requirements of the research and for consistency between different data sources. MODIS reprojection tool (MRT) was used to pre-process the NDVI data. Analysis included reprojection, stitching, and format conversion, while the original HDF format was transformed into the ArcGIS editable TIFF file, and its coordinate system and projection were defined. Furthermore, MATLAB was used to remove outliers through the maximum value synthesis method to obtain the optimal NDVI raster. The meteorological data was interpolated to monthly scale spatial data using Albers projection and Kriging interpolation, before being synthesized into annual precipitation and annual mean temperature spatial data at the monthly scale using Raster Calculator. Finally, the data were masked using the administrative boundaries of the Inner Mongolia region. We also used livestock data from 2010 [59], gross domestic product (GDP) data, and human population (POP) data to represent different anthropogenic factors. Specific information is shown in Table 1 and Figure 2.

Spatial Trend Analysis Methods
The combination of the Theil-Sen trend analysis and the Mann-Kendall test does not require the data to follow a certain distribution when analyzing trends over a time series. Furthermore, the combination is more resistant to data errors and the results from the significance level test are more scientifically reliable than other research methods [29].
Theil-Sen trend analysis is a very reliable non-parametric statistical calculation method that is computationally efficient, insensitive to measurement errors and outlier data, and is often used in the trend analysis of long time series. This is calculated as follows: where NDVIj and NDVIi are the NDVI values of a certain pixel in years j and i, respectively. In this study, 2019 ≥ j ≥ i ≥ 2000. Additionally, the median represented here is the median of the requested series. When β is greater than 0, an increasing trend is observed in the NDVI and when β is less than 0, the NDVI presents a decreasing trend. When β is equal to 0, the NDVI is stable.
Remote Sens. 2021, 13, x FOR PEER REVIEW 7 of 28 Figure 2. NDVI drivers and spatial distribution in Inner Mongolia.

Spatial Trend Analysis Methods
The combination of the Theil-Sen trend analysis and the Mann-Kendall test does not require the data to follow a certain distribution when analyzing trends over a time series. Furthermore, the combination is more resistant to data errors and the results from the significance level test are more scientifically reliable than other research methods [29].
Theil-Sen trend analysis is a very reliable non-parametric statistical calculation method that is computationally efficient, insensitive to measurement errors and outlier data, and is often used in the trend analysis of long time series. This is calculated as follows: The Mann-Kendall method is a non-parametric statistical test that was originally proposed by Mann in 1945 and further refined by Kendall and Sneyers. It is advantageous as it does not require normally distributed measurements or a linear trend and remains unaffected by missing values or outliers. It is widely used for the trend significance testing of long-series data [61]. The Mann-Kendall statistic S is calculated for the time series {NDVIi} as follows: where n is the length of the time series, and sgn is the sign function, defined as: When n ≥ 8, S follows a normal distribution with a mean of zero and a variance of: The standardization of S is as follows: where Z c is the statistic normalized by the Mann-Kendall test and follows a normal distribution. When β = 0, there is no monotonic trend in the study series when the null hypothesis holds and is combined with Theil-Sen trend analysis. Conversely, when the null hypothesis does not hold, there is significant variation in the analyzed series at the alpha level if the absolute value of Z c is greater than the standard normal distribution. This study tested the significance of trends in NDVI time series at a confidence level of α = 0.05. Additionally, the NDVI values for each pixel in the study area during 2000-2019 were subjected to Theil-Sen trend analysis to obtain the spatial distribution of the trend values. If β > 0, the annual average NDVI of the pixel shows an increasing trend, whereas if β < 0, the annual average NDVI of the pixel shows a decreasing trend. Because there is essentially no region where β is strictly equal to 0, the interval between −0.001 and 0.001 was classified as stable and constant in this research, based on the actual situation of β. The results of the Mann-Kendall test, at a confidence level of 0.05, were classified as significant (|Z c | > 1.96) or insignificant changes (|Z c | ≤ 1.96).
In this study, areas with a multi-year average NDVI of less than 0.08 were considered to be unvegetated and were excluded [62].

Coefficient of Variation
The coefficient of variation (C v ) is a statistical index that describes the ratio of a variable relative to its mean, and is regarded as a useful indicator of interannual variability in an ecosystem. Specifically, a lower C v value indicates a lower fluctuation and greater stability of the vegetation, while a higher C v value indicates less stable vegetation. The calculation is as follows: where NDVI i is the NDVI value in year i and NDV I is the average NDVI value from February 2000 to December 2019.
In this study, the stability of the interannual variation of the NDVI in Inner Mongolia was expressed using C v to indicate the fluctuation of NDVI in each region over a 20 year period. The C v of the NDVI in Inner Mongolia during 2000-2019 was calculated pixel-by-pixel.

Trend Estimation Based on Seasonal Trend Model
The STM is based on the BFAST and the classical additive decomposition model, which constructs estimation methods for trends and breakpoints [63,64]. A segmented Remote Sens. 2021, 13, 3357 8 of 25 linear trend and seasonal regression model of NDVI during 2000-2019 to explain the NDVI value y at moment t can be expressed as: where α 1 is the intercept, α 2 is the slope of the trend, γ is the amplitude, δ is the phase of the k harmonic term, ε is the residual error, and f is the number of observations per year. The parameters α 1 and α 2 were estimated using an ordinary least squares regression, wherein the derived time series segments were considered as categorical interaction terms with a trend slope of α 2 . The significance of the trend in each segment was estimated using a t-test for the interaction parameter of the regression between the time series segment and α 2 .
The presence of breakpoints can, under certain circumstances, indicate a shift in the mechanism of influence on the time series. Therefore, if breakpoints are not considered, the results of the trend analysis may lead to a misjudgment of the factors that influence vegetation changes in Inner Mongolia [65,66]. Therefore, in this study, the times at which breakpoints occurred were used as time splitting points, at which the time series was split into two subseries; before and after the breakpoint. Then, the time trends and significance levels of the subseries were monitored separately to obtain accurate conclusions regarding the fluctuations in vegetation trends in Inner Mongolia during 2000-2019.

Hurst Exponent and R/S Analysis
The Hurst exponent and R/S analysis have developed a fractal theory for studying time series and have been widely used in the research of urban environmental change, climate change, population, and economic development [67,68]. R/S analysis can be used to predict future trends in vegetation cover, wherein the main principle is to construct a time series that defines an average series { ξ(t), t = 1, 2 · · · ·} , for any positive integer τ ≥ 1. First, the mean sequence must be defined: Cumulative deviation: The range: Standard deviation: there is a relationship as follows: where H is the Hurst index. According to the value of (τ, R/S), H can be obtained using least squares fitting (Equation (8)) in a double logarithmic coordinate system. The H value provides information on whether the sequence is completely random or if it follows a trend. Furthermore, trends can be expressed as persistent or anti-persistent. Overall, H= 0.5 indicates that vegetation change is random, 0.5 < H < 1 indicates that vegetation change is sustainable, and 0 < H < 0.5 indicates that future vegetation trend is expected to be the opposite of that in the past.

GeoDetector Model
GeoDetector is a statistical method used to explore spatial differentiation and quantify drivers [69]. In this study, four detector tools of the geographical detector model were used to quantify the drivers of NDVI in Inner Mongolia.
The dependent variable Y is the NDVI value for Inner Mongolia and the independent variable X includes natural and social drivers.
(1) The factor detector focuses on exploring the spatial heterogeneity of the dependent variable Y and the extent to which the independent variable X explains the spatial heterogeneity of Y. The expressions are: where q represents the explanatory power of the metric factors on Y. Specifically, a larger q indicates that the explanatory power of each metric factor is stronger, indicating that it has more influence on the spatial distribution of the NDVI, and the range of q values is {0,1}. Note that h = 1, 2 . . . , L is the classification or partition of Y or X; N h and N represent the number of cells in layer h and the entire area, respectively; δ 2 h and δ 2 are the variances of layer h and the entire area, respectively; and SSW and SST are the sum of the within-layer variance and total variance of the entire area, respectively.
(2) An ecological detector was used to determine if a significant difference existed between the effects of two drivers on the spatial distribution of the dependent variable, as measured by the F-statistic: where N x1 and N x2 denote the sample sizes of the two drivers, and SSW x1 and SSW x2 denote the sum of the within-layer variance for the stratum formed by the two drivers.
(3) An interaction detector was used to identify interactions between the different drivers by assessing whether the explanatory power of factors X1 and X2 increases or decreases for the dependent variable Y when they function together [70]. It is also possible that the effects of these factors on the dependent variables were independent. The results of this assessment were obtained by comparing the q-values q(X1) and q(X2) for the two factors on Y, with the q-values on Y after the interaction between the two factors.
(4) A risk detector was used to detect the appropriate extent or type of vegetation cover impact by different drivers, and uses the t-statistic as follows: where Y h is the mean value of the NDVI attributes in sub-region h, n h is the number of samples in sub-region h, and Var denotes variance. GeoDetector was used to analyze variables that must be discretized for continuous types of data. Additionally, the natural intermittency grading method was used to maximize the differences between the groups. This method considers the large latitudinal span of Inner Mongolia, the diverse climate, temperature, soil, etc., and examines the representativeness, typicality, and scientificity of various driving factors. Nine representative and available datasets were selected as drivers for exploring the distribution and changes in the NDVI in Inner Mongolia, including five natural factors: Precipitation, temperature, elevation, soil type, and slope; and four social factors: Population, GDP, land use type, and livestock density. Based on the closest distance method, the average NDVI and nine drivers during 2000-2019 were resampled using ArcGIS 10.6, and the spatial resolution of the raster data was unified to 10 × 10 km. The natural breakpoint method was used to discretize the continuous driver data, and the spatial distribution is shown in Figure 2.

NDVI Temporal Annual Variation Characteristics
To study the annual changes in the NDVI

Spatial Distribution and Trend Analysis of NDVI
The remaining regions were classified into five levels using the equal spacing method in ArcGIS 10.6, into low (0.08-0.25), medium-low (0.25-0.41), medium (0.41-0.58), medium-high (0.58-0.75), and high (0.75-0.92) vegetation cover (Figure 4a) [29]. From the spatial distribution of NDVI levels, the high vegetation cover (0.75-0.92) area was the largest, accounting for 24.95% of the total area. Influenced by the large east-west span of Inner Mongolia, diverse and complex topography, large vertical differences in vegetation growth, and the complex spatial distribution, the NDVI gradually decreased from east to west from a regional perspective. Hulunbuir, the Xing'an League, Tongliao, and Chifeng, locations with forests and typical steppes, recorded high NDVI values, while areas with NDVI values less than 0.41 were mainly in the western Alxa League, Ordos, and Bayannur areas, where the vegetation cover types largely consist of desert and desert steppe. Figure 4b shows the results of the Sen's slope and Mann-Kendall test for annual NDVI. Referring to the grading method of Lu Hao et al. [71], we identified three trend

Spatial Distribution and Trend Analysis of NDVI
The remaining regions were classified into five levels using the equal spacing method in ArcGIS 10.6, into low (0.08-0.25), medium-low (0.25-0.41), medium (0.41-0.58), mediumhigh (0.58-0.75), and high (0.75-0.92) vegetation cover (Figure 4a) [29]. From the spatial distribution of NDVI levels, the high vegetation cover (0.75-0.92) area was the largest, accounting for 24.95% of the total area. Influenced by the large east-west span of Inner Mongolia, diverse and complex topography, large vertical differences in vegetation growth, and the complex spatial distribution, the NDVI gradually decreased from east to west from a regional perspective. Hulunbuir, the Xing'an League, Tongliao, and Chifeng, locations with forests and typical steppes, recorded high NDVI values, while areas with NDVI values less than 0.41 were mainly in the western Alxa League, Ordos, and Bayannur areas, where the vegetation cover types largely consist of desert and desert steppe.

NDVI Stability Analysis
In terms of the area occupied by each type of variation level, low fluctuation areas > medium-low fluctuation areas > high fluctuation areas > medium fluctuation areas > medium-high fluctuation areas. As shown in Figure 5, areas with large fluctuations in NDVI accounted for 27.2% or the entire area, and were mainly distributed in Xilingol league, parts of Ulanqab and Baotou, and Alxa league. This is due to the sparse precipitation and poor soils in these areas. The corresponding vegetation types are mainly grasslands, deserts, or transition areas between grasslands and deserts. These have fragile ecological environments and vegetation cover that is vulnerable to external influences. Areas with low NDVI fluctuations accounted for 61.7% of the total area, mainly distributed in Hulunbuir, the Xing'an League, Tongliao, and Chifeng, with the vegetation types of forest, typical steppe, and meadow steppe. Fluctuation of NDVI in the Xilingol League, where typical steppes are distributed, was also low, which may be related to the good ecological environmental conditions of its forest and grassland areas; these are not easily influenced by external forces [1]. Additionally, the amplitude of fluctuations in grassland areas was greater than that in forest areas, which may be attributed to the grassland areas being under high pressure from population and livestock, and more strongly disturbed by human activities than forest areas.  Referring to the grading method of Lu Hao et al. [71], we identified three trend types, "increased," "stable and constant," and "decreased." Increased and decreased represent areas where the rising and falling NDVI trends pass the 0.05 significance level in the Mann-Kendall test, respectively, and we classified the remaining areas with NDVI trends that did not pass the 0.05 significance level into stable and constant. Area increasing NDVI were larger than those with decreasing NDVI, which is 34.65%, while 33.71% of the regions showed stable and constant NDVI. Increased and decreased NDVI mostly occur in eastern Hulunbeier, Xing'an Meng, Tongliao, southern Erdos, and southeastern Chifeng. The areas with stable and constant NDVI were mainly located in Alxa League. Finally, areas where NDVI decreased significantly were concentrated in Baotou, Hohhot, Ulanqab, Xilingol league, and the western part of Hulunbuir.

NDVI Stability Analysis
In terms of the area occupied by each type of variation level, low fluctuation areas > medium-low fluctuation areas > high fluctuation areas > medium fluctuation areas > medium-high fluctuation areas. As shown in Figure 5, areas with large fluctuations in NDVI accounted for 27.2% or the entire area, and were mainly distributed in Xilingol league, parts of Ulanqab and Baotou, and Alxa league. This is due to the sparse precipitation and poor soils in these areas. The corresponding vegetation types are mainly grasslands, deserts, or transition areas between grasslands and deserts. These have fragile ecological environments and vegetation cover that is vulnerable to external influences. Areas with low NDVI fluctuations accounted for 61.7% of the total area, mainly distributed in Hulunbuir, the Xing'an League, Tongliao, and Chifeng, with the vegetation types of forest, typical steppe, and meadow steppe. Fluctuation of NDVI in the Xilingol League, where typical steppes are distributed, was also low, which may be related to the good ecological environmental conditions of its forest and grassland areas; these are not easily influenced by external forces [1]. Additionally, the amplitude of fluctuations in grassland areas was greater than that in forest areas, which may be attributed to the grassland areas being under high pressure from population and livestock, and more strongly disturbed by human activities than forest areas.

Evaluation of Estimated Breakpoints and Trends
STM was applied to the NDVI sequences in Inner Mongolia, and the number of breakpoints and trend changes before and after the mutations and corresponding signifi-

Evaluation of Estimated Breakpoints and Trends
STM was applied to the NDVI sequences in Inner Mongolia, and the number of breakpoints and trend changes before and after the mutations and corresponding significance levels were calculated (Figures 6 and 7). Considering the variations, approximately 16

Future NDVI Trends
In this study, the Hurst index values for annual NDVI data were calculated on a pixel-by-pixel basis using MATLAB (Figure 8a). The Hurst index of NDVI in Inner Mongolia is low, with maximum, minimum, and mean values of 0.93, 0.09 and 0.44, respectively. The proportion of areas with a Hurst index less than 0.5 was 74.91%, while that of areas greater than 0.5 was 25.08%, indicating that the reversal of the vegetation trend in the future is likely in most areas of Inner Mongolia. This indicates that future changes will be the opposite to those in the present. Furthermore, the Hurst values are low in the border areas between Inner Mongolia and Mongolia and relatively high in other areas. Additionally, the areas with a Hurst index greater than 0.5 are concentrated in the south of Chifeng, Xilingol League, and Alxa League.
The Sen's slope and Mann-Kendall test can gradually quantify trends in the NDVI, and the Hurst index can be used to qualitatively predict whether future trends will be similar to those at present. However, neither of these methods can predict future increases or decreases in trends. The NDVI trends were coupled with the Hurst index classification results to obtain the spatial distribution of future changes of vegetation cover in Inner Mongolia (Figure 8b), and the future change trends were classified into the following four levels: "increased," "constant," "decreased," and "unable to determine," which were used to provide reasonable predictions of future vegetation change trends. Where "increased" and "decreased" represent predicted increases and decreases in future NDVI trends passing the 0.05 significance level in the Mann-Kendall test, respectively. The remaining NDVI trends that did not pass the 0.05 significance level, as well as the areas that were unchanged were notes as constant. The results showed that the area with stable and constant NDVI was the largest. Areas with increased NDVI in Inner Mongolia are larger than areas with decreased NDVI. Increased areas account for approximately 24.62% of the total area and are mainly located in Hulunbeier, Xing'an League, eastern Tongliao, and south-central Alxa League, which may be due to the large-scale ecological management and soil and water conservation and afforestation policies in these areas. Meanwhile, decreased areas, which account for about 10.56% of the total area of Inner Mongolia, are concentrated According to the trends and the significance of changes before and after the mutations, we divided the trends into three levels of "increased," "stable," and "decreased." The vegetation changes in the southeastern part of Alxa League, northwestern part of Bayannur, and most parts of Ordos before the mutation showed a significant upward trend. An obvious decreasing trend was observed in the south-central and eastern part of Alxa Left Banner, the south-central part of Ulanqab, Erguna, and Genhe of Hulunbeier, while the vegetation in the areas along the northwestern border of Baotou and Ulanqab showed a stable NDVI. Meanwhile, the areas with significantly increasing vegetation trends after the mutation occurred in the northern part of Alxa League, the northwestern part of Bayannur, and the border areas of Baotou and Ulanqab. The areas with significantly decreasing trends were in the central and eastern parts of Alxa League. Northwestern Hulunbeier, northern Hohhot, southern Ordos, and northwestern Bayannur, as well as small areas in Xing'an League, Tongliao, and Chifeng showed stable NDVI trends after mutation.
The spatial distribution of NDVI trends before and after the mutation shows obvious variability. For example, most pixels in the south-central region of Alxa League showed a significant increasing trend before an abrupt change to a decreasing trend after the mutation. Meanwhile, the NDVI values at the borders of Baotou and Ulanqab changed from a stable trend to an obvious increasing trend after the mutation, and the NDVI trends in most of the central-eastern parts of Inner Mongolia, which originally showed obvious decreases, gradually stabilized after the mutation. The occurrence of breakpoint years and the differences in the trends before, and after, the mutation may indicate changes in the potential influence mechanisms of NDVI in Inner Mongolia.

Future NDVI Trends
In this study, the Hurst index values for annual NDVI data were calculated on a pixelby-pixel basis using MATLAB (Figure 8a). The Hurst index of NDVI in Inner Mongolia is low, with maximum, minimum, and mean values of 0.93, 0.09 and 0.44, respectively. The proportion of areas with a Hurst index less than 0.5 was 74.91%, while that of areas greater than 0.5 was 25.08%, indicating that the reversal of the vegetation trend in the future is likely in most areas of Inner Mongolia. This indicates that future changes will be the opposite to those in the present. Furthermore, the Hurst values are low in the border areas between Inner Mongolia and Mongolia and relatively high in other areas. Additionally, the areas with a Hurst index greater than 0.5 are concentrated in the south of Chifeng, Xilingol League, and Alxa League.
in the southern part of Inner Mongolia, which is mainly arable land and highly influenced by human activities. In addition, the Wuliang Suhai area in Bayannur, located in the lower reaches of the Yellow River, has experienced vegetation degradation due to increasing ecological pollution [72]. The areas with stable and constant future trends of NDVI accounted for 64.81% of the total area and were distributed throughout the region. In the areas with stable and constant NDVI levels, areas with positive Sen's slope values and p > 0.05 are the largest. This indicates that areas with increasing NDVI in Inner Mongolia will increase more areas with decreasing NDVI in the future.

Driving Factor Impact
The results of the factor detector reflect the influence of the magnitude of each driver on the NDVI in Inner Mongolia, as expressed by the q value of each driver ( Figure 9). Specifically, higher q values indicate a greater influence on the NDVI. As shown in Figure 9, the q value of the total annual precipitation is the largest, reaching 0.81, followed by that of soil type, mean annual temperature, and land use type, all of which have a significant influence greater than 0.5. This indicates that total annual precipitation was the dominant factor of NDVI change in Inner Mongolia during 2000-2019, followed by soil type. Among the factors with an influence below 0.5, the influence of natural factors is higher than that of social factors, in the following order: DEM > slope > population > GDP > livestock density. DEM has a significant influence of 0.34, and the slope has an average influence of 0.13 on the NDVI. The other three drivers had q-values below 0.1, indicating that they did not have a substantial influence on the NDVI. The Sen's slope and Mann-Kendall test can gradually quantify trends in the NDVI, and the Hurst index can be used to qualitatively predict whether future trends will be similar to those at present. However, neither of these methods can predict future increases or decreases in trends. The NDVI trends were coupled with the Hurst index classification results to obtain the spatial distribution of future changes of vegetation cover in Inner Mongolia (Figure 8b), and the future change trends were classified into the following four levels: "increased," "constant," "decreased," and "unable to determine," which were used to provide reasonable predictions of future vegetation change trends. Where "increased" and "decreased" represent predicted increases and decreases in future NDVI trends passing the 0.05 significance level in the Mann-Kendall test, respectively. The remaining NDVI trends that did not pass the 0.05 significance level, as well as the areas that were unchanged were notes as constant. The results showed that the area with stable and constant NDVI was the largest. Areas with increased NDVI in Inner Mongolia are larger than areas with decreased NDVI. Increased areas account for approximately 24.62% of the total area and are mainly located in Hulunbeier, Xing'an League, eastern Tongliao, and south-central Alxa League, which may be due to the large-scale ecological management and soil and water conservation and afforestation policies in these areas. Meanwhile, decreased areas, which account for about 10.56% of the total area of Inner Mongolia, are concentrated in the southern part of Inner Mongolia, which is mainly arable land and highly influenced by human activities. In addition, the Wuliang Suhai area in Bayannur, located in the lower reaches of the Yellow River, has experienced vegetation degradation due to increasing ecological pollution [72]. The areas with stable and constant future trends of NDVI accounted for 64.81% of the total area and were distributed throughout the region. In the areas with stable and constant NDVI levels, areas with positive Sen's slope values and p > 0.05 are the largest. This indicates that areas with increasing NDVI in Inner Mongolia will increase more areas with decreasing NDVI in the future.

Driving Factor Impact
The results of the factor detector reflect the influence of the magnitude of each driver on the NDVI in Inner Mongolia, as expressed by the q value of each driver ( Figure 9). Specifically, higher q values indicate a greater influence on the NDVI. As shown in Figure 9, the q value of the total annual precipitation is the largest, reaching 0.81, followed by that of soil type, mean annual temperature, and land use type, all of which have a significant influence greater than 0.5. This indicates that total annual precipitation was the dominant factor of NDVI change in Inner Mongolia during 2000-2019, followed by soil type. Among the factors with an influence below 0.5, the influence of natural factors is higher than that of social factors, in the following order: DEM > slope > population > GDP > livestock density. DEM has a significant influence of 0.34, and the slope has an average influence of 0.13 on the NDVI. The other three drivers had q-values below 0.1, indicating that they did not have a substantial influence on the NDVI.

Appropriate Range of Different Driving Factors
The most suitable range or type of vegetation growth in Inner Mongolia for factor was analyzed using the risk detector method, and the results of each f passed the statistical test at 95% confidence ( Figure 10). Among the findings, annual precipitation, which was the most influential driver of the NDVI in Mongolia, showed a positive correlation with the NDVI, and the mean value o NDVI reached a maximum when the total annual precipitation was 452-535 mm areas with the highest precipitation are mainly located in Hulunbuir and the Xin League in northeastern Inner Mongolia, where forests and typical steppe are pre these areas had the highest vegetation cover in Inner Mongolia. This indicates arid and semi-arid areas are more sensitive to changes in precipitation condi and that increased precipitation has a more pronounced effect on the growth o ests and typical steppe landscapes [73,74]. As grassland forest areas with high tation cover in Inner Mongolia are widely distributed with luvisols, which have soil fertility and are suitable for forestry, the mean NDVI values were the high areas with luvisols. Moreover, an increase in temperature stimulates the grow vegetation, but also leads to a decrease in soil moisture, negatively impacting tation growth. Therefore, the annual mean temperature in Inner Mongolia is tively correlated with NDVI, and the overall annual mean temperature in Inner

Appropriate Range of Different Driving Factors
The most suitable range or type of vegetation growth in Inner Mongolia for each factor was analyzed using the risk detector method, and the results of each factor passed the statistical test at 95% confidence ( Figure 10). Among the findings, total annual precipitation, which was the most influential driver of the NDVI in Inner Mongolia, showed a positive correlation with the NDVI, and the mean value of the NDVI reached a maximum when the total annual precipitation was 452-535 mm. The areas with the highest precipitation are mainly located in Hulunbuir and the Xing'an League in northeastern Inner Mongolia, where forests and typical steppe are present; these areas had the highest vegetation cover in Inner Mongolia. This indicates that arid and semi-arid areas are more sensitive to changes in precipitation conditions, and that increased precipitation has a more pronounced effect on the growth of forests and typical steppe landscapes [73,74]. As grassland forest areas with high vegetation cover in Inner Mongolia are widely distributed with luvisols, which have high soil fertility and are suitable for forestry, the mean NDVI values were the highest in areas with luvisols. Moreover, an increase in temperature stimulates the growth of vegetation, but also leads to a decrease in soil moisture, negatively impacting vegetation growth. Therefore, the annual mean temperature in Inner Mongolia is negatively correlated with NDVI, and the overall annual mean temperature in Inner Mongolia ranges from −3.7 to 11.2 • C, while the most suitable mean temperature range for vegetation growth is within this range. The mean NDVI value is the highest in Inner Mongolia. The mean NDVI of vegetation fluctuates with DEM and slope, wherein the highest mean NDVI values occur at a DEM of 505-840 m and a slope of 15-21 • . Among social drivers, the implementation of ecological restoration projects in recent years played an important role in the recovery of vegetation in Inner Mongolia, while population, GDP, and livestock are unevenly distributed in a few areas, occupying a small proportion of the area, thereby having a low overall impact on vegetation growth.

Driving Factor Interactions
The spatial pattern of the NDVI changes in Inner Mongolia is not influenced by a single factor, but rather by a combination of drivers. Therefore, interaction detectors were used to detect the interrelationships between different drivers that affect NDVI changes in vegetation ( Figure 11). The results show that the interactions between different drivers are all bi-variable or non-linear enhancements, with most being bi-variable enhancements. There were no cases wherein the drivers were independent or weakened. As previously discussed, soil moisture and temperature have an important compounding effect on vegetation growth through a strong coupling of land and atmosphere [75]. Therefore, the coupling of precipitation and temperature plays an important role in vegetation growth because of the increased evapotranspiration caused by global warming [76]. The results of the interaction detectors in this study also reflect the influence of the interactions between the factors, whereof total annual precipitation, mean annual temperature, and soil type had greater influences, and

Driving Factor Interactions
The spatial pattern of the NDVI changes in Inner Mongolia is not influenced by a single factor, but rather by a combination of drivers. Therefore, interaction detectors were used to detect the interrelationships between different drivers that affect NDVI changes in vegetation (Figure 11). The results show that the interactions between different drivers are all bi-variable or non-linear enhancements, with most being bi-variable enhancements. There were no cases wherein the drivers were independent or weakened. As previously discussed, soil moisture and temperature have an important compounding effect on vegetation growth through a strong coupling of land and atmosphere [75]. Therefore, the coupling of precipitation and temperature plays an important role in vegetation growth because of the increased evapotranspiration caused by global warming [76]. The results of the interaction detectors in this study also reflect the influence of the interactions between the factors, whereof total annual precipitation, mean annual temperature, and soil type had greater influences, and most exhibited bi-variable enhancements. The interaction between mean annual temperature and total annual precipitation (Q value = 0.881) specifically had the greatest influence on the spatial distribution of the NDVI, for which the influence of the social and natural factors, especially land use type, livestock density, and soil type factors, increased significantly. The influence of social and natural factors, especially land use type and livestock density, increased significantly with the soil type factor, whereas the population and GDP factors, which were the least influential, showed mostly non-linear enhancements with the natural factors.
Remote Sens. 2021, 13, x FOR PEER REVIEW 21 of 28 the population and GDP factors, which were the least influential, showed mostly non-linear enhancements with the natural factors.

Differences in Influence for Drivers on the NDVI
The ecological detector can reflect significant differences in the effect of factors on the spatial distribution of NDVI. This can further verify the dominant influencing factors and evaluate the variability of their mechanisms of action. The results of the ecological probe ( Figure 12) show that meteorological factors, land use type, elevation, and slope all had significantly different effects, indicating that although these factors play dominant roles in the NDVI, they have different mechanisms for their effects on vegetation growth. Specifically, precipitation and soil type exhibit significant differences in their effects on the NDVI, wherein the effect of soil type on vegetation growth is significant in areas where precipitation is the main driver. Other factors, such as temperature and soil type, do not show significant differences in their effect on vegetation, because while increased temperature enhances photosynthesis, it also leads to faster transpiration by leaves and increased evapotranspiration, thereby reducing soil moisture and inhibiting vegetation growth [75,77]. The NDVI values within the different landscape types did not vary significantly according to landscape changes, and the elevation and slope were not the main factors contributing to the spatial variation in NDVI.

Differences in Influence for Drivers on the NDVI
The ecological detector can reflect significant differences in the effect of factors on the spatial distribution of NDVI. This can further verify the dominant influencing factors and evaluate the variability of their mechanisms of action. The results of the ecological probe ( Figure 12) show that meteorological factors, land use type, elevation, and slope all had significantly different effects, indicating that although these factors play dominant roles in the NDVI, they have different mechanisms for their effects on vegetation growth. Specifically, precipitation and soil type exhibit significant differences in their effects on the NDVI, wherein the effect of soil type on vegetation growth is significant in areas where precipitation is the main driver. Other factors, such as temperature and soil type, do not show significant differences in their effect on vegetation, because while increased temperature enhances photosynthesis, it also leads to faster transpiration by leaves and increased evapotranspiration, thereby reducing soil moisture and inhibiting vegetation growth [75,77]. The NDVI values within the different landscape types did not vary significantly according to landscape changes, and the elevation and slope were not the main factors contributing to the spatial variation in NDVI. Remote Sens. 2021, 13, x FOR PEER REVIEW 22 of 28

Discussion
The NDVI values in Inner Mongolia indicated a generally decreasing trend before 2009, and then increasing, which is consistent with the results of existing studies [78]. Additionally, the sudden change in NDVI in the central part of Genhe city in Hulunbeier between 2002 and 2005 may be attributed to the frequent grassland fire events in the forest areas of Inner Mongolia during this period [79]. Most of the remaining areas north of central Inner Mongolia mutated between 2008 and 2014, and vegetation showed an increased trend after the mutations. This may be attributed to the positive effect of the "Fifth Phase of the Three-Northern Shelter Forest Program" [80]. The breakpoints in the northern parts of Baotou and Wulanchabu occurred between 2014 and 2017. This may be attributed to the government of Baotou that began artificial restoration and the development of heavily polluted mining areas in 2011, for which physical separation and remediation methods were used to remove heavy metals in the soil to comply with the National Environmental Quality Standards for Soil in China [81,82]. Additionally, vegetation in the Siziwang Banner area increased because of the combined effect of increased precipitation and reduced anthropogenic grazing [44]. In the northwestern regions of Ordos, Bayannur, and the Alxa League, the vegetation deteriorated to varying degrees owing to a combination of natural and human factors, such as reduced precipitation, overgrazing, and human reclamation of land. Vegetation around water bodies showed a small improvement. For example, Moon Lake and the surrounding areas of the Alxa League transitioned from a significant decrease to an increase in NDVI before, and after the mutation, respectively [83].
The average temperature in Inner Mongolia increased and precipitation decreased from 2000 to 2009, leading to increased drought in Inner Mongolia and inhibiting vegetation growth [49]. The influence of precipitation on vegetation in Inner Mongolia mainly occurred in areas where agricultural vegetation and grasslands are distributed. Con-

Discussion
The NDVI values in Inner Mongolia indicated a generally decreasing trend before 2009, and then increasing, which is consistent with the results of existing studies [78]. Additionally, the sudden change in NDVI in the central part of Genhe city in Hulunbeier between 2002 and 2005 may be attributed to the frequent grassland fire events in the forest areas of Inner Mongolia during this period [79]. Most of the remaining areas north of central Inner Mongolia mutated between 2008 and 2014, and vegetation showed an increased trend after the mutations. This may be attributed to the positive effect of the "Fifth Phase of the Three-Northern Shelter Forest Program" [80]. The breakpoints in the northern parts of Baotou and Wulanchabu occurred between 2014 and 2017. This may be attributed to the government of Baotou that began artificial restoration and the development of heavily polluted mining areas in 2011, for which physical separation and remediation methods were used to remove heavy metals in the soil to comply with the National Environmental Quality Standards for Soil in China [81,82]. Additionally, vegetation in the Siziwang Banner area increased because of the combined effect of increased precipitation and reduced anthropogenic grazing [44]. In the northwestern regions of Ordos, Bayannur, and the Alxa League, the vegetation deteriorated to varying degrees owing to a combination of natural and human factors, such as reduced precipitation, overgrazing, and human reclamation of land. Vegetation around water bodies showed a small improvement. For example, Moon Lake and the surrounding areas of the Alxa League transitioned from a significant decrease to an increase in NDVI before, and after the mutation, respectively [83].
The average temperature in Inner Mongolia increased and precipitation decreased from 2000 to 2009, leading to increased drought in Inner Mongolia and inhibiting vegetation growth [49]. The influence of precipitation on vegetation in Inner Mongolia mainly occurred in areas where agricultural vegetation and grasslands are distributed. Conversely, the effect of precipitation is less than that of temperature in forests, probably because agricultural vegetation and grasslands are more sensitive to changes in precipitation. The influence of natural factors is also reflected in the vegetation breakpoint years and the changes in their trends before and after the breakpoints. The cumulative effect of precipitation and soil moisture on desert steppe and desert areas in the northwest is more pronounced and lasts longer than that in other areas because of their fragile ecological environment. This makes them less adaptable to changes in external natural conditions. Additionally, the increased frequency and intensity of extreme climatic events, such as droughts and high-temperature heat waves, can lead to significant changes in vegetation cover [77,84]. Soil types in Inner Mongolia are distributed from northeast to southwest in the order of black soil, dark brown loam, black calcium soil, chestnut calcium soil, brown loam, black kiln soil, gray calcium soil, wind-sand soil, and gray-brown desert soil, and the stability of vegetation growth varies greatly among different soil types. Black and blackish-chestnut soils in the central and northeastern parts of the country have high contents of organic matter, strong water retention properties, and good natural conditions, which are conducive to the growth of vegetation. The stability of vegetation in this area is correspondingly high, with the occurrence of fewer breakpoints. Meanwhile, the eolian sandy soil in the northwest gray desert zone has low soil fertility and poor soil texture, making it susceptible to changes in external conditions. Therefore, this is the main area where sudden changes in vegetation occur and the vegetation is less stable [10].
Inner Mongolia is one of the regions currently experiencing significant environmental changes as a result of tremendous economic development, including changes and intensification of grazing systems and increased intensity of external anthropogenic activities, such as agricultural production. In this study, we highlighted the contribution of anthropogenic factors to vegetation growth, and found that land use types strongly influence vegetation in Inner Mongolia. The impact of human activities is more spatially concentrated, with human activities leading to changes in land-use types and shifts in land-use types affecting the growth of vegetation. This leads to the spatial differentiation of vegetation, as evidenced by the effects of policy implementation, such as returning farmland to forest and grazing bans. Since the late 1990s, the Chinese central and local governments have implemented and introduced large-scale ecological restoration projects, including the "Returning Farmland to Forests Project," the "Returning Grassland to Grass Ecological Project in Northern China," and the "National Ecological Protection." These programs have had a positive impact on vegetation growth in Inner Mongolia owing to the implementation of specific measures to protect grasslands and vegetation through rotational grazing, grazing bans, and converting pasturelands to croplands.
NDVI is the most widely used observation method for vegetation. NDVI is a good proxy for vegetation density parameters such as leaf area index (LAI), vegetation cover (FVC), and absorbed photosynthetically active radiation (fAPAR). However, NDVI has two major limitations in characterizing biomass and productivity. First, the relationship between NDVI and green biomass is non-linear and can be saturated in areas with high vegetation cover. The second limitation is that NDVI mainly reflects vegetation greenness rather than photosynthesis itself. However, total primary productivity (GPP) can decline without any reduction in LAI or chlorophyll. Combining the shortcomings of NDVI, a new vegetation index, kNDVI, was proposed by applying the theory of nuclear methods to NDVI, using machine learning. The kNDVI correlated with GPP similarly or better than other indices globally. kNDVI correlated with SIF better than other indices in general and in all biomes, especially in deciduous broadleaf forests and for herbaceous and cultivated crops. The correlations of kNDVI were higher in nearly all cases (e.g., Spearman correlation, distance correlation), thus confirming the advantage of kNDVI over other indices. Using the kNDVI index in future studies on vegetation will increase the significance of geomonitoring and terrestrial biosphere studies [85].
The use of MODIS NDVI data as the only input in this study may reduce the accuracy of the vegetation breakpoint results. This limitation can be addressed in future studies by incorporating data sources with higher temporal resolution and finer spatial resolution.
For example, the data provided by the GF-1 and Sentinel-2 satellite multispectral scanners could be combined into more comprehensive driving force data. The spatial pattern driving mechanism of vegetation change should be further explored in conjunction with climate change and ecological protection policies, to provide more scientific data support for improving the ecological environment in Inner Mongolia.

Conclusions
Considering Inner Mongolia as the study area, and MODIS-NDVI, total annual precipitation, DEM, mean annual temperature, and data on nine other drivers, this study explored the spatiotemporal characteristics and drivers of NDVI change in Inner Mongolia using spatiotemporal abrupt change and trend analyses, significance tests, the Hurst index, and the GeoDetector model. The spatiotemporal evolution of the vegetation in Inner Mongolia was characterized, and the trend of breakpoints was used to predict the possible future direction of NDVI change, revealing the driving mechanisms and dominant drivers of the NDVI. The detailed conclusions are as follows: Precipitation, soil type, temperature, and land use type were the main driving factors of NDVI changes in Inner Mongolia, with precipitation being the most influential factor. Therefore, areas with the highest total annual precipitation had the most significant impact on NDVI changes, and the soil and land use types with the largest areas in Inner Mongolia had significant impacts on the NDVI. Nevertheless, the interactions between factors exhibited mostly bi-variable enhancements. Although meteorological, land use type, and topography factors dominantly influenced vegetation growth, their mechanisms of influence showed significant differences.