Revealing the Eco-Environmental Quality of the Yellow River Basin: Trends and Drivers

Abstract

The Yellow River Basin (YB) acts as a key barrier to ecological security and is an important experimental region for high-quality development in China. There is a growing demand to assess the ecological status in order to promote the sustainable development of the YB. The eco-environmental quality (EEQ) of the YB was assessed at both the regional and provincial scales utilizing the remote sensing-based ecological index (RSEI) with Landsat images from 2000 to 2020. Then, the Theil–Sen (T-S) estimator and Mann–Kendall (M-K) test were utilized to evaluate its variation trend. Next, the optimal parameter-based geodetector (OPGD) model was used to examine the drivers influencing the EEQ in the YB. Finally, the geographically weighted regression (GWR) model was utilized to further explore the responses of the drivers to RSEI changes. The results suggest that (1) a lower RSEI value was found in the north, while a higher RSEI value was found in the south of the YB. Sichuan (SC) and Inner Mongolia (IM) had the highest and the lowest EEQ, respectively, among the YB provinces. (2) Throughout the research period, the EEQ of the YB improved, whereas it deteriorated in both Henan (HA) and Shandong (SD) provinces. (3) The soil-available water content (AWC), annual precipitation (PRE), and distance from impervious surfaces (IMD) were the main factors affecting the spatial differentiation of RSEI in the YB. (4) The influence of meteorological factors (PRE and TMP) on RSEI changes was greater than that of IMD, and the influence of IMD on RSEI changes showed a significant increasing trend. The research results provide valuable information for application in local ecological construction and regional development planning.

1. Introduction

Sustainable development will have an impact on future generations [1]. The Chinese government regards sustainable development as a major strategy to promote China’s national development. The Yellow River Basin (YB) is a crucial ecological and economic hub in China, contributing significantly to regional development and ecological civilization construction. The monitoring and quantitative evaluation of EEQ in the YB can provide a reference for environmental planning and ecological protection policy formulation and is of great importance in realizing regional sustainable development goals.

Remote sensing technology can efficiently and accurately objectively gather extensive ground feature information and has been widely applied in assessing EEQ [2,3]. The eco-environment is comprehensively affected by multiple factors, and the evaluation results of a single index can only represent changes in a certain aspect of the terrestrial ecological status [4,5]. As a result, scholars have started focusing on all-encompassing measures to better assess the regional EEQ comprehensively and effectively. The most commonly utilized EEQ models at present are the analytic hierarchy process (AHP) [6], the ecological index (EI) [7], the “pressure-state-response” conceptual framework (PSR) [8] and the remote sensing-based ecological index (RSEI) [7].

The AHP is a comprehensive decision-making method combining qualitative and quantitative methods, and since factor weights are relatively fixed in the analysis process, it is unable to deal with possible dynamic changes in indicators [9]. The EI was developed by China’s Environmental Protection Administration in 2006 and can be used for annual EEQ evaluation in areas above the county level [7]. During the process of EI construction, index weight selection is affected by subjectivity, thus affecting the accuracy of EEQ evaluation results. The PSR framework, proposed by the OECD, chose metrics from the three levels of pressure, state, and response, combining social, economic, and environmental elements [10]. According to the PSR framework, the RSEI couples greenness, humidity, heat, and dryness to represent the overall ecological status of the area [4]. The RSEI is based entirely on remote sensing image inversion and has strong comparability with the EI [11]. It has been widely applied to EEQ evaluation at various spatial scales in cities [12], mining areas [13], nature reserves [14], basins [15], and countries [16]. Compared with EI and AHP methods, RSEI is more conducive to the evaluation of the ecological status of uninterrupted land cover [11].

Ecological conditions are affected by a variety of factors, which can be divided into natural factors such as terrain, soil, and climate, as well as human factors including social economy, among others [17]. The analysis methods for driver identification mainly include correlation analysis, principal component analysis, linear regression analysis, geographical detectors, and spatial regression analysis [18,19,20,21]. Geographical detectors can detect both numerical and continuous data and can avoid the influence of multivariable collinearity [22]. The geographically weighted regression (GWR) model is capable of establishing the spatial heterogeneity of parameters across different regions. Both geographical detectors and GWR take into account the spatial effects of data and have been widely used in driver analysis [21,22,23].

