Grassland Health in Xilin Gol League from the Perspective of Machine Learning—Analysis of Grazing Intensity on Grassland Sustainability

Abstract

As one of the most widespread and important types of terrestrial vegetation in the world, grasslands play an irreplaceable role in global climate change. The grasslands of Inner Mongolia, represented by the Xilin Gol League, are typical of Eurasian grasslands and have an important ecological status in the world. In this paper, taking the grassland of Xilin Gol League as the research object, based on the machine learning method, we mainly carry out two aspects of work: the prediction of grassland soil health and evaluation of grassland sustainable development. To address the issue of predicting soil health in grasslands, we focus on an important indicator in grasslands: soil moisture. By analyzing the characteristics of soil moisture time series values and related influencing factors, based on a NAR neural network model, three important factors of soil moisture were predicted: soil evaporation data, average air temperature, and precipitation. Subsequently, the corresponding soil moisture calculation model was trained using regression models based on hyperparameter optimization, and the final predicted soil moisture values were obtained for different months and depths in 2023 and 2024. To evaluate the sustainability of grassland development, we developed a model for the degree of grassland desertification based on the kernel principal component analysis, focusing on three dimensions: environmental factors, surface factors, and human factors. Based on this, a quantitative definition of soil denudation is given by analyzing the main influencing factors of grassland soil degradation. At the same time, a prediction model for the evaluation of soil slumping was established based on a fuzzy comprehensive evaluation matrix, and the evaluation weights of each major factor were given and analyzed. Based on the above research, this paper suggests a reasonable grazing strategy for the grassland areas of the Xilin Gol League: when the grazing intensity is medium and the total number of grazing days is [85, 104] days in a year, the degree of soil slumping and soil desertification in the pastures is minimized. The research results of this paper are useful for the future maintenance and management of the grasslands of Xilin Gol League and other similar areas.

1. Introduction

As one of the most widely distributed and important types of terrestrial vegetation in the world, grasslands are an important functional component of terrestrial ecosystems and play an irreplaceable role in climate change globally, supporting more than 1 billion people worldwide and serving as shelter for a broad variety of animals and plants [1]. The grasslands of Inner Mongolia in China, represented by the Xilin Gol league [2], are typical of Eurasian grasslands and enjoy a significant ecological status in the world [3]. In recent years, however, with the intensification of extreme climatic phenomena represented by soil erosion and dust storms, rational grazing policies have received increasingly widespread attention [4], and research into a range of optimization issues related to grazing will provide a scientific basis for the country to make rational grassland management decisions [5]. In practical grassland management and production, there are often instances where grazing plans are sketchy, production details are vague, and feeding efficiency is low [6]. It is, therefore, necessary to develop a more systematic, scientific, precise, and detailed model for the sustainable use of grassland grazing and ecological conservation that is universally applicable to the typical Eurasian grassland region represented by the Xilin Gol league.

The level of sustainable grassland development is determined by a combination of factors [7], among which soil moisture is one of the more important indicators and a reliable basis for judging the ecological health of a region [8], as previous studies have proved that precipitation, evapotranspiration (ET), and soil moisture are the key controls for the productivity and functioning of temperate grassland ecosystems in Inner Mongolia, northern China [9]. Precipitation dominated the water balances in the grassland, and the balance between ET and precipitation explained the seasonal changes in soil water content in different soil layers [10]. Therefore, the prediction of future values of soil moisture in grasslands [11] plays a very important role in this study. Other parameters such as soil nutrients, ambient temperature, grazing intensity, and intergrass competitions also play key roles in the reflection of grassland health, however, to highlight the significant role of precipitation in the sustainable development of grasslands in typical Eurasian steppe regions and to demonstrate the feasibility of a machine-learning-based approach to grassland ecological monitoring, this paper focuses on one important indicator in predicting soil health in grasslands: soil moisture.

However, soil moisture time series have the characteristics of nonlinearity, nonsmoothness, large noise, and large variation in associated risk, which increases the difficulties of soil moisture time series prediction [12,13]. In recent years, a variety of prediction methods have been studied to solve the soil moisture time series forecast problem. Ref. [13] utilizes a Long Short-Term Memory Model to harness newly available satellite-sensed data for long-term soil moisture projections. Ref. [14] applies three imputation methods, linear interpolation, multiple linear regression, and extended ARMA models with exogenous climatic variables for soil moisture time series forecast in the Kulunda Steppe in the Altai Krai, Russia. Ref. [15] adapts the Autoregressive Integrated Moving Arima (ARIMA) model for soil moisture forecasting in Mongolia, using SMAP and MODIS satellite data. However, there remain a few drawbacks to the above methods, such as overfitting problems. Based on this, there is an urgent need to study models and methods for predicting key indicators for typical Eurasian steppe regions.

On the other hand, under the influence of global climate change and human activities, grassland ecosystems are undergoing severe degradation, undermining their capacity to support biodiversity, ecosystem services, and human well-being, which is manifested by the homogenization of grassland plant populations, the reduction in biodiversity, soil desertification, and denudation [16]. Desertification and denudation, as extreme forms of grassland degradation, will not only lead to the formation of sandy landscapes but also the coarsening and impoverishment of grassland soils and the reduction in grassland production potential [17,18]. Based on this, it is necessary to provide indicators for the evaluation of the sustainability of grassland development [19,20], represented by soil consolidation and grassland desertification, to provide a scientific basis for the designation of grazing strategies and the monitoring of grassland health [21].

The main contributions of this paper are summarized as follows.

(i)

Based on dynamic time series prediction and support vector regression, using soil moisture, soil evaporation, and precipitation data for the grasslands in the Xilin Gol region from 2012 to 2022, a model was developed to predict soil moisture at different depths in 2023 and 2024, while keeping the current grazing strategy unchanged. Firstly, a time series neural network was used to predict the soil evaporation data, average temperature, and precipitation in 2022 and 2023, and then a Gaussian process regression was used to train the corresponding soil moisture calculation model to finally obtain soil moisture data for different months in 2022 and 2023.

(ii)

