Effects of multi human disturbance on soil organic carbon: mediation by landscape patterns and forest community attributes in forest villages of the Henan section, Yellow River Basin

National Forest Villages (NFVs) represent a policy-relevant rural greening initiative with potential contributions to soil organic carbon (SOC) maintenance, yet how multi human disturbances are associated with forest SOC remains poorly understood. We investigated all 99 NFVs in the Henan section of the Yellow River Basin by integrating field surveys, remote sensing data, and statistical approaches. Forest SOC concentration (SOCc), SOC density (SOCs), and forest-to-cropland enrichment ratios (RSOCc and RSOCs) were quantified using paired forest–cropland plots. Human disturbance was characterized by land use human disturbance (LUHD) and multiscale socioeconomic human disturbance (SEHD), while landscape patterns and forest community attributes were quantified independently. After controlling for geoclimatic variation, partial correlations, partial redundancy analysis, variation partitioning, and piecewise structural equation modeling were used to identify their associations with SOC. Weak LUHD NFVs generally maintained 25%–35% higher SOCc and relative SOC enrichment than severely disturbed NFVs, whereas SOCs and landscape diversity were highest under moderate LUHD. Increasing LUHD was also associated with pronounced changes in forest landscape configuration, including reduced forest dominance and greater forest isolation under severe disturbance. Human disturbance and landscape patterns independently explained 13.0% and 7.9% of SOC variation, respectively, whereas the unique contribution of forest community attributes was smaller and non-significant (3.8%). Piecewise SEM further showed that landscape configuration had the largest total association with overall SOC variation. LUHD and urban socioeconomic activity exhibited opposing direct and landscape-mediated associations with SOC. These findings highlight landscape configuration as a key linkage between human disturbance and forest SOC and suggest that maintaining forest dominance, limiting excessive impervious surface expansion, and conserving tree diversity are important for coordinating rural development with SOC management in NFVs.

As global warming is projected to exceed the 1.5 °C threshold, the world is approaching its first major climate tipping point (Lenton et al., 2025). Excessive CO 2 emissions are widely recognized as one of the primary drivers of global warming (IPCC, 2021). To mitigate climate change, carbon reduction targets and mitigation strategies have been widely implemented worldwide (Fawzy et al., 2020). China has committed to peaking carbon emissions before 2030 and achieving carbon neutrality by 2060, commonly referred to as the "dual-carbon" goals. Against the backdrop of increasing carbon emissions and ongoing socioeconomic transformation, both policymakers and researchers are actively exploring effective pathways toward carbon peaking and carbon neutrality (Zhao et al., 2022). Carbon neutrality can generally be achieved through two complementary approaches: reducing carbon emissions and enhancing carbon sequestration (Zhan et al., 2024). Carbon sequestration refers to the process of removing atmospheric carbon and storing it in long term reservoirs, particularly soils (Anjum et al., 2022). Globally, soils store approximately two to three times more carbon than terrestrial vegetation and nearly twice as much as the atmosphere (Davidson et al., 2000). As the most significant carbon reservoir, soil organic carbon (SOC) is of paramount importance for achieving carbon neutrality (Swift, 2001; Sykes et al., 2019). Since topsoil is extremely sensitive to environmental changes, including human disturbances and climate change (Yang et al., 2014), accurately quantifying topsoil SOC is essential. The most commonly used SOC stock parameters are SOC concentration and SOC density, which are referred to as SOC storage (Li et al., 2024). The increase in SOC stocks is referred to as SOC accrual, which is related to C sequestration in soils (Don et al., 2024). Thus, besides absolute SOC concentration and SOC density, introducing parameters characterizing the relative increase in SOC stocks should be a comparatively accurate approach to quantifying SOC storage.

Human disturbance, a concept first introduced by Jalas in 1955, describes the influence of human activities on terrestrial ecosystems. Commonly adopted indices for quantifying human disturbance include the human activity intensity index, the landscape development intensity index, and the Modis Global Disturbance Index (Yang and Song, 2021). In urban ecosystems, the proportion of impervious surfaces and the proportion of built-up areas are typically used to quantify human activity intensity or urbanization intensity. When considering various ecosystems in a region, all land cover types are employed to assess human disturbance. Moreover, socioeconomic factors are strong drivers of land use related human disturbance (Liu et al., 2023; Wang et al., 2022). Numerous previous studies have reported that land-use-related human disturbance significantly leads to the decline of the SOC storage (Lv et al., 2018; Pyles et al., 2022). These studies have rarely addressed the impacts of social and economic development. Given that socioeconomic factors have great potential to reflect human activities and have been frequently used in carbon emission estimation (Kii et al., 2024; Zhao et al., 2024), they offer an opportunity to explore the relationships between SOC and socioeconomic factors. Generally, land use human disturbance (LUHD) characterizes anthropogenic pressure associated with local land use composition and intensity (Ai et al., 2024; Wang et al., 2022; Zhou Y. et al., 2018), whereas socioeconomic human disturbance (SEHD) represents broader pressures associated with population, infrastructure, transportation, and economic activity (Liu et al., 2023). Considering SEHD alongside LUHD provides an additional socioeconomic dimension for characterizing anthropogenic pressure. In addition, the influences of climate and terrain factors on soil properties cannot be overlooked and have been utilized in soil mapping (Sothe et al., 2022a; Sothe et al., 2022b). Therefore, when evaluating the net effect of human disturbance on SOC storage, appropriate scientific research methods should be employed to statistically control for major geoclimatic variation, especially at large scales. The relevant research results will help fill the gap in making large scale predictions of SOC storage based on multiscale human disturbance using big data.