The YB is a key belt connecting China’s eastern, central, and western regions. It is an important task in regional coordinated development in the YB to evaluate the EEQ and its changing patterns. In this study, we examined the spatial–temporal differences in the YB’s EEQ from 2000 to 2020 at the overall and provincial levels and identified the driving factors affecting EEQ. The following were the main focuses of this paper: (1) We performed quantitative analysis of EEQ based on the RSEI model. (2) The Theil–Sen (T-S) estimator and Mann–Kendall (M-K) method were used to analyze the spatiotemporal variation of RSEI. (3) The main factors affecting the spatial heterogeneity of RSEI were obtained by using the optimal parameter-based geographical detector (OPGD) model. (4) The GWR model was used to examine the responses of driving factors to RSEI changes.

2. Study Area and Datasets

2.1. YB Region

The YB in the northwest of China covers Qinghai (QH), Sichuan (SC), Gansu (GS), Ningxia (NX), Inner Mongolia (IM), Shaanxi (SN), Shanxi (SX), Henan (HA) and Shandong (SD) provinces (Figure 1) and lies between the latitudes 32°10′N and 41°50′N and the longitudes 95°53′E and 119°05E. The YB consists of four geomorphological units, with varying temperatures and precipitation across the landscape. The majority of areas in the upper and middle sections of the YB are situated in dry and semi-dry regions, with ecosystems that are delicate and easily affected. This is especially the case in central and the southern IM, northern SN, southeastern GS, and southern NX provinces [24].

2.2. Datasets