Based on the kernel principal component analysis, and concerning the green cover, vegetation index, and runoff volume of the Xilin Gol League in the past ten years, a model of the degree of grassland desertification was developed with three dimensions: environmental factors, surface factors, and human factors. Based on this, a quantitative definition of soil denudation was given by analyzing the main influencing factors of grassland soil degradation; meanwhile, a prediction model for soil denudation evaluation was established based on the fuzzy comprehensive evaluation matrix, and the evaluation weights of each main factor were given and analyzed.

(iii)

Based on the above analysis, a comprehensive analysis of the grassland ecological data and herders’ livestock data in the Xilin Gol League region in the past ten years was conducted to give a more reasonable grazing strategy in the future: when the grazing intensity is medium and the total number of grazing days is within the range of $[85,104]$ days per year, the degree of soil slumping and soil desertification in the pasture is minimized and the indexes of soil moisture are satisfied.

(iv)

Based on machine learning methods and fuzzy evaluation matrices, this paper makes preliminary predictions and analyses of the future health and sustainability of the Eurasian steppe region represented by the Xilin Gol league.

(v)

The research results of this paper have certain guiding significance and reference value for the future sustainable development of grassland ecology and livestock policy formulation in the Xilin Gol league grassland and similar areas.

2. Material and Methods

2.1. Study Site

The research was carried out from January 2013 to December 2022 in the Xilin Gol Meng region of Inner Mongolia, China. As shown in Figure 1, the coordinates of the Xilin Gol grassland range from 110°50′ E to 119°58′ E and 41°30′ N to 46°45′ N. The average annual precipitation is 340 mm, with most of the annual rainfall concentrated between June and August. It is a representative and typical grassland in the temperate zone, which belongs to the northern temperate continental climate. The main climatic characteristics are windy, dry, and cold, with an average annual temperature of 0–3 °C, a freezing period of up to 5 months, and a cold period of up to 7 months, with the lowest temperatures in January and the highest temperatures in July. The specific month-by-month data on climatic characteristics (precipitation and average temperature) for the past 5 years are illustrated in Figure 2 and Figure 3.

The number of grazing days in the region varies from 70 to 300 days. The region is one of the most concentrated and severely damaged areas of desertification and sandy land in China, with the four deserts of Badangilin, Tengri, Ulanbuhe, and Kubuqi and the four grains of sand of Mawusu, Hunsandak, Horqin, and Hulunbuir. Desertified land accounts for 51.50% of the total area of the region, and sandy land accounts for 34.48% of the total area of the region, with a very fragile ecological base. In the late 1980s, the area of sandy land in the region accounted for 10.12% of the total land area of the region, accounting for 51.47% of the region’s desert and sandy land area. Among them, mobile sand, semifixed sand, and fixed sand accounted for 18.51%, 35.01%, and 46.48% of the total sand area, respectively. According to the definition of desertification set out in the United Nations Convention to Combat Desertification, the potential area of land desertification in the region is 701,100 square kilometers, accounting for 59.30% of the total land area of the autonomous region. Of these, the arid, semiarid, and dry subhumid zones cover an area of 170,100 square kilometers, 267,300 square kilometers, and 263,700 square kilometers, respectively. The overall trend of desertification in the region is gradually increasing from east to west, with lightly desertified land mainly distributed in the subhumid zone in the east, moderately desertified land mainly distributed in the semiarid zone in the center, and very heavily desertified land mainly distributed in the western part of the autonomous region. From 2012 to 2022, an average of 12 million mu of land was completed annually in the region for sand prevention and control.

2.2. Description of the Grazing Strategy Selected for this Study

There are currently five main types of grazing in typical Eurasian grasslands: continuous grazing (CG), grazing prohibition (NG), rotational grazing (RG), light grazing (LG), and seasonal rest grazing (SRG) [22,23]. A diagrammatic representation of the different grazing methods is shown in Figure 4. The strategy of RG is to divide the area into several zones of approximately equal size and to alternate grazing within these zones, with each zone being grazed for approximately one month. Previous studies have shown that grassland communities will maintain high biodiversity under moderate grazing intensity, while with changes in grazing practices and increasing grazing intensity, grasslands will gradually degrade [24,25,26]. The vegetation biomass and Shannon–Wiener diversity indices of grassland plant communities are the highest in zoned rotational grazing and lowest in CG [27]. NG also leads to a homogenization of the community structure, which is not conducive to the maintenance of plant diversity [28].

Furthermore, multiple studies show that managing grazing for high livestock density has certain effects on plant and livestock productivity, or may even promote biodiversity in the region [29,30]. Ref. [31] conducted a study in the northern slopes of Mount Kilimajaro in order to examine the impact of grazing on natural regeneration of the grazed vegetation, pointing out that overgrazing has ecological ramifications which lead to degradation of the ecosystem, which has become a problem in many parts of Tanzania. Refs. [32,33,34] proved that grasslands in Tanzania/Kenya with high-density grazing and 1 year resting could show high grassland biodiversity. Ref. [35] established phenological profiles of differing vegetation in the Kenya–Tanzania borderlands and characterized the trends and drought resilience across the landscape through analyzing MODIS vegetation indices in 8-day time steps from 2002 to 2012.

In conclusion, RG is a necessary condition for the sustainable use of grassland resources and the healthy development of animal husbandry. Therefore, according to the actual grassland animal husbandry and the needs of the research problem, in this study, the grazing strategy of RG is mainly considered.

2.3. Data Collection

Basic data on the soils and climate of the Xilin Gol grasslands in Inner Mongolia are monitored and provided by specialized institutions, especially core data such as soil moisture, soil evaporation, and vegetation index NDVI. In addition, in recent years, a large number of experiments have been conducted on representative grasslands in Xilin Gol league, Inner Mongolia Autonomous Region, and data have been collected on sample community surveys of different herding households’ ecological animal husbandry patterns; livestock numbers of different model herding households; soil carbon and nitrogen monitoring data of different grazing intensities; community structure monitoring data of rotational grazing sample plots; and livestock numbers of different model herding households in Xilin Gol league. The data were also obtained.

