The Yellow River Basin (32°10′-41°50′N, 95°53′-119°05′E; Fig. 1) covers an area greater than 7.5×10⁵ km² (Lin et al., 2023). The basin is composed of four geomorphic units: the Inner Mongolia Plateau, the Loess Plateau, the Qinghai-Xizang Plateau, and the Huanghuai Plain. The terrain decreases from the west to the east, showing three major ladders, with a relief of more than 4000 m (He et al., 2023a). The climate is characterized by large temperature differences, unevenly distributed precipitation, large interannual precipitation changes, low air humidity, and significant evaporation. Vegetation is mostly distributed in the middle and lower reaches of the basin and the main vegetation types comprise forest, grassland, and farmland.

The data used in the present study were obtained from the National Aeronautics and Space Administration (NASA) of the United States (https://www.nasa.gov/). The data sources were Moderate Resolution Imaging Spectroradiometer (MODIS) product datasets (MOD09A1, MOD11A2, MOD13Q1, MOD15A2H, MOD17A3, and MCD43A3) from June to October of 2000-2022. The above datasets were preprocessed using MODIS Reprojection Tool (MRT), including stitching, projection (Social Impact Network (SIN) projection to Krasovsky- 1940-Albers), and resampling (all the above MODIS products were resampled into grids with the spatial resolution of 1000 m utilizing the bilinear interpolation method).

Nighttime light data from 2000 to 2022, corrected using Defense Meteorological Satellite Program/Operational Linescan System (DMSP/OLI) and National Polar-orbiting Operational Environmental Satellite System Preparatory Project/Visible Infrared Imaging Radiometer Suite (NPP/VIIRS) sensors, were obtained from Harvard Dataverse (Wei et al., 2023; Xu et al., 2023). As mentioned, the lighting data were obtained from Harvard Dataverse (https://dataverse.harvard.edu/dataset.xhtml?persistentId=doi:10.7910/DVN/GIYGJU), whereby the data from the two nighttime light sensors, NPP and DMSP, were fitted. The data of each period represent the average value for that year.

In the following study, nine driving factors were selected, namely elevation, slope, soil type (black soil, loess soil, red soil, sandy soil, and clay soil), temperature, precipitation, vegetation coverage, land use type, gross domestic product (GDP) density, and population density (Table 1). These driving factors were then divided into natural environmental factors and socio-economic factors. The selection of these factors involves a comprehensive consideration of natural and economic factors; first, they can cover the impact of the Yellow River Basin's ecological fragility to the greatest degree possible, and second, they are not only more conventional but also easier to obtain.

Natural environmental factors: digital elevation model (DEM) data with a spatial resolution of 500 m were obtained from General Bathymetric Chart of the Oceans (GEBCO). Slope was calculated in ArcGIS v.10.7 software (Esri, Redlands, California, USA) using DEM data. The soil type data were retrieved from Environmental and Scientific Data Center of the Chinese Academy of Sciences (https://www.resdc.cn). The temperature and precipitation data were downloaded from China Meteorological Data Service Center (https://data.cma.cn). Utilizing Kriging interpolation method and national meteorological station data, the grids of temperature and precipitation were obtained using ArcGIS v.10.7. The land use type and vegetation coverage data were obtained from MOD13Q1 and MCD12Q1 product data from NASA.

Socio-economic factors: GDP density data were derived from Figshare (https://doi.org/10.6084/m9.figshare.17004523.v1), which were calculated using the nighttime light data. Population density data were downloaded from Oak Ridge National Laboratory (ORNL) archive (https://landscan.ornl.gov).

PCA was used to transform the original multiple indicators into a few or several comprehensive indicators. The comprehensive indicators could be regarded as a linear combination of the initial indicators; therefore, these comprehensive indicators had no relationship with each other and could reflect the main information of the initial indicators (Ankush et al., 2023; Krzyśko et al., 2024).

The gravity center refers to the action point of the resultant force generated by the gravity of each part of an object. The gravity center could indicate the change deviation of ecological fragility. The transfer trajectory of gravity center could reflect the uneven variation in ecological vulnerability. The equations used in the present study are as follows (Wang et al., 2023b):

where Z_i is the attribute value of ecological vulnerability in the i^th year; x_j and y_j are the latitude and longitude coordinates of the j^th point, respectively; and

\overline{x}

and

\overline{y}

are the latitude and longitude coordinates of the gravity center of ecological vulnerability after the calculation of the n points, respectively.

The gravity center migration distance and direction could be utilized to describe the degree of distribution imbalance (Xia et al., 2023). The equations of direction and distance are as follows:

where θ is the gravity center deviation direction (°); d_m is the gravity center deviation distance (m); y_t+m and y_t are the latitude coordinates in time t+m and time t, respectively; and x_t+m and x_t are the longitude coordinates in time t+m and time t, respectively.

In the following study, the mechanism driving ecological vulnerability changes was analyzed using geodetector. Based on geospatial differentiation concept, researchers used geodetector to detect the spatial distribution consistency of dependent variables and independent variables (Shan and Wang, 2023). In this study, the dependent variable Y represents the ecological vulnerability value for each grid, and the independent variable X represents the value of nine indicators such as DEM, slope, temperature, etc. The formulas used are as follows:

where q is the explanatory power of a factor, which is typically applied to measure the factor's contribution to the target variable; h is the layer or partition of the dependent variable Y and the independent variable X; N_h is the number of units in the layer h; N is the number of units;

{\sigma }_h^2

is the variance of the Y value in the layer h; σ² is the variance in the Y value; W_ss is the total variance within each layer; and T_ss is the total variance of the entire area. For q∈[0,1], the greater the q value, the stronger the explanatory power of X to Y and vice versa.

Based on the above six indicators of humidity, heat, greenness, dryness, land degradation, and social economy, we applied PCA method to obtain the final remote sensing ecological vulnerability index for each year from 2000 to 2022. The PCA method can be used to not only reflect the specific ecological vulnerability index of each year but also better reflect the dynamic changes in the impact of different driving factors on ecological vulnerability.

where RSEVI is the remote sensing ecological vulnerability index; a_x is the contribution rate of the x^th principal component; and P_x is the x^th principal component.

In order to verify the accuracy of remote sensing ecological vulnerability index obtained via inversion, we selected 276 sample points using uniform point distribution method in 2022 (Table 2). We divided the ecological vulnerability into five levels based on the natural break method of ArcGIS v.10.7 to more intuitively reveal the spatial and temporal distribution characteristics of different levels of ecological vulnerability (Hou et al., 2020). Remote sensing ecological vulnerability index <0.80 indicated slight vulnerability, 0.80-0.95 indicated mild vulnerability, 0.95-1.10 indicated moderate vulnerability, 1.10-1.25 indicated intensive vulnerability, and >1.25 indicated severe vulnerability. The error matrix was constructed and calculated using the sample points provided by Google Earth (https://www.google.cn/intl/zh-en/earth/) and the evaluation results. The overall accuracy of remote sensing ecological vulnerability index was 86.36%. The accuracy of moderate vulnerability was the highest, with an accuracy rate of 91.53% being recorded, followed by intensive vulnerability and slight vulnerability, with accuracy rate of 90.48% and 90.38%, respectively, while the accuracy rate of mild vulnerability was the worst, at 84.78%. In general, the remote sensing ecological vulnerability index has high applicability.

Based on the series MODIS products and other datasets during 2000-2022, we calculated the remote sensing ecological vulnerability index for each year, and then utilized the raster calculator of ArcGIS v.10.7 to obtain the average remote sensing ecological vulnerability index for the Yellow River Basin (Fig. 2).

From 2000 to 2022, the average remote sensing ecological vulnerability index of the Yellow River Basin was 1.03, representing moderate vulnerability level. The spatial distribution patterns of ecological vulnerability at different levels differed greatly (Fig. 3). The area of slight vulnerability covered roughly 9.72×10⁴ km², accounting for 12.40% of the total area. This particular area was mostly located in the eastern part of Qinghai Province, the northern part of Sichuan Province, and the northwestern area of Henan Province, spanning the east-west direction. The mild vulnerability area was mostly located in the western part of Shandong Province and the central and eastern parts of Qinghai Province, with an area of roughly 1.65×10⁵ km², accounting for 20.98% of the total area. The moderate vulnerability covered an area of 1.66×10⁵ km², accounting for 21.18% of the total area, and was mainly distributed in the southern part of Gansu Province and the central part of Shaanxi Province. The intensive vulnerability area was the largest (1.99×10⁴ km²), accounting for 25.39% of the total area and being mainly distributed in the northern part of Shaanxi Province and the eastern part of Shanxi Province. The severe vulnerability covered an area of 1.57×10⁴ km², accounting for 20.05% of the total area, mainly located in the northern Shaanxi Province and the central Inner Mongolia Autonomous Region.

In this study, we used ArcGIS v.10.7 to determine the change rate of ecological vulnerability from 2000 to 2022, which could better display the change pattern and law of ecological vulnerability. The natural discontinuity method can reasonably divide the ecological vulnerability index into different intervals. In this study, the difference between the discontinuous values of different year intervals was not large, the error was less than 0.02, which was within a reasonable controllable range. The change rate of ecological vulnerability (CR) was divided into five levels: CR< -0.20 as severe improvement, -0.20<CR< -0.05 as moderate improvement, -0.05<CR<0.05 as stable, 0.05<CR<0.20 as moderate intensification, and CR>0.20 as severe intensification (Fig. 4).

From 2000 to 2022, the stable area constituted the largest area, accounting for 58.43% of the entire study area, and was mainly distributed in the central and southern parts of Gansu Province, the southern part of Inner Mongolia Autonomous Region, and the northern and southern parts of Shaanxi Province. Severe intensification areas were mainly distributed in the southeastern part of Qinghai Province and the northern of Sichuan Province, primarily due to their geographic location and unique climate types, which leaded to frequent natural disasters such as droughts, wind and sand, landslides, and debris flows, significantly impacting the ecological environment. Moderate intensification areas were mainly found in the western part of Qinghai Province, the central of Ningxia Hui Autonomous Region, and the western part of Shaanxi Province, including western of Tongren City and northwestern of Guyuan City. Moderate improvement areas were primarily in the northern part of Shaanxi Province and the western part of Shanxi Province, such as eastern of Yan'an City, Weinan City, eastern of Yulin City, and Linfen City. Severe improvement areas were concentrated in the eastern part of Lishi City and the western part of Yuncheng City in Shanxi Province. In summary, the overall ecological vulnerability showed a stable trend over the past 23 a, and the vulnerability of the central and eastern regions showed an obvious improvement trend.

The gravity center migration trajectory of ecological vulnerability can better reflect the spatial imbalance of ecological vulnerability in the study area. Using ArcGIS v.10.7, we calculated and analyzed the gravity center of ecological vulnerability at a 5-a scale (2000-2004, 2005-2009, 2010-2014, 2015-2019, and 2020-2022), and the trajectory change of the gravity center was analyzed.

The center of gravity of ecological vulnerability was very concentrated, with it being mainly distributed in Huanxian County, Henan Province, which is located in the central of the basin (Fig. 5). The change degree of the gravity centers in the northeastern and southwestern orientations was much greater than that in the other directions. Compared with that of 2000-2004, the average gravity center of 2005-2009 moved northwestward, indicating that the intensification degree in the northwest parts was greater than that of the southeast parts. The average gravity center of 2010-2014 moved southwest compared with that of 2005-2009, indicating that the intensification degree in the southwest parts was greater than that of the northeastern parts. Moreover, the average gravity center of 2015-2019 continued to move southwest, whereas that of 2020-2022 moved northeast. In summary, the gravity center of ecological vulnerability showed a trend of moving southwestward from 2000 to 2022, indicating that the intensification degree of ecological vulnerability in the southwest regions was greater than that of the northeast regions.

To explore the spatial distribution imbalance of ecological vulnerability and the migration trajectory of its gravity center, this study utilized the mean center tool in ArcGIS v.10.7 to calculate the gravity centers of ecological vulnerability of the Yellow River Basin in the periods of 2000-2004, 2005-2009, 2010-2014, 2015-2019, and 2000-2022. The gravity centers of ecological vulnerability from 2000 to 2022 were mostly located in the third quadrant (Fig. 6). From 2000 to 2004, the gravity center of ecological vulnerability was in the southwestern direction of 62°27′26′′N, and the migration distance was 5.50 km. The migration distance of the gravity center from 2005 to 2009 was the shortest, with 3.39 km in the northwestern direction of 33°15′03′′N. From 2010 to 2014, the gravity center of ecological vulnerability migrated 10.59 km in the southwestern direction of 29°47′15′′N. The migration distance from 2015 to 2019 was the longest, with 16.49 km in the southwestern direction of 25°29′59′′N, whereas that moved 10.85 km in the northeastern direction of 45°42′04′′N from 2020 to 2022.

Due to the large terrain fluctuation in the Yellow River Basin, there were also great differences in climate, land use type, and vegetation coverage. According to the above geographical spatial differentiation patterns, this study divided the Yellow River Basin into three sub-regions (Fig. 7): Region I (plateau semi-arid area), Region II (mid-temperate arid area), and Region III (warm temperate semi-humid area) (Zhang et al., 2021). Using geodetector, the dominant single factors and interaction factors in the different sub-regions in 2000, 2010, and 2022 were determined using q-value (Table 3).

The dominant single factor of ecological vulnerability in the Yellow River Basin in 2000 was NDVI (q=0.784), and the dominant interaction factors were temperature∩NDVI (q=0.884) (Fig. 8). This study revealed that in 2000, the dominant single factor varied across different regions, but the dominant interaction factors were consistently temperature∩NDVI, with high q-values for all three sub-regions and the whole basin. In 2010, the dominant single factor across different regions was NDVI, and the dominant interaction factors were temperature∩NDVI. In 2022, the dominant single factor across different regions was temperature, the dominant interaction factors were temperature∩precipitation except Region Ⅰ.

At present, scholars primarily use MODIS data to invert ecological remote sensing index based on multiple indicators to reflect the comprehensive condition of regional ecological environment. However, the impacts of human activities on ecological environment are often ignored (Duo et al., 2023). In this study, the socio-economic indicator indicated by nighttime light data was introduced to reflect the interference from human activities. The previous studies have utilized NDVI or fraction vegetation coverage (FVC) to reflect the vegetation condition of ecological environment (Yu et al., 2022b). Although NDVI or FVC was the most effective parameter for characterizing vegetation changes and could better reflect changes in vegetation greenness, it could not indicate the core canopy structure of vegetation growth season status (Wang et al., 2024b). Vegetation NPP could reflect the regional vegetation production capacity in the natural environment. In this study, NDVI, LAI, and NPP were introduced to the proposed greenness index to comprehensively reflect vegetation ecosystem condition (Cheng et al., 2022). Moreover, under the combined actions of global change and urbanization, land degradation has become an increasingly serious ecological problem, which should not be ignored in the evaluation of ecological vulnerability (Li et al., 2022). In addition, the soil erosion was still a serious ecological and environmental problem, which greatly affected the health of regional ecosystems (Zhao, 2023). Therefore, the Albedo, TGSI, and SI were introduced to indicate the condition of salinization, desertification, and soil erosion. The evaluated accuracy in this paper was higher than that of most previous studies in arid and semi-arid areas, which just considered the vegetation, land surface temperature, soil moisture, and surface albedo (Li et al., 2022; Yu et al., 2022). In addition to the traditional remote sensing evaluation index system, this study introduced the regional typical land degradation monitoring indicators (soil erosion, desertification, and soil salinization) to enhance the evaluation accuracy.

The spatial distribution of different levels of ecological vulnerability differed greatly. The slight and mild vulnerability areas were predominantly located in the middle and eastern parts of Qinghai Province, the northern part of Sichuan Province, the middle of Shaanxi Province, and the northwestern part of Henan Province (Niu et al., 2024). The primary explanation for this finding is that the western part of the vulnerable area is located in the plateau, with this area being characterized by a high altitude, low temperature, less precipitation, and being rarely disturbed by human activities (Niu et al., 2022). During the past decades, more ecological protection measures, namely ecological protection zone, returning farmland to forests and grasslands, return grazing land to grassland have been implemented to improve the ecological conditions in the upper reaches of the Yellow River Basin. Moreover, with the climate becoming warmer and more humid, regional vegetation ecosystem had been rapidly restored (Sun and Wang, 2022). In terms of the eastern part of the basin, the level of precipitation is sufficient, which is conducive to vegetation growth (Guo et al., 2024).

The moderate vulnerability areas were mainly concentrated in the southern part of Gansu Province and the northern part of Shaanxi Province, which was consistent with findings of Wang et al. (2023b) that the junction zone of Gansu Province and Shaanxi Province were covered by mild and moderate vulnerability. In these regions, the precipitation and temperature vary greatly with seasons, with significant evaporation and frequent floods and droughts (Yang et al., 2024). With the rapid expansion of urbanization, regional vegetation ecosystem has been greatly destroyed, which intensified the ecological vulnerability of these regions (Wen et al., 2023a).

The severe vulnerability areas were mostly concentrated in the northern part of Ningxia Hui Autonomous Region, the northern part of Gansu Province, and the southern part of Inner Mongolia Autonomous Region, which is similar with the discoveries from Zhang et al. (2022) that the intensive and severe vulnerability areas were mainly distributed in the junction zone of Ningxia Hui Autonomous Region-Gansu Province-Inner Mongolia Autonomous Region. In these regions, the hydrothermal conditions are poor, which is not conducive to vegetation growth (Wang et al., 2024c). Before 2000, under the stress of drought and unreasonable human development, namely overgrazing and disorderly development of mineral resources, the grassland and forests were greatly destroyed. Due to lack of precipitation, it is difficult to restore vegetation ecosystems (Wang et al., 2023b). In addition, the Loess Plateau is located in these regions, where the intensive and severe soil erosion are widely distributed (Niu et al., 2022). In the northeastern part of the basin, the precipitation was scarcity, the dryness is relatively high, and the vegetation is sparse, which leads to severe desertification in these regions combined with the unreasonable human activities, such as overgrazing (Wang et al., 2024c). It finally exacerbated the fragility of the ecological environment (Wen et al., 2023b).

From 2000 to 2022, natural environmental factors were the key driving factors of changes in ecological vulnerability. The dominant single factor in 2000 and 2010 was NDVI, while that of 2022 was temperature. The primary explanation for this finding is that vegetation greatly contributed to regulating ecosystem functions, including the hydrothermal cycle and exchange (Kong et al., 2024). Vegetation has the functions of windbreak, sand fixation, water conservation, and soil stabilization, which could mitigate the land degradation process, such as soil erosion and desertification (Sun and Wang, 2022). In 2022, temperature became the dominant factor of ecological environment vulnerability, due to the gradual increase in average temperature with global warming. The increasing temperature would lead to the instability of water supply, increasing evaporation in the basin, frequent droughts in the upper reaches, a decline in water level in the lower reaches, significant sediment deposition, and an irreversible impact on the environment (Bai et al., 2023). In addition, in alpine regions, the rising temperature is conductive to the photosynthesis of vegetation, which would promote the restoration of vegetation ecosystems (Kong et al., 2024).

The dominant interaction factors in 2000 and 2010 were temperature∩NDVI and the dominant interaction factors in 2022 were temperature∩precipitation. The influence of temperature on the ecological environment ran through the whole process. The main difference was that the influence of precipitation on ecology gradually increased under the influence of global warming. Due to the increase in evaporation, the water level of the basin decreased, and drought and shortage of vegetation and crops were significant. An increase in precipitation would lessen the basin's vulnerability (Wu et al., 2022). Under the combined actions of temperature, precipitation became the driving interaction factor with the greatest explanatory power. The hydrothermal resources were the core factor that affected the processes of regional ecosystems in the semi-arid and arid areas, while in the humid eastern parts of the basin, the water resources also played an important role in the regional economic development, such as agriculture and industry (Bai et al., 2023).

From 2000 to 2022, the single factor explanatory power of population density remained relatively stable, while its interaction with other driving factors exhibited nonlinear enhancement. This nonlinear enhancement may be attributed to the complex geographical and socio-economic characteristics of the Yellow River Basin (Wang et al., 2023b). As population density increases, the complexity of resource demands, environmental pressures, and socio-economic activities also grows, potentially leading to nonlinear effects (Niu et al., 2022; Wu et al., 2022). For instance, high population density may cause interactions with factors such as land use changes, pollution emissions, and water resource utilization to become more pronounced, resulting in these interactions becoming nonlinear rather than linear (Sun and Wang, 2022). This nonlinear enhancement reflected the increased complexity and significance of interactions among driving factors in areas with higher population density. The rapid increase in population would have a double-edged sword effect on the protection of regional ecological environment. On one hand, the disorderly mining, overgrazing, steep slope cultivation, and rapid urbanization process would greatly destroy the regional ecosystem, resulting in serious soil erosion, vegetation damage, water pollution, and desertification (Duo et al., 2023). On the other hand, with the rapid development of economy, more environmental policies and measures have been implemented, such as afforestation, returning farmland to forests and grasslands, which have greatly promoted the restoration of regional ecosystems (Bai et al., 2023; Kong et al., 2024).