With the increasing availability of open-access big data, a variety of social and economic factors have become easily accessible and mutually complementary. Socioeconomic factors have been widely employed to characterize human disturbance and to investigate their impacts on terrestrial ecosystems, as well as the associated characteristics and functions (Ai et al., 2024; Yang and Song, 2021). Firstly, socioeconomic conditions serve as key drivers of land use variation (Wang et al., 2022). Land use human disturbance (LUHD) exerts significant influences on landscape patterns, leading to increased landscape fragmentation and a greater diversity of landscape types (Wu et al., 2021). Over the past decade, landscape patterns have been consistently demonstrated to regulate forest community attributes (Li et al., 2020; Schindler et al., 2013; Yang et al., 2019). Meanwhile, the effects of forest community attributes on SOC have been a focus of research for decades (Olsson et al., 1996; Tian et al., 2023; Wei et al., 2021a; Wei et al., 2021b). Furthermore, the impacts of socioeconomic human disturbance (SEHD) on forest community attributes (Yang et al., 2022), LUHD on forest community attributes (Beche et al., 2022), and landscape pattern on SOC (Liu et al., 2024; Takoutsing and Rodríguez Martín, 2025; Wei et al., 2025) have also been elucidated. To clarify the driving pathway of multi human disturbance (SEHD and LUHD) on SOC through landscape patterns and forest community attributes, well-designed research that integrates socioeconomic factors, land use human disturbance, landscape patterns, forest community attributes, and SOC is essential.

With rapid social economy development, land use related human disturbance has increasingly shaped rural landscape patterns. On average, more than 40% of landscape fragmentation is caused by the growth of built-up land (Li G. et al., 2021). Under the “rural revitalization” target, China's central government launched the National Forest Villages (NFVs) initiative from 2019 to 2020 (Gong et al., 2022). The aim was to simultaneously safeguard the living space of villages under external pressures, restore degraded ecological environments, and improve residents' quality of life. Among the 7,586 designated National Forest Villages, Henan Province hosts 503 (7% of the total), ranking first nationwide. For a long time, Henan Province has been one of the major grain producing areas in China. In 2021, the grain yield in Henan Province reached 65.42 million tons, accounting for almost 10% of the national total grain output, ranking second among all provinces in China (Wang Z. et al., 2023). Given Henan Province’s plans to construct additional forest villages in the future, it means more forest would be afforested. It is imperative to evaluate the forest ecological function of villages on the list. Specifically, comparing forest SOC storage with adjacent long term cultivated cropland is an effective approach to determine whether forests store more SOC. Moreover, understanding the relationships between forest SOC storage and biotic and abiotic factors is urgently needed to inform the management and construction of forest villages.

As a vital link between land and sea, rivers are pivotal in geochemical cycles, land–ocean interactions, and human activities (Liu et al., 2021). The Yellow River, the world’s fifth longest and China’s second longest (5,464 km), drains a catchment area of 795,000 km 2. However, rapid population growth and urban expansion have driven continuous ecological degradation, making it one of China’s most severely eroded regions (Xu et al., 2023). The Henan section of the Yellow River Basin (YRB) is a key population and economic hub and plays a central role in the basin’s ecological protection and high-quality development strategy (Jiang et al., 2021). Yet, climate change and human activities have recently posed major challenges, including a fragile ecological environment and imbalanced regional development (Fan et al., 2024). From 1990 to 2018, socioeconomic development in the lower YRB (Henan section) intensified gradually: water bodies, paddy fields, and grasslands were characterized by outflow, while urban land, rural settlements, and other construction land expanded dominantly via inflow (Xiao et al., 2020). Clarifying the pathways through which socioeconomic factors affect SOC is therefore critically needed, providing a scientific foundation for sustainable urbanization in this region.

This study hypothesizes that SOC, forest community attributes, and landscape patterns vary across land use human disturbance (LUHD) intensity gradients, and that there are one-way coupling relationships among them, ultimately driven by socioeconomic human disturbance (SEHD). Specifically, this study aims to address the following questions: (1) Do SOC, forest community attributes, and landscape patterns differ significantly among LUHD intensity gradients? (2) What are the potential relationships between SOC and those biotic and abiotic factors? (3) What are the potential pathways through which SEHD influences SOC via LUHD, forest community attributes and landscape patterns? The findings of this study will provide a scientific basis for enhancing SOC storage in the YRB and offer support for achieving future carbon neutrality through optimized construction of NFVs.

Study area