At the same time, data on the community structure of typical grassland rotational grazing sites in the Xilin Gol League were also obtained. These data were collected from 12 grazing plots, and each plot was sampled in five 1m X 1m squares according to the initial sampling zone. The experimental design was a randomized grouping with three replications of each grazing intensity. The grazing experiment began in 2013 with four rotations per year, and at the end of each grazing period, 12 grazing plots were sampled, with five 1 m × 1 m samples collected from each plot according to the initial sampling strips. The above-ground plant parts were harvested flush using the harvesting method, and the collected plants were weighed fresh on the same night and placed in a blast oven at 65 °C. After 24 h, the dry weight was weighed to determine the above-ground plant biomass.

The above data are publicly available at http://www.ncdc.ac.cn/portal/ (accessed on 1 January 2023). and https://ladsweb.modaps.eosdis.nasa.gov/search/ (accessed on 4 January 2023). The data sets used for this study are available on request from the research team.

2.4. Assumptions and Notations

From the perspective of mechanism analysis, based on data processing, model building, and model prediction, the following assumptions need to be made in this paper:

(i)

The existence of grassland requires annual precipitation greater than a certain threshold, which is set at 300 mm in this paper;

(ii)

There is a correspondence between grass production and annual precipitation;

(iii)

Overgrazing leads to desertification;

(iv)

The duration of grazing and grazing rest is a full month;

(v)

Moderate grazing means that the amount of pasture grazed is in a dynamic balance with the carrying capacity of the pasture, thus maintaining the normal production of livestock and the healthy and sustainable development of the grassland in the region;

(vi)

A standard sheep unit is a unit of calculation for livestock. A sheep unit is 1.8 kg of hay with 14% water and 10% crude protein per day, weighing 50 kg for a lambing ewe. One sheep equals one sheep unit, while one cow equals five sheep units, and one horse, donkey, and mule each equals five sheep units.

The main symbols used in this paper are illustrated below in

Table 1. Main symbols used in this study.

Symbols	Description
CG	Continuous grazing
NG	Grazing Prohibition
RG	Rotational Grazing
LG	Light Grazing
SRG	Seasonal Rest Grazing
CK	Zero grazing intensity
LGI	Light grazing intensity
MGI	Moderate grazing intensity
HGI	Heavy grazing intensity

.

2.5. General Framework of This Study

The overall framework of this study is shown in Figure 5. To address the issue of predicting soil health in grasslands, we focus on an important indicator in grasslands: soil moisture. By analyzing the characteristics of soil moisture time series values and related influencing factors, based on a NAR neural network model, three important factors of soil moisture were predicted: soil evaporation data, average air temperature, and precipitation. Subsequently, the corresponding soil moisture calculation model was trained using Gaussian process regression, and the final predicted soil moisture values were obtained for different months and depths in 2023 and 2024. To evaluate the sustainability of grassland development, we developed a model for the degree of grassland desertification based on the kernel principal component analysis, focusing on three dimensions: environmental factors, surface factors, and human factors. Based on this, a quantitative definition of soil denudation is given by analyzing the main influencing factors of grassland soil degradation. At the same time, a prediction model for the evaluation of soil slumping was established based on a fuzzy comprehensive evaluation matrix, and the evaluation weights of each major factor were given and analyzed. Based on the above research, this paper suggests a reasonable grazing strategy for the grassland areas of the Xilin Gol League.

2.6. Soil Moisture Prediction in the Grasslands of Xilin Gol League Based on NAR and SVM

In this study, based on the soil evapotranspiration, mean temperature, and precipitation for the previous 10 years in the Xilin Gol League region, a time series neural network model was selected to predict the three indicators for 2023 and 2024. The prediction can be understood as a nonlinear autoregressive model with no external input, so the NAR (Nonlinear Auto Regressive) neural network model was chosen for the prediction.

In a time-series-based dynamic neural network, the output of the model at each moment is based on the synthesis of the dynamic results of the system before the current moment and, therefore, has a feedback and memory function. This prediction model can be seen as a special case of nonlinear system identification. Let $\left\{y_i^k\mid i=1,\dots ,N\right\},k=1,2,3$ be a time series of local environmental data, where N is the length of the time series and k is the number of the different environmental data. Using the variable $\left\{y_i^k\mid i=1,\dots ,N\right\}$, a nonlinear autoregressive prediction model can be obtained as follows:

$$y_t=f\left(y_{t-1},y_{t-2},\dots ,y_{t-m}\right)$$

(1)

where the function f was estimated by the NAR, with order m determined from the experimental data. After estimating the models, they were used in a simulation model in which past predictions were used and future predictions were calculated using the estimated model f, as shown below.

$${\widehat{y}}_t=\widehat{f}\left({\widehat{y}}_{t-1},{\widehat{y}}_{t-2},\dots ,{\widehat{y}}_{t-m}\right)$$

(2)

In this paper, a NAR model is constructed using the technique of OLS, a simultaneous structure selection parameter estimation algorithm for a variety of nonlinear system identification models. For a given time series output, the polynomial of the NAR model is expressed as

$$y\left(t\right)=\sum _{m=1}^{n_p}{P}_m{\theta }_m+\xi \left(t\right)$$

(3)

where $n_p$ denotes the number of polynomial expansions, $P_m$ is the mth regression term, which consists of lagged combinations of the output terms, and ${\theta }_m$ is the mth regression coefficient. In this paper, the number of lagged terms is set according to the characteristics of each data set, as temperature, soil evaporation, and precipitation are usually seasonal events that repeats each year and are therefore designed according to the different tracks.

In matrix form, identification involves the formulation and solution of the least squares problems as follows:

$$P\theta +\zeta =y$$

(4)

where P is the $n\times m$ order regression matrix, $\theta$ is the $m\times l$ dimensional coefficient vector, and y is the $n\times l$ dimensional actual observation vector. The OLS technique will perform both structure selection and parameter estimation based on the QR decomposition as follows:

$$P=QR$$

(5)

Define g as

$$g=Q^{T}y$$

(6)

Thus, it can be derived as