Landsat 5 and Landsat 8 tier-1 surface reflectance (SR) data obtained from the GEE platform were used, with a 30 m spatial resolution. Landsat images captured from June to September during the target years 2000, 2005, 2010, 2015, and 2020, and each year before and after, were selected to synthesize cloud-free mean data. DEM data, which were published by NASA and NIMA, had a spatial resolution of 30 m and were used to extract elevation and slope measurements in the study area (https://developers.google.com/earth-engine/datasets/catalog/USGS_SRTMGL1_003, accessed on 10 November 2023). Sand, clay, silt, soil organic carbon content, and soil available water content were obtained from the soilgrid250m dataset [25], available at https://soilgrids.org/, accessed on 26 May 2023. Soil data (clay, sand, silt, and soil organic carbon) were used to calculate the soil erodibility factor. More details are available in the Supplementary Materials. Precipitation and temperature datasets have a spatial resolution of 1000 m and are available from https://poles.tpdc.ac.cn/zh-hans/, accessed on 16 April 2023. The land cover data [26] with a 30 m spatial resolution were released by Huang Xin’s team at Wuhan University (http://irsip.whu.edu.cn/resources/CLCD.php, accessed on 7 April 2023). The GDP data with a 1000 m spatial resolution were obtained from https://www.resdc.cn/DOI/, accessed on 31 July 2023. The distance from impervious surfaces (IMD) was calculated using Euclidean distance, with impervious data obtained from land cover data. The geographic reference system for all datasets was unified into the WGS-84 coordinate system and UTM projection.

3. Methods

This study’s framework comprised three parts (Figure 2). The first part was the calculation of the RSEI. The second part evaluated the EEQ of the study area through the spatial distribution analysis and spatial–temporal change analysis of the RSEI. The third part employed the OPDG model and GWR model to explore the factors influencing the EEQ in the YB.

3.1. RSEI Model

Principal component analysis (PCA) was used to derive the RSEI model [7,11] from NDVI, WET, LST, and NDBSI. Before conducting PCA, a water body mask needed to be applied to each index using the modified normalized difference water index proposed by Xu (MNDWI) [27].

The following formula was used to calculate RSEI:

$$RSE{I}_0=1-PC1\left(NDVI,WET,LST,NDBSI\right)$$

(1)

$$RSEI=Nor\left(RSE{I}_0\right)$$

(2)

where Nor( ) represents the min-max normalization. Before conducting PCA, four components were normalized by means of min–max normalization methods.

3.2. T-S Estimator and M-K Test

The T-S estimator [28,29] and the M-K test [30,31] are commonly utilized in analyzing ecological factors to detect trend changes in time series data and minimize the influence of data inaccuracies.

The T-S estimator calculation formula is as follows:

$$\beta =median\left(\frac{RSE{I}_j-RSE{I}_i}{j-i}\right),\forall j>i$$

(3)

where RSEI_i and RSEI_j denote the RSEI values for year i and year j, respectively. When β > 0, RSEI displays a growing trend; otherwise, it displays a decreasing trend.

The M-K test calculation formula is as follows:

$$Z=\left\{\begin{array}{c}(S-1)/\sqrt{\mathrm{var}(S)},S>0\\ 0,S=0\\ (S+1)/\sqrt{\mathrm{var}(S)},S<0\end{array}\right.$$

(4)

$$S={\displaystyle {\sum }_{i=1}^{n-1}{\displaystyle {\sum }_{j=i+1}^n\mathrm{sgn}\left(RSE{I}_j-RSE{I}_i\right)}}$$

(5)

$$\mathrm{sgn}\left(RSE{I}_j-RSE{I}_i\right)=\left\{\begin{array}{c}1,RSE{I}_j-RSE{I}_i>0\\ 0,RSE{I}_j-RSE{I}_i=0\\ -1,RSE{I}_j-RSE{I}_i<0\end{array}\right.$$

(6)

$$\mathrm{var}(S)=\mathrm{n}\left(n-1,2n+5\right)/18$$

(7)

where n is the length of the time series, sgn( ) is the symbolic function, and Z is the significance test statistic with values in the range (−∞, +∞).

The trend categories of the T-S estimator and M-K test are shown in

Table 1. T-S estimator and M-K test trend categories.

β	Z	Trend Type	Trend Features
<0	≤−1.96	SD	Significant decrease
<0	−1.96–1.96	NSD	No significant decrease
0	−1.96–1.96	NC	No change
>0	−1.96–1.96	NSI	No significant increase
>0	≥1.96	SI	Significant increase

NOTE: β is the variation trend of the RSEI sequence calculated by the T-S estimator; Z is the standardized test statistics calculated by the M-K test.

.

3.3. OPGD Model

Based on the OPGD model [32], the best discretization scheme for the continuous variables was selected and geographical exploration analysis was subsequently conducted. The classification methods used were based on equal intervals, Jenks natural breaks (Jenks), and quantiles, and the classification classes were 4–10. The factor detector and interaction detector were used for the experiment (R version 4.1.3). Factor detection can identify the q value, representing the impact of the drivers on the RSEI. The q can have a value between 0 and 1. Greater q values indicate increased explanatory strength. Interaction detection can detect the influence of the interaction between two drivers on the explanatory power of the RSEI.

In this paper, elevation (ELE), slope (SLO), soil erodibility (SE), soil available water content (AWC), annual precipitation (PRE), annual mean temperature (TMP), land cover types (LCT), GDP and distance from impervious surfaces (IMD) were selected. The spatial distribution of driving factors is shown in Figure 3. These nine factors were divided into four categories based on their content and attributes: (1) topographic factors (ELE and SLO); (2) soil factors (SE and AWC); (3) meteorological factors (PRE and TMP); and (4) human factors (LCT, GDP, and IMD). Considering the area of YB and satisfying the requirements of spatial heterogeneity, a 10 km × 10 km grid was generated. LCT data were obtained for the types that occupy the largest proportion in each grid, and other factors were calculated for the grid mean data.

3.4. Geographically Weighted Regression

The geographically weighted regression (GWR) model is a widely used local spatial regression model, which can effectively reflect the spatial heterogeneity of regression relationships [33].

The general expression of GWR is as follows:

$${\mathrm{y}}_i={\beta }_0(u_i,v_i)+{\displaystyle \sum {\beta }_k(u_i,v_i)}{x}_{ik}+{\epsilon }_i$$

(8)

where y_i is the RSEI fitting value of sample i, β_k is the band estimation coefficient of sample i, (u_i,v_i) is the central coordinate of the spatial unit of sample i, β₀(u_i, v_i) is the estimated value of the constant term of sample i, and ε_i is the independent identically distributed error term.

The GWR model was used to analyze the influence of driving factors on RSEI changes. Since topographic and soil factors can be regarded as constants during the study period, meteorological factors and human factors are selected as regression variables. Excluding numerical data (LCT) and local collinearity data (GDP), the final factors are selected as PRE, TMP, and IMD. The fitting results of the GWR model are shown in

Table 2. Comparison of the results of OLS and GWR models.

Period	Model	AICc	Sigma	R2	Adjust R2
2000–2005	OLS	24,085.113	1312.632	0.047	0.047
	GWR	10,723.586	8225.790	0.848	0.826
2005–2010	OLS	24,076.407	8217.495	0.048	0.047
	GWR	11,776.483	1469.866	0.830	0.830
2010–2015	OLS	21,448.259	6059.903	0.298	0.297
	GWR	9117.758	1099.883	0.873	0.855
2015–2020	OLS	23,851.793	8006.352	0.072	0.072
	GWR	9843.700	1139.525	0.868	0.846
2000–2020	OLS	23,281.808	7494.584	0.131	0.131
	GWR	10,390.218	1229.309	0.858	0.835

.

4. Results

4.1. Spatial Distribution and Change Characteristics of EEQ

The RSEI was categorized into five classes via the Jenks method, namely Poor, Fair, Moderate, Good, and Excellent. The RSEI classes were lower in the north and higher in the south (Figure 4). The Poor and Fair entries were concentrated in the upper and middle reaches of the YB, including Lanzhou and Baiyin in GS, IM, and the arid zone of central NX. The Moderate and Good entries were primarily found in the Yellow River’s origin areas and the YB’s lower reaches. The Excellent entries were primarily found in the Ruoergai-Maqu Ecological Function Reserve, the Qilian Mountains, and the Qinling Mountains.

At both the provincial and the land cover levels, the RSEI classes in the YB had regional characteristics (Figure 5a). The ecological condition of NX and IM was not good, with more than 60% of this area being categorized as Poor or Fair. The ecological condition of other provinces was higher, with the majority of regions categorized as Moderate and above (an area exceeding 70% of the whole). SC had the highest categorized ecological condition, with more than 90% of its total area classified as Good or Excellent. The distribution of RSEI classes varied across different land cover types (Figure 5b). The distribution of RSEI classes in forestland, shrubland, and wetland was mainly Good and Excellent, while barren land was almost entirely Poor. The ecological condition of cropland, grassland, and impervious land was dominated by multiple RSEI classes.

The T-S estimator and M-K test were used to perform a statistical analysis of EEQ changes in the YB during the period of 2000–2020 (Figure 6,

Table 3. Grading table of RSEI change amplitude in YB.

Class	Trend Type	Area/km2		Proportion/%
Degeneration	SD	21,375	34,371	2.63	4.27
	NSD	13,356		1.64
No change	NC	599,800	599,800	73.87	73.87
Improvement	NSI	43,267	177,488	5.33	21.86
	SI	134,221		16.53

). A total of 21.86% of EEQ improvement areas were mainly distributed in northern QH, SC, eastern GS, southern NX, northeastern IM, western SN, and SX. Of these, 16.53% were significant improvement areas, mainly located in ecological function protection areas on the Loess Plateau. The portion of land that was degraded was 4.27%, while the portion that was significantly degraded was 2.63%. Most of the deteriorated regions were located in the Ningxia Plain urban agglomeration, Hetao Plain in IM, the Guanzhong Plain urban agglomeration (centered on Xi’an in SN), the Jinzhong urban agglomeration in SX, HA, and SD.

Table 4. The grading scale of EEQ change in YB provinces (%).

Trend Type	QH	SC	GS	NX	IM	SN	SX	HA	SD
SD	1.61	1.81	1.56	2.99	0.93	3.74	2.31	10.91	9.96
NSD	0.74	1.15	0.97	1.41	1.08	2.07	1.19	6.85	8.31
NC	81.76	86.50	75.71	79.49	79.77	62.75	60.58	72.91	70.42
NSI	7.37	4.83	5.15	3.78	6.00	4.86	3.39	4.28	5.08
SI	8.53	5.73	16.61	12.33	12.21	26.59	32.54	5.06	6.22

shows the RSEI trend change types in the provinces of the YB from 2000 to 2020. Throughout this timeframe, aside from HA and SD, the ratio of land improvement areas in the other seven provinces exceeded that of the land degradation area. The change ratio of land ecological status in the nine provinces in descending order was SX > SN > SD > HA > GS > NX > IM > QH > SC. SX and SN showed the most obvious improvement in EEQ, with area improvement rates reaching 35.93% and 31.45%, respectively—far greater than the land degradation area. The degradation of EEQ was most obvious in HA and SD, with area degradation ratios of 17.76% and 18.27%, respectively.

Table 5. The grading scale of EEQ change in different land cover types (%).

Trend Type	FO	SH	WE	CR	GR	IM	BA	CR-GR
SD	0.60	1.03	1.27	5.68	1.01	13.16	0.01	5.50
NSD	0.20	0.51	0.62	4.44	0.60	4.59	0.17	2.53
NC	93.06	81.62	76.36	61.27	78.02	70.28	98.10	62.43
NSI	1.60	4.33	6.03	7.03	5.73	3.13	0.70	5.62
SI	4.53	12.51	15.71	21.58	14.63	8.84	1.01	23.92
Trend type	BA-GR	GR-CR	CR-IM	GR-FO	GR-BA	CR-FO	GR-IM	BA-CR
SD	0.01	0.43	42.29	0.07	1.07	0.26	10.02	0.03
NSD	0.03	0.45	12.37	0.08	0.67	0.19	3.43	0.18
NC	90.99	49.23	40.42	41.48	96.01	55.93	77.40	48.15
NSI	0.53	7.04	2.23	3.83	0.81	4.96	2.69	3.05
SI	8.44	42.85	2.70	54.54	1.43	38.66	6.45	48.59

shows the RSEI trend change types in different land cover change types of the YB from 2000 to 2020. The seven types of land cover in the area without conversion and the top nine types of land cover change were selected. In areas where the land cover has not been changed, the EEQ of cropland and impervious land has changed greatly, with the proportion of change being more than 38%. Only the EEQ of impervious areas has been degraded, and the EEQ of the other land cover has been improved. Changes in land cover change are often accompanied by changes in ecological conditions. Except for BA-GR and GR-BA, the ecological area change rates of other transfer types ranged from 22.6% to 59.58%. In general, the RSEI index can reflect the difference in EEQ among different land cover types and the change in EEQ brought about by type conversion.

4.2. Analysis of Driving Factors of RSEI Spatial Differentiation

The explanatory power of driving factors was as follows: AWC > PRE > IMD > LCT > SLO > TMP > ELE > GDP > SE (

Table 6. q statistics of detection factors (95% confidence level).

Year	ELE	SLO	SE	AWC	PRE	TMP	LCT	GDP	IMD
2000	0.289	0.280	0.150	0.818	0.663	0.328	0.318	0.186	0.453
2005	0.229	0.296	0.121	0.791	0.615	0.303	0.344	0.158	0.430
2010	0.243	0.313	0.116	0.797	0.633	0.305	0.332	0.132	0.417
2015	0.176	0.332	0.084	0.725	0.602	0.247	0.350	0.110	0.357
2020	0.229	0.348	0.100	0.803	0.555	0.303	0.306	0.160	0.447

). AWC, PRE, and IMD were the dominant factors of RSEI, with 5-year mean q values of 0.787, 0.614, and 0.421, respectively. LCT, SLO, TMP, and ELE also have an important impact on EEQ, with five-year average q values of 0.330, 0.314, 0.297, and 0.233, respectively. The explanatory power of GDP and SE is relatively low.

The interactive detection analysis results clearly show that the RSEI in the YB is influenced by both bilinear enhancement and nonlinear enhancement (Figure 7). This suggests that the combined effect of these two factors is more significant than that of a signal factor. The explanatory power of AWC∩ELE and AWC∩LCT for the interaction of factors is significantly stronger than that of other factors. Although GDP and SE may not individually have a strong impact, their combined influence on the spatial differentiation features of the RSEI is evident in interactive detection.

Figure 8 shows the factor detection results of provinces in the YB. In QH, SC, and GS, natural factors are the main influencing factors. In NX, IM, SN, SX, HA, and SD, both natural and human factors show a strong influence. The interactive detection results show that the factors selected in this study are reasonable, and the interaction between factors can explain the spatial differentiation characteristics of RSEI well (

Table 7. The main interaction factors of nine provinces in YB (95% confidence level).

Name	2000	2005	2010	2015	2020
QH	AWC∩TMP (0.898)	AWC∩TMP (0.898)	AWC∩TMP (0.874)	ELE∩AWC (0.867)	AWC∩TMP (0.884)
	AWC∩IMD (0.858)	AWC∩IMD (0.878)	AWC∩IMD (0.870)	AWC∩TMP (0.870)	AWC∩IMD (0.875)
SC	ELE∩AWC (0.703)	ELE∩SE (0.624)	SE∩TMP (0.642)	ELE∩PRE (0.600)	ELE∩PRE (0.624)
	PRE∩IMD (0.782)	PRE∩IMD (0.676)	AWC∩TMP (0.614)	AWC∩PRE (0.595)	AWC∩PRE (0.647)
GS	ELE∩AWC (0.915)	ELE∩AWC (0.919)	ELE∩AWC (0.922)	ELE∩AWC (0.908)	ELE∩AWC (0.925)
	AWC∩LCT (0.929)	TMP∩AWC (0.916)	TMP∩AWC (0.917)	TMP∩AWC (0.899)	TMP∩AWC (0.920)
NX	ELE∩AWC (0.884)	ELE∩AWC (0.884)	ELE∩AWC (0.880)	SLO∩AWC (0.878)	SLO∩AWC (0.895)
	SLO∩AWC (0.880)	SLO∩AWC (0.879)	SLO∩AWC (0.889)	AWC∩LCT (0.882)	AWC∩LCT (0.896)
IM	ELE∩AWC (0.817)	AWC∩PRE (0.839)	AWC∩PRE (0.808)	AWC∩PRE (0.798)	AWC∩PRE (0.813)
	AWC∩LCT (0.824)	AWC∩GDP (0.843)	AWC∩GDP (0.817)	AWC∩GDP (0.823)	AWC∩GDP (0.828)
SN	PRE∩AWC (0.932)	SLO∩AWC (0.927)	ELE∩AWC (0.903)	ELE∩AWC (0.910)	ELE∩AWC (0.946)
	AWC∩IMD (0.931)	AWC∩IMD (0.926)	TMP∩AWC (0.902)	TMP∩AWC (0.909)	SLO∩AWC (0.948)
SX	ELE∩AWC (0.840)	ELE∩AWC (0.862)	TMP∩AWC (0.854)	ELE∩AWC (0.861)	ELE∩AWC (0.855)
	TMP∩AWC (0.837)	TMP∩AWC (0.859)	AWC∩IMD (0.856)	TMP∩AWC (0.867)	TMP∩AWC (0.855)
HA	SLO∩PRE (0.806)	SLO∩PRE (0.846)	ELE∩SE (0.873)	ELE∩SE (0.845)	SLO∩TMP (0.902)
	PRE∩IMD (0.784)	PRE∩IMD (0.829)	SE∩IMD (0.876)	SLO∩GDP (0.841)	SE∩IMD (0.890)
SD	ELE∩AWC (0.822)	ELE∩AWC (0.794)	ELE∩AWC (0.853)	ELE∩IMD (0.815)	ELE∩IMD (0.840)
	ELE∩SE (0.796)	ELE∩SE (0.770)	ELE∩SE (0.830)	ELE∩SE (0.855)	ELE∩SE (0.786)

).

In conclusion, topographic factors, soil factors, meteorological factors, and human factors jointly influence the formation of spatial patterns of regional RSEI. Natural factors can explain most of the spatial differentiation characteristics of RSEI. The influence of human factors on RSEI can be strongly demonstrated when they interact with other factors.

4.3. Response of the Driving Factors to EEQ Changes

To explore the influence of changes in driving factors on the change in RSEI, the GWR model was performed with the grid difference of PRE, TMP, and IMD as the independent variable and the grid difference of RSEI as the dependent variable (

Table 8. Statistical analysis of regression coefficients of the GWR model.

Factor	2000–2005		2005–2010		2010–2015		2015–2020		2000–2020
	N+	N−	N+	N−	N+	N−	N+	N−	N+	N−
PRE	1892	1711	1723	1652	1054	868	2063	1621	2703	1867
TMP	1448	1428	1819	1652	2189	1905	927	828	792	749
IMD	1485	665	970	595	1906	707	1929	1261	1854	664

NOTE: N+: the number of grids of positive influence; N−: the number of grids of negative influence.

). By comparing the absolute value of the regression coefficient, the dominant driving factor of each grid was obtained, and the results are shown in Figure 9 and Figure 10.

According to the number of positive and negative regression coefficients, PRE and TMP have a slight positive influence on RSEI change in the YB, while IMD has a significant positive influence on RSEI change. From left to right, the dominant role of meteorological factors (PRE and TMP) on RSEI changes gradually decreases, while the dominant role of human factors on RSEI changes gradually increases. Meteorological factors almost dominate the RSEI changes in QH and SC, while IMD almost dominates the RSEI change in SD, with the proportion of dominant grids exceeding 90% from 2000 to 2020. From 2000 to 2005, 2005 to 2010, 2010 to 2015, and 2015 to 2020, the leading role of driving factors in different provinces of the YB was judged by summing the number of grids dominated by driving factors. The provinces dominated by PRE are QH and SC, the provinces dominated by TMP are GS, NX, IM, and SN, and the provinces dominated by IMD are HA and SD. It is worth noting that the number of grids dominated by IMD showed an increasing trend during the study period, indicating that the impact of human activities on EEQ change is becoming more and more obvious.

5. Discussion

5.1. Computational and Analytical Problems of RSEI

Landsat SR products on GEE contain pixel quality information, atmospheric opacity information, and radiometric saturation information. This information can be used to synthesize the best available data for target years in the YB. However, there were two anomalies in Landsat5 SR data: Christmas tree anomalies (Figure 11a) and caterpillar tracks (Figure 11b). For Christmas tree anomalies, we can use the establishment of an internal buffer zone of 3000 m to remove them [34]. Caterpillar tracks are found throughout the image, and the larger outliers in the ST_B6 band can be used as screening criteria. To save computation, we set up an inner ring buffer of 500 m for the image set with the Christmas tree anomaly removed and obtain the maximum value of ST_B6 band in this buffer. If the maximum value is greater than 60,000, we remove the image. The annual image can then be generated by combining the mean values.

Before RSEI calculation, annual image synthesis should be performed. Images with similar dates in different years are selected to calculate RSEI values [12,35]. This mostly occurs in small research areas; that is, one remote sensing image can completely cover the research area. The development of cloud platforms has also made the computation and processing of a large number of images very straightforward. In larger research areas, images covering two to four months with a high NDVI in the research year are selected to synthesize the images [36,37]. The YB was selected for the research, encompassing over 800,000 km². To balance the phenological differences in different periods and realize the coverage of the whole basin, we chose to synthesize the mean image using three years of Landsat data, from June to September.

Some scholars replace component indicators with alternative indicators to construct new RSEI models to monitor EEQ in specific areas (arid deserts, alpine regions, and so on). In the case of the YB, however, its vast area and diverse geological and geomorphological types make it difficult to find suitable component indexes. Moreover, creating new indicators through substituting component indicators [38] appears to introduce more subjectivity and challenges in the consistent application of the model. The effectiveness of the initial RSEI model has been demonstrated in different settings [39], leading to its application in monitoring EEQ in the YB. The average correlation coefficients of RSEI with NDVI, WET, LST, and NDBSI are 0.9758, 0.8726, −0.7964, and −0.9550, as shown in Figure 12, further validating the comprehensive representativeness of RSEI.

This study found that the EEQ of the YB showed an improving trend from 2000 to 2020, which is consistent with the research of Yang et al. [40]. The significant improvement area is mainly located in the Loess Plateau protection area. Zhang et al. [41] analyzed the eco-environment of the Loess Plateau and found that 84.51% of the EEQ in the study area had improved. The variations in RSEI are affected by many factors. The selection of factors in this paper is not sufficiently comprehensive. Natural factors such as wind speed, sunshine hours, and environmental protection policy-related factors were not selected. Future work will involve selecting a wider range of influencing factors, taking into account the current state of the research field.

5.2. Relationship between RSEI and Land Cover Change

Land cover is an important factor affecting the spatial differentiation of RSEI (Section 4.2). The distribution characteristics of RSEI differ significantly among different land cover types, and EEQ improvement or degradation areas tend to occur in areas where land type conversion occurs (Section 4.1). A comparative experiment was conducted to further explore the relationship between land cover types and RSEI. Figure 13 shows that the order of EEQ of different land cover types from high to low was forest, shrub, wetland, cropland, grassland, impervious land, and barren land.

Area changes of different land cover types in the YB for the period 2000–2020 were calculated (Figure 14). During this period, the area of cropland and grassland in QH decreased, while the area of barren land increased by 1177.93 km², and the grassland was severely degraded [42]. Problems with reduced grassland areas are also present in SC and GS. According to relevant studies, increased grassland degradation and areas of desertification have been caused by overgrazing in QH, SC, and GS [43], with the problem of grassland degradation remaining prominent. China’s prevention and control efforts regarding desertification have achieved a double reduction, in both desertification land and desertified land area. Indeed, the bare land areas of NX and IM reduced by 1141.14 km² and 6653.82 km², respectively, during the study period. A series of soil and water conservation measures, such as the comprehensive watershed management, the construction of the Three-North Shelterbelt, and the return of cropland to forest (grass), have partially improved the ecological status of the region. In addition, the ecological degradation of the YB is also related to the large-scale expansion of impervious surface, which has increased by 10,844.71 km² in the past 20 years.

In recent decades, China has attached great importance to the sustainable development of the YB and has launched and implemented a large number of major ecological protection, construction, and restoration projects, which are the main reasons for the improvement of the EEQ. Certain unresolved issues remain, however, such as the vulnerability and systematic degradation of the YB. The ecological environmental management of the YB is a long-term project, which heightens the requirements for ecological environmental quality evaluation. More accurate data support and enhanced analysis, potentially combined with expanded aspects of evaluation such as ecosystem services, are future research directions.

6. Conclusions

This research examines the patterns and drivers of RSEI in the YB region from 2000 to 2020. The distribution of EEQ in the YB area has regional characteristics, and the overall distribution was low in the north and high in the south. SC and IM had the highest and lowest EEQ, respectively. Throughout the research period, the improved area made up 21.86%, while the degraded area made up 4.27%. In the YB provinces, the degraded area of HA and SD is greater than the improved area. Soil-available water content (AWC), annual precipitation (PRE), and distance from impervious surfaces (IMD) were the main factors affecting the spatial distribution of RSEI. Precipitation, temperature, and IMD have important effects on RSEI variation in the YB, whereas the area dominated by precipitation and temperature fluctuated continuously during the study period, but the area dominated by IMD showed a significant increasing trend. The dominant factors vary among provinces. From west to east, IMD has an increasing influence on EEQ. Future research will focus on analyzing the changes and factors influencing EEQ in the YB at an annual scale with a more comprehensive approach, aiming to offer valuable insights for decision making related to ecological development and regional growth.