The Henan section of the Yellow River Basin (YRB), located in northern central Henan Province, covers the middle and lower reaches of the YRB (Figure 1). The study area extends from 33°41′ to 36°06′ N and from 110°21′ to 116°06′ E, and lies within the transitional zone between the second and third geomorphic terraces of China (Figure 1). The region is characterized by diverse and complex landforms, with mountainous and hilly terrain mainly distributed in the western part, whereas plains dominate the eastern region (Figure 1). Administratively, the study area encompasses eight cities: Sanmenxia, Jiyuan, Luoyang, Zhengzhou, Xinxiang, Kaifeng, Anyang, and Puyang. Together, these cities cover approximately 36,200 km 2 (Figure 1). The number of NFVs in each city is presented in Table 1.

The study area has a temperate continental monsoon climate, with most precipitation occurring between May and September. Based on the 1980–2020 gridded climate dataset (Table 2) used in this study, mean annual precipitation and mean annual temperature across the study area were 661.8 mm and 13.8 °C, respectively. Vegetation patterns vary by topography: broadleaf forests dominate mountainous regions, while croplands prevail in the east-central area (Fan et al., 2024).

Field survey and soil sampling

Field investigations and soil sampling were carried out from July to September 2023 in all 99 National Forest Villages (NFVs) within the study area (Figure 1). Within each NFV, the largest contiguous forest patch dominated by tree species was selected, and a 20 m × 20 m quadrat was established. Tree species, diameter at breast height (DBH; 1.3 m), tree height, shrub species, shrub basal diameter, shrub canopy dimensions (measured as the maximum crown diameter and the perpendicular crown diameter), and shrub height were recorded for diversity and biomass estimation. The geographic position and altitude were recorded at the center of each quadrat using a handheld GPS device (Jisibao G128BD, UniStrong, Beijing, China).

To compare forest SOC with a locally matched agricultural reference, an adjacent cropland plot under a corn–wheat rotation system and continuously cultivated for more than 15 years was selected in close proximity (mean distance < 200 m) to each forest quadrat. Four random soil-sampling points were established in both the forest quadrat and its paired cropland plot. At each point, the surface litter layer was removed, and soil was collected from a depth of 0–20 cm using a 100 cm 3 cutting ring; samples from the four points were composited into one representative sample (400 cm 3) for forest and paired cropland, respectively. Forest and cropland soil samples were collected using the same procedure, air-dried to constant mass, and sieved through 2-mm and 0.25-mm meshes. Soil bulk density was measured using the cutting-ring method. Visible roots and rock fragments were manually removed before the soil samples were air-dried and weighed.

Determination of SOC parameters

Soil bulk density (BD, g·cm −3) was calculated as:

where M is the air-dried soil mass (g), and 400 cm 3 is the volume of the soil core.

Soil organic carbon concentration (SOCc, g·kg −1) was measured using the heated dichromate/titration method (Bao, 2000). Forest SOC density (SOCs, t·ha − 1) for the 0–20 cm soil layer was calculated as:

where D is soil depth (20 cm), and 0.1 is the unit conversion factor from g·cm − 2 to t·ha −1.

Two dimensionless forest-to-cropland enrichment ratios were calculated as:

and

Thus, SOCc and SOCs represent the absolute SOC status of forest soils, whereas RSOCc and RSOCs represent the relative SOC enrichment of forest soils compared with their paired cropland references.

Determination of forest community attributes

Tree richness and shrub richness were quantified as the number of species. Simpson’s index, Shannon–Wiener’s index, and Pielou’s evenness index were calculated for tree and shrub layers. Aboveground biomass of trees and shrubs was estimated using species-specific allometric equations developed by Zhou G. et al. (2018) and Xie et al. (2018). Biomass was converted to carbon using a conversion coefficient of 0.5 (Lv et al., 2018).

Land use human disturbance (LUHD)

To characterize land-use-related anthropogenic pressure, we quantified LUHD from the proportional areas of eight land use types within circular buffers of 1–10 km surrounding each NFV and their corresponding normalized disturbance coefficients (Ai et al., 2024; Wang et al., 2022). LUHD was treated separately from socioeconomic human disturbance (SEHD) because it represents the realized local land use state rather than broader demographic, economic, and infrastructure pressures:

where n is the number of land use types (n = 8), S i is the area of the ith land use type within the buffer, S is the total buffer area, and D i is the corresponding normalized disturbance coefficient. The coefficients were bare land = 0, cropland = 0.67, forest = 0, grassland = 0.33, impervious surface = 1, shrubland = 0.33, water = 0.33, and wetland = 0.33. LUHD was evaluated at 1 km intervals from 1 to 10 km. Because LUHD is bounded within the theoretical range 0–1, this range was divided into three equal-width intervals to define disturbance grades: weak (0 ≤ LUHD ≤ 1/3), moderate (1/3 < LUHD ≤ 2/3), and severe (2/3 < LUHD ≤ 1). The same fixed thresholds were applied at every buffer distance; the cutoffs were not recalculated separately by scale. Thus, the number of NFVs assigned to each grade could vary among buffers because LUHD values changed with the surrounding land use composition, not because the classification criterion changed. Notably, LUHD represents a land-use-based proxy for anthropogenic pressure rather than a measure of temporal land-use/land-cover change.

Landscape metrics calculation

Based on the 10 m land use raster data used in our previous work (Wei et al., 2025), 13 landscape metrics were calculated within a 10 km radius around each NFV using FRAGSTATS 4.2 (McGarigal, 2014). The metrics represented overall landscape diversity and the density, dominance, shape complexity, and isolation of cropland, forest, and impervious surface patches. Definitions and ecological interpretations are summarized in Table 2.