$$R\theta =g$$

(7)

Each column in P is calculated as a measure of $ERR$, and the $ERR$ value represents the contribution of each regression column in reducing the error variance of output $\widehat{y}$ and observation y. Based on g, the ERR value for each i row yields

$$ER{R}_i=\frac{g_i^2}{y^{T}y}$$

(8)

The model structure is determined by selecting regression columns until the cumulative $ERR$ value reaches a sufficiently high value. In the experiments in this paper, $0.98$ is used as the threshold, which means that the chosen regression volume should account for at least 98 percent of the reduction in the error variance. The parameter value E for the reduced P matrix ($PR$) can be re-estimated again using the $QR$ decomposition (Equation (5) to Equation (7)) to obtain the final model.

Finally, by rearranging and solving Equation (7), the value of ${\theta }_R$ can be estimated as

$${\theta }_R=R_R^T{g}_R$$

(9)

The obtained candidate models were validated and analyzed to determine the best model. One-step ahead (OSA), correlation, and histogram tests were used to select the best model that met the validation criteria. Meanwhile, parameters such as precipitation were used to solve for soil moisture and, to pursue maximum rendezvous coverage, the target model for this problem was defined by Equation (10) as follows, which is a regression fitting problem, and three indicators such as soil evaporation rate, average air temperature, and precipitation were selected as model input features:

$$g=g(h,i,\Omega )$$

(10)

where g indicates soil moisture.

Support vector machines (SVM) are mathematical models and powerful general-purpose approximators that can be used for both classification and regression tasks and are recognized as one of the best classifiers that can be implemented effectively to date. Assuming a data set ${\left\{x_i,y_j\right\}}_{i=1,\dots ,n}$, where $x_i\in {\Re }^{\mathrm{D}}$ denotes the input features, d is the number of features in the training data set, and $y_i\in \Re$ is its corresponding actual output, the main goal of the SVM algorithm is to draw the linear decision function shown as follows:

$$f\left(x\right)=<w,{\varphi }_i\left(x\right)>+b$$

(11)

where w and b are the weights and constants, respectively, that must be estimated from the data set. $\varphi$ is a nonlinear function that maps the input features to a higher feature space. $<.,.>$ denotes the dot product in ${\Re }^D$. This problem can be expressed as minimizing the following function:

$$R\left(C\right)=\frac{C}{n}\sum _{i=1}^n{L}_{\epsilon }\left(f\left(x_i\right),y_i\right)+\frac{1}{2}\Vert w{\Vert }^2$$

(12)

where $L_g\left(f\left(x_i\right),y_i\right)$ is the $\epsilon$ density loss function, and such optimization problems with specified constraints can be treated as quadratic optimization problems using Lagrange multipliers.

An optimizable SVM model is a regression model based on hyperparameter optimization, which is trained by adjusting the hyperparameters to optimize the hyperparameters of the ensemble to achieve better training results. By choosing a cross-validation method, this method can better prevent overfitting and is suitable for small- and medium-sized data sets, which is the case in this study. To compare the advantages and disadvantages of the models used, the steps shown in Algorithm 1 were used for training.

Algorithm 1: Regression model training process based on hyperparameter optimization

Step 0: Prevent overfitting by dividing the track data set into folds and estimating the accuracy of each fold. The number of folds in this model is chosen to be 5 because of the small amount of data in this study. Step 1: Select a general predefined training model to facilitate later comparison with the hyperparametric optimization model. Step 2: Initialize the parameters, such as the basis functions. Step 3: Perform model training. Step 4: Calculate fitness: the standard SVM model is trained on the new training set $trai{n}_{s}et*$, and the validation result on the validation set $validatio{n}_{s}et*$ is taken as the fitness value, with the global optimal fitness value as $globa{l}_{f}itness$, the individual optimal fitness value as $loca{l}_{f}itness$, the global optimal particle position as $globa{l}_x$, the individual optimal particle position as $globa{l}_x$, and the individual historical optimal position as $loca{l}_x$. Step 5: Compare the predefined training model with the model based on the adjusted hyperparameters. Step 6: Obtain the optimal model and use it for computation.

2.7. Desertification and Slabbing Analysis of Grassland Soils Based on KPCA and Fuzzy Evaluation Matrix

The prediction models for the degree of desertification index and the degree of the slabbing index were first developed by using kernel principal component analysis to obtain the weight coefficients of the different factors in this model. After the data were collated, the modeling was carried out for historical data on nine main factors, including environment, NAVI index, runoff, soil moisture, the population of Xilin Gol Union, livestock, and income of urban and rural people in the past years, respectively, for the 2013–2022 period. In selecting the data for quantitative analysis, the soil moisture data were mainly selected for the 10 cm depth of soil, as soil slabbing is mainly related to surface soil moisture. As the sampling period for the different data in the raw data is not consistent, the data are averaged by year when making predictions.

The PCA method usually completes the comprehensive evaluation of the model performance based on the individual performance indexes of each dimension. The contribution rate of each principal component reflects the amount of information extracted from the initial individual index, the larger the contribution rate, the more original information contained in the principal component. However, when the contribution rate of the first principal component is less than 85 percent, the data of individual performance indexes are not sufficiently compressed, in which case it is difficult to extract the nonlinear characteristics of individual performance indexes and, therefore, the KPCA method is applied to the comprehensive evaluation of model performance. The initial single-indicator data can be mapped to a high-dimensional space by a nonlinear transformation, while effective dimensionality reduction is calculated by principal component analysis.

The main idea of KPCA is to replace the original data multidimensional variables with some new independent principal components based on linear transformation with minimal loss of information. The nature of the original sample can be represented by principal components with some significance. For example, changes in the original data are reflected by the first principal component, and other characteristics of the sample can be reflected by other principal components. The main steps of the principal component analysis used in this study are illustrated in Algorithm 2.

Algorithm 2: Main steps of the principal component analysis