Data statistics and analysis

LUHD grade comparisons across spatial scales

One-way ANOVA followed by Fisher’s LSD test was used to compare SOC and forest community attributes among weak, moderate, and severe LUHD grades at each of the 1–10 km buffers. Complete multiscale results are provided in Supplementary Figures S2–S5. Landscape metrics, which were calculated within a 10 km radius, were compared among the three LUHD grades at the 10 km scale (Figure 2). Because the overall differences in SOC and forest community attributes among LUHD grades were broadly consistent across the 1–10 km buffers, LUHD at 10 km was retained for subsequent multivariate analyses to match the spatial scale of the landscape metrics.

Partial Pearson correlations and scale selection for scalable SEHD indicators

To account for geographic and climatic variation that may influence spatial patterns of SOC, GeoClimate_PC1 was extracted by PCA from longitude, latitude, altitude, annual mean temperature, and annual precipitation, explaining 77.9% of the total variance. Partial Pearson correlations, controlling for GeoClimate_PC1, were calculated between SOC and forest community attributes, landscape metrics, LUHD, and SEHD indicators.

Socioeconomic human disturbance (SEHD) was characterized using eight indicators: distance from the city center, road length, traffic, nighttime light, drainage-canal length, GDP density, potential crop yield, and population density. Except for distance from the city center, SEHD indicators were evaluated across 10–50 km buffers. For each scalable SEHD indicator, the buffer with the highest mean absolute partial correlation across the four SOC variables was retained. Selected scales were 10 km for road length and population density, 40 km for traffic and nighttime light, and 50 km for drainage-canal length, GDP density, and potential crop yield. Details and raw results are provided in Supplementary Section S1 and Supplementary Figures S6–S9.

Multivariate reduction and predictor selection

PCA was conducted separately for SEHD indicators, landscape metrics, forest community attributes, and SOC variables, with principal components retained until the cumulative explained variance reached ≥85%. All variables were standardized prior to PCA. SOC_PC2 was sign-oriented so that higher scores represented stronger relative forest-to-cropland SOC enrichment.

Partial redundancy analysis (pRDA), controlling for GeoClimate_PC1, was then used to quantify the associations between the retained predictor components and SOC_PC1–PC2. Variation partitioning was further performed to separate the unique and shared contributions of human activity (LUHD and SEHD PCs), landscape pattern, and forest community attributes. Additional details are provided in Supplementary Section S2.

Piecewise structural equation modeling and robustness assessment

Before conducting piecewise structural equation modeling (Piecewise SEM), the small sample corrected Akaike information criterion (AICc) was used to select parsimonious terminal SOC component models before their incorporation into the piecewise SEM. Details and results can be found in Supplementary Section S3 and Supplementary Table S2.

Piecewise structural equation models were constructed using an a priori hierarchical framework linking SEHD, LUHD_10 km, landscape pattern, forest community attributes, and SOC. GeoClimate_PC1 was included as a control covariate in all component equations. Model fit was evaluated using Fisher’s C statistic, with P > 0.05 indicating acceptable global fit. Standardized path coefficients and R 2 values were reported, and bootstrap resampling was used to assess the robustness of pathway estimates. Details of the bootstrap resampling process are provided in the Appendix.

ANOVAs and Fisher’s LSD post hoc tests were performed using IBM SPSS Statistics for Windows, version 26.0 (IBM Corp., 2019). All other statistical analyses were conducted in R version 4.4.1 (R Core Team, 2024). The R packages used included ppcor v1.1 (Kim, 2015), stats (R Core Team, 2024), vegan v2.7-5 (Oksanen et al., 2026), lmtest v0.9-40 (Zeileis and Hothorn, 2002), performance v0.13.0 (Lüdecke et al., 2021), corrplot v0.95 (Wei and Simko, 2024), MuMIn v1.48.4 (Bartoń, 2024), and piecewiseSEM v2.3.0.1 (Lefcheck, 2016).

Differences in landscape metrics, forest community attributes, and SOC among LUHD grades

Land use human disturbance (LUHD) intensity within individual NFVs remained relatively stable across the 10 buffer zones (Supplementary Figure S1). In the 3–4 km and 7–10 km buffers, forest soil organic carbon concentration (SOCc) in weak LUHD NFVs was approximately 1.25 times that in severe LUHD NFVs (P < 0.05) (Figure 3 and Supplementary Figure S2). Similarly, in the 4 km and 7–10 km buffers, SOC density (SOCs) in NFVs with moderate LUHD was 1.16 times higher than in those with severe LUHD (P < 0.05). Across all buffer zones, RSOCc in NFVs with weak LUHD was 1.35 times greater than in those with severe LUHD (P < 0.05), and RSOCs in weak LUHD NFVs was 1.27 times greater than in severe LUHD NFVs (P < 0.05) (Figure 3 and Supplementary Figure S2).

Tree richness differed significantly among LUHD grades across most 1–10 km buffers, except at 7 and 8 km (P < 0.05; Figure 4 and Supplementary Figure S3). Specifically, tree Simpson diversity differed at 5 and 6 km, whereas tree Shannon–Wiener diversity and Pielou evenness did not differ significantly across the full 1–10 km range. Overall, tree diversity (encompassing richness and evenness) was higher in NFVs with weak LUHD than in those with severe LUHD. For shrub communities, both richness and Shannon diversity were consistently higher under weak than moderate or severe LUHD across all 1–10 km buffers (P < 0.05; Supplementary Figure S4). Tree aboveground biomass carbon was higher under severe than weak LUHD, whereas shrub aboveground biomass carbon showed the opposite pattern (P < 0.05; Supplementary Figure S5).

Landscape metrics varied significantly among the three LUHD grades (weak, moderate, and severe) (Figure 2). The moderate LUHD NFVs exhibited the highest landscape diversity, which was 1.4- and 1.9-fold higher than that in the weak and severe LUHD NFVs, respectively (P < 0.05). Cropland patch density and landscape shape index followed the order: weak > moderate > severe (P < 0.05). In contrast, the cropland largest patch index showed the opposite trend: severe (60.1%) > moderate (40.4%) > weak (4.2%) (P < 0.05). The Euclidean Nearest−Neighbor Distance of cropland patches in weak LUHD NFVs was 1.3-fold higher than that in moderate and severe LUHD NFVs (P < 0.05). Forest patch density and landscape shape index were highest in moderate LUHD NFVs, followed by weak and severe LUHD NFVs (P < 0.05). The Euclidean Nearest-Neighbor Distance of forest patches in severe LUHD NFVs was 5- to 7.7-fold higher than that in weak and moderate LUHD NFVs (P < 0.05). The largest patch index of forest patches in weak LUHD NFVs was 60.3%, significantly higher than that in moderate (7%) and severe (0.08%) LUHD NFVs (P < 0.05). Impervious patch density, largest patch index, and landscape shape index were lower in weak LUHD NFVs than in moderate and severe LUHD NFVs (P < 0.05). The Euclidean Nearest– Neighbor Distance of impervious surface patches in weak LUHD NFVs was the highest among all three grades (P < 0.05).

Partial correlations between SOC and biotic and abiotic factors

After controlling for GeoClimate_PC1, tree richness showed the most consistent positive relationships with SOCc and SOCs, whereas shrub richness was negatively related to SOCs and RSOCs (Figure 5, Supplementary Figure S6). Among landscape metrics, SHDI, cropland patch density, cropland landscape shape index, forest patch density, and forest landscape shape index were generally positively associated with SOC, while forest nearest neighbor distance showed negative associations (Figure 5, Supplementary Figure S7). Correlations between LUHD and the four SOC parameters were generally weak, indicating that the overall LUHD–SOC relationship could not be adequately characterized by a simple bivariate association (Figure 5, Supplementary Figure S8). SEHD indicators were mainly related to SOCc and SOCs, and these associations were scale-dependent (Figure 5, Supplementary Figure S8). The mean-absolute-partial-r criterion retained road length at 10 km, traffic at 40 km, nighttime light at 40 km, drainage canal length at 50 km, GDP at 50 km, potential crop yield at 50 km, and population density at 10 km; city distance was retained as a non-buffered proximity variable (Supplementary Figure S9). Distance from city center was positively related to RSOCc (P < 0.05), and nighttime light index and GDP were positively related to SOCc and SOCs (P < 0.05). Road length and potential crop yield were negatively related to SOCs, whereas drainage canal length was negatively related to both SOCc and SOCs (P < 0.05) (Figure 5). The main text heatmap summarizes the retained variables, while full partial correlation matrices are provided in Supplementary Figures S6–S8.

Partial RDA and variation partitioning of the two SOC components

PCA was performed to reduce the dimensionality of the raw data and create new variables for subsequent analyses, including partial RDA and variation partitioning. The new variables, major contributors, and ecological meanings of high scores are shown in Table 3 and Supplementary Table S1.

After controlling for GeoClimate_PC1, partial RDA revealed that biotic and abiotic factors significantly explained the variation in forest SOC (overall adj. R 2 = 0.182, P < 0.001; Figure 6a). SOC_PC1 and SOC_PC2 were significantly explained by biotic and abiotic factors, with adjusted explanatory rates of 24.4% (P < 0.001) and 12.1% (P = 0.017), respectively. Among individual predictors, Landscape_PC1, SEHD_PC2, and SEHD_PC3 showed significant independent effects on SOC variation, explaining 3.24% (F = 4.50, P = 0.04), 4.59% (F = 5.96, P = 0.03), and 3.36% (F = 4.64, P = 0.04) of the variation, respectively (Figure 6a).

Landscape_PC1 was positively associated with SOC_PC1, indicating that areas characterized by greater forest dominance and cropland fragmentation tended to have higher overall forest SOC (Figure 6a). SEHD_PC2 showed a positive relationship with SOC_PC1 and a negative relationship with SOC_PC2. In contrast, SEHD_PC3 was positively associated with SOC_PC2, indicating that increasing remoteness and reduced human accessibility were associated with stronger relative forest-to-cropland SOC enrichment (Figure 6a). Other biotic and abiotic factors, including forest community variables and remaining landscape or human disturbance components, showed no significant independent contributions (Figure 6a).