Step 0: The raw data are standardized to eliminate the adverse effects caused by different dimensions, and the raw data can be expressed as follows: $$X={\left(x_{ij}\right)}_{n\times p}=\left(X_1,X_2,\cdots ,X_p\right)$$ (13) The standardization formula is shown below. $$Z_{ij}=\left(x_{ij}-\overline{x_j}\right)/S_j,i=1,2,\dots ,n;j=1,2,\dots ,p$$ (14) Step 1: Calculate the correlation coefficient matrix as follows: $$R=\frac{1}{n-1}{Z}^{\prime }Z$$ (15) where R is an $n\times n$ symmetric matrix with data on the diagonal all being 1; $Z^{\prime }$ is the transpose matrix of the matrix Z. Step 2: Calculate the eigenvalues and eigenvectors of the correlation coefficient matrix. The eigenvalue ${\lambda }_i$ is obtained from $|\lambda E-R|=0$ and classifies ${\lambda }_i$ according to its size, ${\lambda }_1\ge {\lambda }_2\ge \cdots \ge {\lambda }_p\ge 0$; the eigenvector $Z{X}_i$ is obtained from $(\lambda E-R)X=0$, i.e. $Z{X}_i={(Z{X}_{i1},Z{X}_{i2},\cdots ,Z{X}_{ip})}^{\prime }$. Step 3: Determine the principal component contribution as $${\alpha }_i=\frac{{\lambda }_i}{\sum {\lambda }_i}$$ (16) The cumulative contribution of the first m principal components is (in general, the first m principal components were selected such that the cumulative contribution $C{V}_i\ge 85\%$) $$C{V}_i=\frac{{\displaystyle \sum _{i=1}^m}{\lambda }_i}{{\displaystyle \sum _{i=1}^p}{\lambda }_i{\alpha }_i}$$ (17) Step 4: Determine the expression of principal components: $$F_i=Z{X}_{i1}\times Z_1+Z{X}_{i2}\times Z_2+\cdots +Z{X}_{ip}\times Z_p$$ (18) Step 5: Determine the comprehensive evaluation function: $$F=\left({\lambda }_1{F}_1+{\lambda }_2{F}_2+\cdots +{\lambda }_m{F}_m\right)/\left({\lambda }_1+{\lambda }_2+\cdots +{\lambda }_m\right)$$ (19)

3. Results

3.1. Soil Moisture Prediction in 2023 and 2024 in the Xilin Gol Region

Based on the soil moisture prediction model in the grasslands of Xilin Gol League based on NAR and SVM, the following neural network training results were obtained for the three inputs for the next two years (2023 and 2024) predictions. As each value has different characteristics, and therefore the same neural network could not be applied, the training and prediction results for the three inputs are shown below. All the plots below on the performance of the prediction models are arranged in the order of soil evaporation rate, average temperature, and precipitation. In particular, Figure 6 shows the variation in MSE values for the three networks.

Figure 7 shows the histograms of the prediction errors of soil evaporation rate, average temperature, and precipitation for the next two years (2023 and 2024). It can be found that the error distribution of the three variables is relatively concentrated, basically all concentrated in the vicinity of $Errors=0$, which can infer that the overall error of the prediction model is small and the prediction effect is satisfactory.

Figure 8 shows the data fit of the training, validation, test, and overall data. It can be found that the fit of all three sets of data is reasonable, proving that the model can be considered to have good training results, which could be used for the prediction of more model variables.

Figure 9, Figure 10 and Figure 11 show the response of the time series of the prediction model and its prediction results, from which it can be seen that the error distribution is relatively uniform and can be used in the prediction of the model. In Figure 11, the precipitation data are highly variable and have a different distribution of regular characteristics, with different climatic characteristics from year to year; for example, some years have a concentration of rainfall in the middle of the year, while others have a concentration of rainfall at the beginning and end of the year.

Figure 12, Figure 13 and Figure 14 show the distribution values of the error autocorrelation of the prediction model.

Figure 15 shows a plot of the historical data for the previous ten years (2013–2022) for the three inputs, as well as the predicted values for 2023 and 2024. It can be seen that the patterns are relatively similar between the years before and after the three sets of data, and therefore the forecasts are considered to be good.

As can be found in the previous Figure 3, although the annual precipitation in the region is below 600 mm in most years, there are still some years in the last decade when the precipitation exceeds 1000 mm. In 2018, it even exceeded 1200 mm. It is with this in mind that, based on thorough research and studies, this paper examines soil moisture after 100 cm and up to 200 cm. In addition, although the soil moisture data above 100 cm remained stable in most years, they still contrasted sharply with the soil moisture data below 100 cm. Therefore, the study of soil moisture data above 100 cm is still of interest in this paper.

Next, the training effect of the SVM model used in this study is shown in Figure 16. In the prediction process for soil moisture in 2022 and 2023, the moisture at a depth of 200 cm was excluded from the prediction range as there was little variation between years. It was thus predicted using a linear regression equation.

Figure 17 shows the trend in predicted soil moisture at the depths of 10 cm, 40 cm, and 100 cm, respectively, with the blue section showing the historical data values for the last 10 years of soil moisture data and the red dashed line showing the predicted values based on the model training. It can be seen that as the soil depth increases, the soil moisture increases, and the trend is more pronounced over time.

Finally, the predicted soil moisture in the Xilin Gol grassland for different months in 2023 and 2024 at soil the depths of 10 cm, 40 cm, 100 cm, and 200 cm, respectively, while keeping the current grazing strategy unchanged, are shown in

Table 2. Predicted soil moisture in the Xilin Gol grassland for different months in 2023 and 2024 at soil depths of 10 cm, 40 cm, 100 cm, and 200 cm.