Variation partitioning further demonstrated that human disturbance and landscape factors were the dominant explanatory groups after removing the effects of geo-climatic variation (Figure 6b). Human disturbance independently explained 13.0% of SOC variation (P < 0.001), followed by landscape factors (7.9%, P = 0.005), whereas forest community characteristics explained 3.8% of the variation but did not reach significance (P = 0.178). The shared contributions among forest community, landscape, and human disturbance were relatively small (0.7%–1.2%) (Figure 6b). Together, these results suggest that forest SOC variation was primarily shaped by landscape configuration and human disturbance gradients.

Piecewise structural equation models for the two SOC components

Structural equation modeling (SEM) was used to quantify the direct and indirect effects of human disturbance, landscape metrics, and forest community attributes on forest SOC variation (Figure 7). The final model showed a good fit to the observed data (Fisher’s C = 0.051, P = 0.975; AICc = 1420.57) and explained 34.4% of the variation in SOC_PC1 (R 2 = 0.344). Landscape_PC1 showed a significant positive direct effect on SOC_PC1 (standardized path coefficient = 1.087, P < 0.001), with an additional positive indirect effect (0.153) through Forcommunity_PC2, resulting in the largest total effect among all predictors (1.241). In contrast, landscape_PC3 showed a negative direct effect on SOC_PC1 (−0.202, P = 0.026). Human disturbance exhibited complex direct and indirect effects. LUHD had a strong positive direct effect on SOC_PC1 (1.096, P < 0.001); however, its indirect effect was negative (−0.923), resulting in a much weaker total effect (0.173). Similarly, SEHD_PC2 showed a positive direct effect on SOC_PC1 (0.222, P = 0.018), but a negative indirect effect (−0.049) through Landscape_PC3, leading to a total effect of 0.173. Forest community attributes showed a weaker but significant positive effect on SOC_PC1, with Forcommunity_PC2 exhibiting a standardized direct effect of 0.183 (P = 0.044) and an equivalent total effect (0.183). SEM was not constructed for SOC_PC2 because model selection retained only one predictor, SEHD_PC2 (R 2 = 0.087, adjusted R 2 = 0.068, Supplementary Table S2).

Overall, SEM indicated that landscape metrics showed the strongest associations with SOC_PC1. LUHD and SEHD_PC2 were associated with SOC_PC1through both direct and landscape-mediated pathways. LUHD was negatively associated with Landscape_PC1, which was strongly positively associated with SOC_PC1, whereas SEHD_PC2 was positively associated with Landscape_PC3, which was negatively associated with SOC_PC1. These results suggest complex associations among human disturbance, landscape configuration, forest community attributes, and SOC.

Effects of human disturbance on the NFVs’ SOC

SOC in the National Forest Villages (NFVs) generally showed lower values under severe land use human disturbance (LUHD) than under weak or moderate LUHD. Specifically, SOCc, RSOCc, and RSOCs in NFVs with weak LUHD were 25%, 35%, and 27% higher than those in NFVs with severe LUHD, respectively, whereas SOCs in NFVs with moderate LUHD was 16% higher than that under severe LUHD (Figure 3). These patterns align with previous studies reporting an obvious decline in SOC content from weak to severe LUHD gradients (Pandey, 2021). At the global scale, Xu et al. (2024) observed a "decline-rise-decline" pattern in SOC stocks across urbanization intensity gradients, attributed to land use transitions—a typical form of land use human disturbance historically. In this study, all land use types were integrated to define land use human disturbance intensity. Severe LUHD NFVs were predominantly located in cropland-concentrated areas, whereas weak LUHD NFVs were primarily found in forest-concentrated regions. Thus, the observed SOC differences among LUHD grades were closely associated with spatial variation in current land-use composition. Extensive field surveys have previously assessed soil organic carbon (SOC) storage capacity under varying degrees of anthropogenic disturbance, though most relied on qualitative methods (Arunrat et al., 2022), unlike the quantitative approach adopted in this study. For instance, deforestation in China resulted in substantial carbon losses of approximately 29.81 Tg C yr −1 during 2000–2018 (Chang et al., 2022). Relative to undisturbed natural ecosystems, reductions in SOC have been widely documented in planted forests, grazed grasslands, and conventionally tilled croplands, largely owing to the rapid depletion of particulate organic carbon (POC) (Liu et al., 2025). These associations may be linked to processes involving soil aggregation, soil enzyme activities, and soil fungal communities (Cai et al., 2023; Duan et al., 2025; Song et al., 2026; Zhang et al., 2024). Consequently, future research on NFVs should investigate whether and how land-use-related anthropogenic pressure is associated with soil aggregation and, in turn, SOC storage and relative enrichment.

When considering socioeconomic human disturbance, multivariate analyses indicated that socioeconomic gradients were not uniformly associated with lower soil organic carbon (SOC) (Figure 6a and Table 3). After geoclimatic control, higher urban socioeconomic activity (SEHD_PC2) including traffic density, nighttime light, and GDP, was negatively associated with relative SOC-enrichment component (SOC_PC2) (Figures 5, 6a). In the final piecewise SEM, SEHD_PC2 had a positive direct effect on overall SOC component (SOC_PC1), partly offset by a small negative indirect effect through Landscape_PC3 (Figure 7). However, SEHD_PC2 had a significant negative effect on relative SOC-enrichment component (Supplementary Table S2). Chien et al. (2024) reported that when human population density exceeded 300 people/km 2 (defined as "urban domains" in the European Commission database), mangrove soil carbon stocks were approximately 37% lower than those in isolated mangrove forests, indicating a negative relationship between population density and mangrove soil carbon. Our findings are broadly consistent with this pattern: higher urban remoteness (SEHD_PC3), characterized by greater distance from the city center and lower population density, was positively associated with the forest/cropland SOC ratio (Figures 5, 6a). Song et al. (2024) found that moderate distances from roadways enhance soil carbon storage compared to areas near roads or undisturbed natural landscapes, partially supporting the moderate disturbance hypothesis. In the present study, SOC density was also higher under moderate than severe LUHD (Figure 3), further indicating that the association between disturbance and SOC was not strictly monotonic (Figure 3). When controlling for geoclimatic factors, the piecewise SEM showed opposing direct and indirect LUHD effects on SOC_PC1, resulting in a much weaker total effect (Figures 5, 7). Nevertheless, the LUHD in the present study should be interpreted as a proxy for land-use-related anthropogenic pressure rather than as a direct causal attribution of individual land cover states or changes to human activities. Some land cover patterns may also be influenced by natural processes, which cannot be fully distinguished using static land use data alone (Pickens et al., 2025; Zhu et al., 2022). Future studies could reduce this uncertainty by integrating multi-temporal land use/cover change with independent indicators of direct anthropogenic pressures and spatial records of major natural disturbances, thereby allowing anthropogenic and naturally induced land cover changes to be more explicitly distinguished. Collectively, these findings highlight the complex influence of human disturbance on SOC variability, providing valuable insights for developing targeted SOC management strategies in the context of NFV construction.

Driving pathways through which human disturbance affects NFV SOC via landscape and forest community attributes

After controlling for GeoClimate_PC1, variation partitioning showed that human disturbance and landscape patterns made the largest unique contributions to SOC variation, whereas forest community attributes contributed less and did not reach significance (Figure 6b). Human disturbance independently explained 13.0% of SOC variation and landscape patterns explained 7.9%, compared with 3.8% for forest community attributes. In the piecewise SEM, LUHD and SEHD_PC2 showed both direct and landscape-mediated effects on SOC_PC1, while Landscape_PC1 had the largest total effect and Forcommunity_PC2 contributed a smaller positive effect (Figure 7). Thus, broader socioeconomic disturbance and local land use disturbance were linked to SOC through both direct and landscape-mediated pathways, although the strengths and significance of individual pathways varied.

Firstly, landscape patterns differed among LUHD grades. Previous studies have demonstrated that intensified human disturbance significantly increased landscape fragmentation (+82%), regularity (+56%), and diversity (+37%) in the Dongting Lake Basin (Hou et al., 2024). In the present study, the highest landscape diversity was observed under moderate human disturbance (Figure 2). By contrast, differences in landscape fragmentation and regularity among the three disturbance grades were landscape-type dependent. As LUHD increased, cropland exhibited increased dominance and regularity, but decreased fragmentation; impervious surface showed increased dominance and fragmentation, yet decreased regularity (Figure 2). Forest patch density and landscape shape index were highest under moderate LUHD, whereas the forest largest patch index was highest under weak LUHD (Figure 2). Wu et al. (2021) reported that with continuous intensification of human activities, landscape types on a given scale tend to homogenize, and patches of the same type coalesce into larger, more integrated units, especially near the cores of major urban agglomerations. In this study, higher SEHD_PC2, including traffic density, nighttime light, and GDP, corresponded to greater impervious surface dominance (Landscape_PC3); higher LUHD was associated with lower forest dominance and cropland fragmentation (Figure 7). Li Y. et al. (2021) reported that > 40% of severe landscape fragmentation was mainly associated with the expansion of built-up land. In contrast to the conventional view, the influence of settlement expansion on landscape fragmentation is relatively weaker in large cities than in small- and medium-size cities. Furthermore, the effects of settlement expansion on landscape fragmentation exhibit pronounced regional variability and clear gradient-dependent patterns (Li G. et al., 2021). Consistent with these findings, our results suggest that the influence of human disturbance on landscape patterns in NFVs varies with landscape type. At the regional scale, focusing exclusively on overall landscape pattern without accounting for specific landscape types represents a one-sided approach.

Secondly, landscape metrics exerted both direct and indirect effects on SOC. Our previous work demonstrated that forest patch density is one of the most influential indicators for explaining SOC content and density (Wei et al., 2025). The present study further showed that forest patch density was generally positively associated with the four SOC metrics (Figure 5). Furthermore, SHDI, cropland patch density, cropland landscape shape index, forest patch density, and forest landscape shape index were generally positively associated with SOC, whereas forest nearest-neighbor distance was negatively associated (Figure 5). Liu et al. (2024) reported that at the class level, configuration indices showed weaker correlations with farmland SOC than at the landscape level, among which the two most important were the Mean Fractal Dimension Index and Number of Patches. All these results suggest that landscape configuration was closely linked to SOC. Hou et al. (2024) reported that intensified human disturbance increased landscape fragmentation, regularity, and diversity, resulting in the greatest decline in soil carbon stocks (0.24 Gt, 51%). In our study, moderate LUHD coincided with the highest landscape diversity (Figure 2), which was positively related to tree diversity (Figure 6a), supporting the moderate disturbance hypothesis. Taken together, landscape configuration represented an important indirect pathway linking human disturbance to SOC, with indirect effects that substantially offset the corresponding direct effects in the SEM.