Date	10 cm Humidity (kg/m2)	40 cm Humidity (kg/m2)	100 cm Humidity (kg/m2)	200 cm Humidity (kg/m2)
01.2023	13.37525757	52.22381559	107.6804601	163.51
02.2023	11.83527446	50.65322371	97.49117155	163.51
03.2023	16.82404328	49.46654864	117.4964342	163.48
04.2023	17.00713772	49.59997287	124.5826714	163.42
05.2023	25.75868244	51.18402058	125.0665499	163.21
06.2023	26.00274637	52.32363466	122.5816529	163.88
07.2023	24.85162543	58.63120507	117.6804601	162.91
08.2023	26.65417775	59.54263952	97.49117155	162.90
09.2023	23.30849528	49.48736328	117.4964342	162.85
10.2023	14.40577771	50.93954168	125.5826714	162.56
11.2023	13.33229795	52.35459894	115.0665499	162.42
12.2023	13.41870855	51.06605715	122.5816529	162.17
01.2024	12.45	52.14	93.45	164.48
02.2024	12.1	52.14	93.45	164.48
03.2024	14.96	52.13	93.44	164.48
04.2024	13.32873198	50.23969834	94.58267136	164.48
05.2024	15.34610466	49.37208641	95.06654987	164.45
06.2024	20.56577209	49.92429651	112.5816529	164.41
07.2024	24.94576191	56.89172347	107.6804601	164.22
08.2024	23.32646163	52.06632159	97.49117155	164.15
09.2024	14.4444616	59.44178985	107.4964342	163.86
10.2024	12.21150864	49.40986908	101.5826714	163.71
11.2024	11.05478603	49.74908927	105.0665499	163.56
12.2024	12.81034713	51.44057073	102.5816529	163.51

.

3.2. Degree of Grassland Desertification and Quantitative Definition of Soil Denudation of Xilin Gol League

The model was developed using principal component analysis, and the total variance interpretation and the rotated component matrix are shown in

Table 3. Total variance interpretation of the soil desertification index in the Xilin Gol grassland.

Components	Initial Eigenvalues			SS of Extracted Loads			SS of Rotational Loads
	Total	Variance (%)	Cumulative (%)	Total	Variance (%)	Cumulative (%)	Total	Variance	Cumulative
1	3.009	33.428	33.428	3.009	33.428	33.428	2.914	32.380	32.380
2	2.170	24.110	57.538	2.170	24.110	57.538	2.072	23.024	55.404
3	1.478	16.422	73.960	1.478	16.422	73.960	1.670	18.556	73.960

and

Table 4. Rotated component matrix of the soil desertification index in the Xilin Gol grassland.

Factor	Ingredients
	1	2	3
Average temperature	0.099	0.795	0.408
Precipitation	0.056	−0.093	−0.481
Average wind speed	−0.121	0.032	0.868
Vegetation cover	0.045	0.924	0.052
Surface water resources	−0.138	0.424	0.455
Groundwater resources	−0.961	−0.018	0.132
Population	0.965	−0.030	−0.061
Livestock	0.208	−0.629	0.533
Income per capita	0.983	−0.019	−0.057

. According to the total variance interpretation, the weight ratio of each component can be obtained, and based on this, combined with the rotated component matrix, the weight coefficients of each factor can be obtained.

The final weight coefficients of each factor are shown in

Table 5. Final weight coefficients of each factor of the soil desertification index in the Xilin Gol grassland.

Environmental factors	0.27512341	Average temperature $W_{c_1}$	0.129213773
		Precipitation $W_{c_2}$	0.05623126
		Average wind speed $W_{c_3}$	0.08967837
Surface factors	0.33628267	Vegetation cover $W_{c_4}$	0.106873421
		Surface water resources $W_{c_5}$	0.098943721
		Groundwater resources $W_{c_6}$	0.130465571
Human factors	0.38859344	Population $W_{c_7}$	0.126242129
		Livestock $W_{c_8}$	0.135366619
		Per capita income $W_{c_9}$	0.126985136

.

As the principal component analysis method uses normalized data when calculating the weights of individual indicators, the normalized standard data are used in the subsequent data entry to calculate the desiccation index:

$$SM=\eta ·\sum _{i=1}^n{S}_{Q_{\mathrm{i}}}=\eta ·\sum _{i=1}^n(Q_i·W_{c_i})$$

(20)

where $i=[1,9]$, $\eta$ is adjusted for different indices of desertification according to [36,37], and SM is based on data from the study area (in this study, the Xilin Gol League), which is divided according to the criteria of [38], as shown in

Table 6. Criteria for classifying the degree and index of desertification in typical steppe areas of Eurasia.

Content of Classification	Type and Value of Desertification
Degree of Desertification	Non	Light	Moderate	Severe	Extremely
Desertification index	$\left[0,0.2\right]$	$\left[0.2,0.4\right]$	$\left[0.4,0.6\right]$	$\left[0.6,0.8\right]$	$\left[0.8,1.0\right]$

.

According to the data in [39], the average values of the desertification index of the typical steppe areas of Eurasia in 2021 and 2022 are $0.25165$ and $0.230493$, respectively. Therefore, by correcting $\eta$ by the desertification index value, we finally obtain $\eta =0.5535$. Thus, the evaluation index of the desertification index of the Xilin Gol League is

$$SM==0.5535·\sum _{i=1}^n(Q_i·W_{c_i})$$

(21)

Based on the data in [40,41,42], the NDVI values at different grazing intensities and their accuracy are shown in

Table 7. Grazing intensity and NDVI values in different grazing areas of Inner Mongolia.

Grazing Plot and Intensity	NDVI (%)	Mean Value	Accuracy (%)
CK1	40.40	38.206	80.57
CK2	42.68		75.76
CK3	31.54		85.21
LGI1	41.26	37.923	85.17
LGI2	31.54		83.65
LGI3	40.97		73.46
MGI1	33.82	37.623	89.83
MGI2	42.93		84.06
MGI3	36.12		91.66
HGI1	39.44	37.703	97.09
HGI2	35.04		78.21
HGI3	38.63		91.05

, based on which the effect of different grazing intensities on the degree of desertification in Xilin Gol League can finally be obtained. In the calculations here, by default, grazing intensity only has an effect on the NDVI values, and the other values are kept as the average for the year (2022).

Based on the mean data in Table 7 and Equation (20), the results of the desertification index for different grazing intensities in the Xilin Gol League region over a certain period in the future are shown in

Table 8. The effect of different grazing intensities on the desertification index over a certain period in the future.

Grazing Intensity	Desertification Index Values
CK	$0.2471$
LGI	$0.2473$
MGI	$0.2475$
HGI	$0.2474$

.

Therefore, according to the principal component analysis, different weights can be obtained for the different components, and the gravel diagram of the variables related to soil slabbing in the Xilin Gol League region [43,44] is shown in Figure 18.

Table 9. Score coefficients for the components of the factors associated with soil slabbing in the Xilin Gol League region.

	Composition
	1	2
SOC	0.018	0.488
SIC	0.031	0.509
STC	0.051	0.729
TN	0.003	0.387
C/N	0.034	0.438
SC (10 cm)	0.255	0.056
SC (20 cm)	0.255	0.056
SC (30 cm)	0.255	0.056
Moisture (10 cm)	0.040	0.266
SC (50 cm)	0.255	0.056

illustrates the matrix of score coefficients for the components of the factors associated with soil slabbing in the Xilin Gol League region.

The resulting component weighting factors F are shown in

Table 10. Soil slabbing weighting factor in the Xilin Gol League region.

Chemical properties	0.634974146	SOC	0.119510496
		SIC	0.127538715
		STC	0.184268742
		TN	0.092134885
		C/N	0.111521307
Physical properties	0.365025854	SC (10 cm)	0.07319175
		SC (20 cm)	0.07319175
		SC (30 cm)	0.07319175
		Moisture (10 cm)	0.072258856
		SC (50 cm)	0.07319175

.

Regarding the evaluation of the soil slumping index, the comprehensive evaluation method in [45,46] was referred to and improved by adding weighting factors, resulting in a comprehensive evaluation matrix U, as shown in

Table 11. A comprehensive evaluation matrix on the soil slabbing index in the Xilin Gol League region.

SOC	SIC	STC	TN	C/N	SC (10 cm)	SC (20 cm)	SC (30 cm)	Moisture (10 cm)	SC (50 cm)
16.25	6.93	23.19	2.14	10.02	1.24	1.30	1.36	13.73	1.38
16.91	2.83	19.74	2.06	9.58	1.24	1.30	1.36	13	1.38
…	…	…	…	…	…	…	…	…	…
14.74	3.35	18.10	1.96	9.23	1.22	1.28	1.34	12.21	1.39
14.78	9.82	24.60	1.71	14.36	1.22	1.28	1.34	11.86	1.39

.

Through the comprehensive evaluation matrix U, the evaluation matrix R can be obtained in

Table 12. An evaluation matrix on the soil slabbing index in the Xilin Gol League region.

0.55212	0.484185	0.502826	0.750302	0.148471	1	1	1	0.401615	0
0.594266	0.184307	0.27705	0.702296	0.111994	1	1	1	0.327952	0
…	…	…	…	…	…	…	…	…
0.455517	0.222891	0.169788	0.645289	0.082929	0	0	0	0.248234	1
0.457607	0.69485	0.595075	0.501695	0.506654	0	0	0	0.212916	1

.

Let X denotes the influence of different grazing intensities on other factors affecting slabbing in the soil. Based on the data in Table 7, the X values corresponding to different grazing intensities (CK, LGI, MGI, and HGI) on soil slaking were $0.304081045$, $0.249280441$, $0.359852477$, and $0.286786036$, respectively, which led to the final quantitative analysis results of soil slaking $R_1$, where $R_1=X\times U\times F$, as shown in

Table 13. Quantitative results of soil slabbing analysis in the Xilin Gol League region.

No.	Soil Slumping Index
1	$0.199600716839267$
2	$0.167708241995365$
$...$	$...$
$n-1$	$0.0985937002601645$
n	$0.159870224975688$

. The higher the $R_1$, the less slabbed the soil.

3.3. A More Reasonable Grazing Strategy of the Xilin Gol League in the Future