Thirdly, plant diversity showed apparent human disturbance-dependent variation. Weak LUHD NFVs generally had higher tree diversity than severe LUHD NFVs (Figure 4 and Supplementary Figure S3). In contrast, urban forest tree diversity in cities showed a strictly monotonic increase with urbanization intensity (Wang et al., 2018; Xiao et al., 2016; Yang et al., 2022). Meanwhile, the widely held view is that woody plant diversity generally exerts positive effects on SOC (Tian et al., 2023). Lee et al. (2024) identified species richness as the dominant biotic factor directly influencing soil organic carbon (SOC) in temperate forests of South Korea. Previous studies have shown that greater variation in tree species composition, together with reduced community functional diversity, can increase the heterogeneity of SOC chemical composition across forest communities (Wang H. et al., 2023). Beche et al. (2022) further reported that the combined effects of forest edge influence and forest degradation are highly context dependent and are largely determined by local forest-use practices. In the present study, tree aboveground biomass carbon was higher in NFVs with severe LUHD than in those with weak LUHD. Shrub richness was consistently lower under moderate and severe LUHD, while shrub biomass was generally lower under stronger LUHD, with some scale-dependent variation (Supplementary Figures S4, S5). This pattern may partly reflect the prevalence of monoculture plantation forests (rather than diverse natural stands) in intensively disturbed NFVs–stands typically composed of fast-growing, high-yield species such as Populus spp., or regionally common species like Platycladus orientalis and Robinia pseudoacacia. These plantations are often situated adjacent to farmlands as shelterbelts, or positioned at village entrances and other focal sites as feng-shui forests to support rural livelihoods and ecosystem functions (Chen et al., 2018). Future NFV construction and management should therefore integrate such woodlands into sustainable rural development strategies to balance forest SOC maintenance, biodiversity conservation, and local socio-ecological needs.

Implications

As a rural greening initiative related to nature-based solutions (NbS) (Pereira et al., 2024), NFVs provide a policy-relevant context for considering climate-change mitigation and human wellbeing. In recent years, NbS has emerged as an innovative paradigm in global ecological governance, offering effective responses to interconnected societal challenges (Fu and Wang, 2025). The core principle of NbS is to respect natural laws, harness natural processes, and design context-specific ecosystem carbon sink enhancement plans based on the ecosystem approach and local conditions (Yu et al., 2022). NFVs provide a useful setting for examining these relationships because village forests are embedded within mosaics of cropland, settlements, roads, and other land uses, allowing local LUHD and broader SEHD to be evaluated within the same rural landscape system. Our results document substantial variation in forest SOC and relative forest-to-cropland SOC enrichment among NFVs, thereby revealing how landscape configuration, socioeconomic development, and forest community attributes are associated with SOC storage within these rural systems. However, because this study did not include non-forest villages as controls, the observed SOC patterns cannot be used to determine whether NFVs possess a unique SOC storage advantage relative to ordinary villages, nor should they be interpreted as evidence that NFV designation itself enhances SOC. Future studies incorporating comparable non-forest villages under similar climatic, topographic, and land use conditions will be needed to quantify the additional SOC benefits, if any, associated with NFV construction.

Future NFV construction should balance rural socioeconomic development with maintenance of forest SOC. The current analyses suggest that management should prioritize landscape configuration, particularly maintain forest dominance and avoid excessive impervious surface dominance, while recognizing that socioeconomic effects are scale- and component-dependent rather than uniformly negative (Figures 6a, 7). The positive association of tree diversity with SOC-related components also supports maintaining tree diversity and, where appropriate, restoring shrub layers in existing forests. Xu et al. (2023) emphasized that conversion of built-up land is a critical pathway to enhance carbon storage, substantially improving carbon sequestration efficiency and capacity in the Yellow River Basin (YRB). Our study provides a framework for jointly considering socioeconomic development, landscape optimization, forest community restoration, and forest SOC management, thereby fostering better coordination of human–nature relationships in NFV construction and management.

This study showed that human disturbance was closely associated with forest SOC variation in NFVs of the Henan section of the Yellow River Basin. Weak LUHD NFVs exhibited higher SOC concentration and relative SOC enrichment, whereas moderate LUHD showed the highest SOC density and landscape diversity. After controlling for geoclimatic effects, human disturbance and landscape patterns explained more SOC variation than forest community attributes. Piecewise SEM further indicated that landscape configuration played a central role in linking LUHD and SEHD with SOC. These findings suggest that the effects of human disturbance on SOC are complex rather than simply linear associations. Future NFV management should therefore prioritize maintaining forest dominance and connectivity, limiting excessive impervious surface expansion, and conserving tree diversity. Because comparable non-NFV controls were not included, further studies are needed to determine whether NFV construction itself provides additional SOC benefits.