The range of grazing intensity for typical areas of Eurasian steppe is now obtained as $[0,8]$ in standard sheep/day/ha, based on a review of [47,48,49]. Based on the basic premise that the best grazing period is from June to September in the Xilin Gol Union, the grazing time range for the grazing method was set at [0, 120] in days/year. The grazing method in this study was fixed to rotational grazing. The above predictions of soil moisture in the grasslands of the Xilin Gol League region in 2023 and 2024 and the results of the evaluation indicators of soil desertification and soil consolidation in the grasslands of the region in a certain period of time in the future were combined to obtain various key indicators of the grasslands of the Xilin Gol League region under different grazing strategies (all standardized, as shown in

Table 14. Vegetation cover in the Xilin Gol Meng region under different grazing strategies.

	Intensity	0	1	2	3	4	5	6	7	8
Time (Days)
1		1	1	1	1	1	1	1	1	1
2		1	0.99838	0.996762	0.995146	0.993531	0.991918	0.990306	0.988696	0.987088
3		1	0.996762	0.993531	0.990306	0.987088	$9.84\times {10}^{-1}$	$9.81\times {10}^{-1}$	$9.77\times {10}^{-1}$	$9.74\times {10}^{-1}$
…		…	…	…	…	…	…	…	…	…
121		1	$8.17\times {10}^{-1}$	$6.54\times {10}^{-1}$	$5.10\times {10}^{-1}$	$3.82\times {10}^{-1}$	$2.69\times {10}^{-1}$	0.1683	0.079107	0

,

Table 15. Organic carbon content in the Xilin Gol Meng region under different grazing strategies.

	Intensity	0	1	2	3	4	5	6	7	8
Time (Days)
1		0	0	0	0	0	0	0	0	0
2		0	$1.33\times {10}^{-10}$	$2.69\times {10}^{-10}$	$4.08\times {10}^{-10}$	$5.49\times {10}^{-10}$	$6.93\times {10}^{-10}$	$8.40\times {10}^{-10}$	$9.90\times {10}^{-10}$	$1.14\times {10}^{-9}$
3		0	$2.69\times {10}^{-10}$	$5.49\times {10}^{-10}$	$8.40\times {10}^{-10}$	$1.14\times {10}^{-9}$	$1.46\times {10}^{-9}$	$1.79\times {10}^{-9}$	$2.13\times {10}^{-9}$	$2.48\times {10}^{-9}$
…		…	…	…	…	…	…	…	…	…
121		0	$6.40\times {10}^{-8}$	$7.73\times {10}^{-7}$	$7.80\times {10}^{-6}$	$8.21\times {10}^{-5}$	$8.62\times {10}^{-4}$	0.009059	0.095179	1

,

Table 16. Nitrogen content in the Xilin Gol Meng region under different grazing strategies.

	Intensity	0	1	2	3	4	5	6	7	8
Time (Days)
1		0	0	0	0	0	0	0	0	0
2		0	$1.67\times {10}^{-4}$	$3.34\times {10}^{-4}$	$5.02\times {10}^{-4}$	$6.70\times {10}^{-4}$	$8.39\times {10}^{-4}$	$1.01\times {10}^{-3}$	$1.18\times {10}^{-3}$	$1.35\times {10}^{-3}$
3		0	$1.34\times {10}^{-4}$	$6.70\times {10}^{-4}$	$1.01\times {10}^{-3}$	$1.35\times {10}^{-3}$	$1.69\times {10}^{-3}$	$2.04\times {10}^{-3}$	$2.38\times {10}^{-3}$	$2.73\times {10}^{-3}$
…		…	…	…	…	…	…	…	…	…
121		0	$2.42\times {10}^{-2}$	$5.93\times {10}^{-2}$	$1.10\times {10}^{-1}$	$1.84\times {10}^{-1}$	$2.91\times {10}^{-1}$	0.447009	0.672661	1

,

Table 17. Carbon and nitrogen content in the Xilin Gol Meng region under different grazing strategies.

	Intensity	0	1	2	3	4	5	6	7	8
Time (Days)
1		1	1	1	1	1	1	1	1	1
2		1	$9.83\times {10}^{-1}$	$9.67\times {10}^{-1}$	$9.51\times {10}^{-1}$	$9.35\times {10}^{-1}$	$9.19\times {10}^{-1}$	$9.04\times {10}^{-1}$	$8.89\times {10}^{-1}$	$8.74\times {10}^{-1}$
3		1	$9.67\times {10}^{-1}$	$9.35\times {10}^{-1}$	$9.04\times {10}^{-1}$	$8.74\times {10}^{-1}$	$8.45\times {10}^{-1}$	$8.17\times {10}^{-1}$	$7.90\times {10}^{-1}$	$7.64\times {10}^{-1}$
…		…	…	…	…	…	…	…	…	…
121		1	$1.33\times {10}^{-1}$	$1.77\times {10}^{-2}$	$2.36\times {10}^{-3}$	$3.15\times {10}^{-4}$	$4.18\times {10}^{-5}$	$5.48\times {10}^{-6}$	$6.44\times {10}^{-7}$	0

and

Table 18. Soil moisture in the Xilin Gol Meng region under different grazing strategies.

	Intensity	0	1	2	3	4	5	6	7	8
Time (Days)
1		0	0.125	0.25	0.375	0.5	0.625	0.75	0.875	1
2		0	$1.25\times {10}^{-1}$	$2.50\times {10}^{-1}$	$3.75\times {10}^{-1}$	$5.00\times {10}^{-1}$	$6.25\times {10}^{-1}$	$7.49\times {10}^{-1}$	$8.74\times {10}^{-1}$	$9.99\times {10}^{-1}$
3		0	$1.25\times {10}^{-1}$	$2.50\times {10}^{-1}$	$3.74\times {10}^{-1}$	$4.99\times {10}^{-1}$	$6.24\times {10}^{-1}$	$7.49\times {10}^{-1}$	$8.74\times {10}^{-1}$	$9.99\times {10}^{-1}$
…		…	…	…	…	…	…	…	…	…
121		0	$1.14\times {10}^{-1}$	$2.29\times {10}^{-1}$	$3.43\times {10}^{-1}$	$4.58\times {10}^{-1}$	$5.72\times {10}^{-1}$	0.686393	0.800792	0.91519

; to save space, only some of them are given).

Based on the model of soil desertification and soil consolidation developed above, the changes in soil consolidation and desertification under different grazing strategies in the Xilin Gol League can be obtained as shown in Figure 19 and Figure 20.

The analysis shows that a smaller soil desertification index represents a lower degree of desertification; a higher soil slab index represents a lower degree of slabbing [50,51]. Based on the data obtained, it can be concluded that when the grazing intensity is within $[3,4]$ and the number of grazing days is $[0,25]$ or $[85,120]$, the degree of desiccation is low and stable throughout the range of values; when the grazing intensity is within $[2,3]$ and the number of grazing days is $[47,104]$, the degree of desiccation is low and the values are reasonable.

For the intersection of the two, combined with the influence of the above key indicators for the region, it can be concluded that when the grazing intensity is 3 and the number of grazing days is $[85,104]$, the degree of soil denudation and soil desertification in the rangelands of the Xilin Gol League is the lowest, and the other indicators are in the best range.

4. Conclusions

In this paper, based on machine learning methods and fuzzy evaluation matrices, preliminary predictions and analyses are made of the future health and sustainability of the Eurasian steppe region, represented by the Xilin Gol League. The main work of this paper is as follows: (i) Based on dynamic time series prediction and support vector regression, using three key indicator data for the grasslands in the Xilin Gol region from 2012 to 2022, a model was developed to predict soil moisture at different depths in 2023 and 2024, while keeping the current grazing strategy unchanged. (ii) Based on the kernel principal component analysis, a model of the degree of grassland desertification was developed with three dimensions. Subsequently, a quantitative definition of soil denudation was given by analyzing the main influencing factors of grassland soil degradation; meanwhile, a prediction model for soil denudation evaluation was established based on the fuzzy comprehensive evaluation matrix. (iii) Based on the above analysis, a comprehensive analysis of the grassland ecological data and herders’ livestock data in the Xilin Gol League region in the past ten years was conducted to give a more reasonable grazing strategy in the future. The research results of this paper have certain guiding significance and reference value for the future sustainable development of grassland ecology and livestock policy formulation in the Xilin Gol League grassland and similar areas.