How Land Use Transitions Contribute to the Soil Organic Carbon Accumulation from 1990 to 2020

: Soil organic carbon stock (SOCS) changes caused by land use changes are still unclear, and understanding this response is essential for many environmental policies and land management practices. In this study, we investigated the temporal–spatial and vertical distribution characteristics of the SOCS in the Western Sichuan Plateau (WSP) using the sparrow search algorithm– random forest regression (SSA-RFR) models with excellent model applicability and accuracy. The temporal–spatial variations in the SOCS were modeled using 1080 soil samples and a set of nine environmental covariates. We analyzed the effect of land use changes on the SOCS in the WSP. The total SOCS increased by 18.03 Tg C from 1990 to 2020. The results of this study confirmed a significant increase in the SOCS in the study area since 2010. There was a 27.88 Tg C increase in the SOCS in 2020 compared to the total SOCS in 2010. We found that the spatial distribution of the SOCS increased from southeast to northwest, and the vertical distribution of the SOCS in the study area decreased with increasing soil depth. Forests and grasslands are the main sources of SOCS the total SOCS in the forest and grassland accounted for 37.53 and 59.39% of the total soil organic carbon (SOC) pool in 2020, respectively. The expansion of the wetlands, forest, and grassland areas could increase the SOCS in the study area. A timely and accurate understanding of the dynamics of SOC is crucial for developing effective land management strategies to enhance carbon sequestration and mitigate land degradation.


Introduction
Soil organic carbon (SOC) is a significant indicator to measure soil fertility and evaluate crop growth and development [1,2].It is related to climate change [3], vegetation types [4], land use change (LUC) [5,6], and human activities [7].As land resource management transitions from quantity-focused to quality-focused management [8], it is crucial to study the effect of this transition on the SOC stock (SOCS).In addition, the response of the SOCS to land use change across different soil layers and regions shows significant differences, and it is essential to examine the impact of land use change on the SOCS in both horizontal and vertical directions [9].SOC, pH, and soil bulk density (SBD) are important indicators for estimating the SOCS [10].The accurate estimation of these parameters is crucial for assessing regional changes in the SOCS and understanding carbon cycles in ecosystems.Studying the temporal-spatial characteristic of LUC plays an important part in the exploration of the change in the SOCS.However, only a few studies have focused on the effects of LUC on the SOCS using a long time series [5,6,[11][12][13].Therefore, it is significant to investigate the effects of LUC on the SOCS.
At present, random forest regression (RFR) has been commonly used to estimate the SOC and SBD based on multi-source remote sensing data, and these models have strong applicability [14][15][16].To further improve the accuracy of the RFR model, the sparrow search algorithm (SSA) can be employed, which is known for its easy implementation, versatility, and fast convergence speed [17,18].In addition, it is vital to identify the major environmental variables affecting the SOC content; these variables were utilized as the input dataset for SOC prediction models.However, the optimal combination of these variables for improving the modeling accuracy has different performance in different regions.
The Qinghai-Tibet Plateau (QTP) is an important region for studying carbon sinks and global energy, carbon, and moisture cycles [10,12].However, the effect of LUC on the increase in the SOCS in the QTP is still unclear.It is essential to qualitatively and quantitatively investigate the SOC effects and potential risk assessment.Some studies have reported that LUC can directly affect the SOCS [19,20], while others have indicated that LUC is affected by climate change, indirectly affecting the SOCS [3,12,21,22].Many scholars have regarded the QTP as their study area to investigate the changes in the SOCS with coarse spatial resolution (e.g., 90 m, 250 m, 500 m, 1000 m) and short time series; their results have some limitations in showing the fine spatial and temporal characteristics in the QTP [12,23,24].To address these limitations, this study focused on the Western Sichuan Plateau (WSP), which is located in the southeastern part of the QTP.The WSP is known for its natural ecosystems, including forests, grassland, and wetlands, which have a prominent carbon reduction effect [25].The WSP is a major pastoral area in China, with grasslands playing a critical role [26].The WSP is a sensitive area, in which the vegetation types and soil conditions are closely related to climate conditions [27].Understanding the effects of LUC on the SOCS is important for developing practices that promote agricultural sustainability.However, only a few studies have been conducted on the impact of LUC on the SOCS in the WSP.Therefore, it is important to further investigate the effect of LUC on the SOCS in the WSP.
In this study, we analyzed the impact of LUC on the temporal-spatial variability of the SOCS in the WSP based on long time series field survey data and SSA-RFR.We measured the SOC content, pH, and SBD of the topsoil (0-5, 5-10, 10-20, and 20-30 cm) using data for 54 sampling sites to estimate the spatial and vertical distribution of the SOC and SBD.We analyzed the temporal-spatial and vertical trends of the SOCS under the LUC from 1990 to 2020.The optimal combination of environmental covariates with the highest model accuracy was determined by the different characteristics of the SOC and SBD in different soil layers.Then, SSA-RFR was used to build SOC and SBD estimation models for the WSP; the optimal models were adopted to upscale the field survey samples to the whole study region.The spatial distribution results of the SOC density (SOCD) were mapped with a high resolution (30 m).The objective of this study is to (1) find the optimal combination for the specific region and (2) explore the dynamic change mechanism of the SBD and SOC from 1990 to 2020.These results can be utilized for developing effective land management strategies to enhance carbon sequestration and mitigate land degradation.

Study Areas
The WSP is located in the southeastern part of the QTP (27°28′15.6″-34°18′46.8″N,97°20′56.4″-104°42′50.4″E), with an average elevation of 4177 m (Figure 1a) and a total area of about 2.98 × 10 5 km 2 .The WSP is influenced by the Qinghai-Tibet Plateau and the Hengduan Mountains, exhibiting a typical plateau monsoon climate with dry and rainless conditions and a large temperature amplitude daily [28,29].The average annual temperature was 7.5 °C and the average annual precipitation was 87.61 mm in 2020 (Figure 1b).The main land use types are grasslands and forests (Figure 1c).The main soil types are dark felty soils, felty soils, and brown earths (Figure 1d).Carex myosuroides Vill.and a variety of grass vegetation mostly grow in the dark felty soil, and it is an important area for livestock products [30].The vegetation in the alpine meadows is mainly Kobresia humilis and Carex alatauensis S. R. Zhang, and the accompanying species include Carex tristachya Thunb.and Polygonum macrophyllum D. Don [31].

Field Data Measurements and Processes
We collected a total of 1080 soil samples from August 2018 to August 2020; the sampling sites were evenly distributed and representative (Figure 1c).Considering that forests and grasslands are the primary types of land use in the study area, we have established numerous sampling sites distributed in these areas.In each sampling area, five samples were collected at equidistant locations to obtain average values.Additionally, two 100 m sample lines were established in each sampling area (Figure 2a).For each sample, we collected undisturbed soil samples at fixed intervals of 0-5, 5-10, 10-20, and 20-30 cm in the topsoil layer using steel bulk density rings (Figure 2b).Each soil sample was selected at a weight of 10 g and mixed with distilled water to form a soil-water suspension.The pH value of each sample was measured using the potentiometric method.These samples were ovendried at 105 °C until a constant mass was achieved to determine the SBD (g cm −3 ).The SOC content of the samples was separately determined using potassium dichromate oxidation spectrophotometry.To estimate the SOCS, we calculated the SOCD in a certain soil layer thickness per unit area (1 m 2 ).The different soil structures of various land use types can lead to changes in the accumulation and decomposition processes of SOC, thereby affecting SOCD.The SOCS in the 0-30 cm layer was calculated using Equation (1) [32][33][34].
where A is the area of the study area (m 2 ), SOCD is the density of SOC (t hm −2 ), H is the soil thickness (cm), SBD is the bulk density of the soil layer (g cm −3 ), G is the gravel content of the soil sample (%), and n is the number of soil layers (n = 4).To ensure the accuracy and reliability of the spatiotemporal analysis of the SOCS, many relevant measured data were obtained from several articles to validate and correct the results of this study (Table 1).Moreover, publicly available soil data (including the SBD and SOC data for the 0-30 cm soil layer(s) in the 2010s) were acquired from the above studies, and some of the soil data (including SBD data for the 0-30 cm soil layer) were obtained from the Harmonized World Soil Database (HWSD).This database provided us with additional SBD data for the 0-30 cm soil layer.By incorporating these diverse sources of data, we aimed to enhance the robustness and comprehensiveness of our analysis.

SSA-RFR Model and Accuracy Evaluation
The random forest regression (RFR) algorithm is a machine learning algorithm to estimate the ecological indicators [39][40][41].In this study, 70% of the soil samples were utilized as the training set and 30% of the soil samples were utilized as the testing set to estimate the SBD and SOC with different land types, respectively.In addition, to ensure the optimal and stable model, the SSA was employed to improve the performance of the RFR models; the workflow chart is shown in Figure 3.The SSA draws inspiration from the foraging behavior of sparrows and applies cooperative and competitive strategies to find the most optimal solutions within the solution space [42].The SSA incorporates randomness and local search mechanisms to enhance its ability to explore the solution space and improve the convergence of the algorithm [43].The performance of the RFR models of the SBD and SOC predictions was evaluated using the coefficient of determination (R 2 , Equation (2)), mean absolute error (MAE, Equation (3)), root mean square error (RMSE, Equation ( 4)), and BIAS (Equation ( 5)) [44,45].Both RFR and the SSA are implemented in Python 3.8.(2) where yi and fi represent the field's measured and estimated SOC or SBD values in the ith sample, respectively,  ̅ and  ̅ are the measured and estimated SOC or SBD averaged over all of the soil samples, and n is the number of soil samples.

Satellite Data and Pre-Processing
A total of 224 Landsat Collection 2 Level-2 surface reflectance products with a spatial resolution of 30 m during 1990-2020 were downloaded from the United States Geological Survey Earth Explorer (USGS EarthExplorer) with an image acquisition time that was close to the field measurement time and a cloud cover of less than 10% [40,46].The images were pre-processed using the ENVI 5.3.The vegetation indices were calculated by band math (Table 2).Moreover, a digital elevation model (DEM) was obtained from Shuttle Radar Topography Mission (SRTM) images (version 003, 30 m), which was used to calculate the slope, aspect, terrain ruggedness index (TRI), and topographic wetness index (TWI).Both the above vegetation indices and terrain indicators were used to construct the dataset for the RFR model.The land use change data were obtained from the first Landsat-derived annual land cover product of China (CLCD) with 30 m spatial resolution, covering four periods: 1990, 2000, 2010, and 2020.The meteorological data were downloaded from the European Centre for Medium-Range Weather Forecasts (ECMWF).The meteorological and LUC data were used to analyze the spatial and temporal characteristics of the SOC and SBD.Note: NDVI: normalized difference vegetation index; EVI: enhanced vegetation index; CIg: chlorophyll index-green; and SAVI: soil-adjusted vegetation index.

Statistical Analysis of SBD and SOC
The vertical distribution of the SOC decreased with increasing soil depth, with an average SOC value of 66.12 g kg −1 in the 0-30 cm layer (Figure 4).The largest decrease in the SOC content occurred between 5-10 cm and 0-5 cm and in adjacent soil layers, with a maximum difference of 117.26 g kg −1 .The average SBD value for the 0-30 cm interval in the study area was 0.83 g cm −3 .Both the pH and SBD exhibited significant negative correlations with the SOC (Figure 5).Soil is rich in microorganisms which are suitable to grow in a neutral soil environment.The study area is an acidic soil, in which the decomposition process of organic carbon is inhibited and can increase the stability of the SOC.SOC has the ability to combine with soil particles, enhancing soil stability and reducing SBD.Additionally, higher SOC in the soil results in increased microbial activity and soil biodiversity, which facilitate the breakdown and decomposition of organic matter in the soil.This process contributes to the formation of soil aggregates and improves soil structure and reduces soil compaction.The average SOCD value in the 0-30 cm range was 1.98 t hm −2 .The SOCD decreased with increasing soil depth.And there were noticeable differences in the SOCD among the different land use types (Figure 6).The SOCD can directly influence the SOCS.Vegetation types and land use types strongly influence SOCD and the SOCS.

Applicability of SSA-RFR Model
SSA-RFR models were used to construct the SBD and SOC estimation models for the study area.The indicator significance results for each model also exhibited significant differences (Figure 7).According to the indicator significance results, the DEM we chose mainly has indicators to estimate SBD and SOC in each different soil layer.This is because different spectral bands reflect different vegetation properties and hence soil properties.Topographic indices can reflect soil properties and compensate for the limitations of remotely sensed data in obtaining surface data.Therefore, the combination groups of these vegetation and topographic indicators are a good description of the soil properties.The optimal SBD and SOC in the different soil depth estimation models were selected to achieve the point-to-region expansion based on the cross-validation.The precision evaluation indicators of the SSA-RFR models are listed in Table 3.The simulated SBD and SOC exhibited strong correlations with the measured SBD and SOC across each model (Figure 8).The estimation accuracy of the SBD prediction model was significantly better than that of the SOC prediction model.Moreover, all of the SSA-RFR models have strong stability and applicability in SBD and SOC estimation.We validated the results of the SBD and SOC in different years using datasets obtained from the relevant studies and the HWSD.Our results have excellent consistency with the validation data.

Spatial Distribution of SOCD in 2020
Using the nine spatial explanatory variables and the optimized SSA-RFR models for the different soil layers (within 0-30 cm), we generated a map of the SOCD in 2020 (Figure 9).The SOCD was higher in the northern and western parts of the study area.In addition, the wetlands had the highest SOCD in the WSP, and the magnitude of the decrease in the SOCD with increasing soil depth was different.Moreover, the response of the SOCD to LUC varies with soil depth.The SOCD in built-up land exhibited the most pronounced reduction as soil depth increased due to the intense impact of human activities.The SBD and SOC content in the 0-5 cm soil layer were strongly influenced by land management patterns.These activities frequently resulted in soil quality degradation and an increased loss of SOC.In contrast, the shrub, barren land, grassland, and forest exhibited a higher soil quality due to the high vegetation cover and better sequestration conditions.The wetlands showed the highest stable SOCD trend due to the higher inputs and fixation of SOC.Combined with the land use data for the WSP, we calculated the SOCS of the different land use types at different soil depths (Table 4).The grassland had the highest SOCS in the WSP.The total SOCS in the study area in 2020 was 232.94 Tg C. The total SOCSs of the seven land use types in topsoil were different.The SOCS of the forest, shrub, and grassland remained stable and can be fixed in the soil for long periods.

Characteristics of Temporal-Spatial Distribution of SOCD
The average SOCD values were 7.21, 7.06, 6.88, and 7.82 t hm −2 in 1990, 2000, 2010, and 2020, respectively (Figure 10).The total SOCD increased by 0.61 t hm −2 from 1990 to 2020.The greatest increasing trend of the SOCD occurred in the northern and western regions of the study area, which were mostly grassland and forest areas.The total SOCD decreased from 1990 to 2010 and increased from 2010 to 2020.

Effects of LUC on SOCS during 1990-2020
We analyzed the spatial distribution of the LUC in the study area from 1990 to 2020.Considering that forests and grasslands were the main sources of the SOCS in the study area, we examined the carbon dynamics and LUC associated with the conversions between forests and grasslands (Figure 11).By comparing the changes in the SOCS in the different land use types, we aimed to better understand the implications of LUC on soil carbon dynamics.The land use types were predominantly grasslands and forests in the WSP.The changes in the built-up land and barren land were concentrated in the northeastern and western parts of the study area, respectively.Different land management patterns can significantly affect the SOC content of the soil at different depths.From 1990 to 2000, there was a significant reduction in wetland, with 367.84 km 2 being converted into grassland mainly due to climate warming and human activities.Between 2000 and 2010, there was a notable expansion in forests, increasing by 1.94 × 10 3 km 2 .In addition, from 2010 to 2020, there was a considerable increase in both wetlands and forests, with wetlands expanding by 238.60 km 2 and forests by 1.11 × 10 3 km 2 .

Impact of Land Management Practices on SOCS
Recently, significant changes in land use patterns have occurred in the WSP, including rapid urbanization [48] and land degradation [49].These changes have had a notable impact on the SOCS in different land use types.Forests and grasslands have played a role in enhancing carbon accumulation (Figure 10, Table 5), which is consistent with previous studies showing the positive effects of these natural ecosystems on SOC.Ploughing and agricultural activities, on the other hand, lead to a decrease in SOC.The total SOCS in the study area has significantly improved since 2010, primarily due to the conversion of grasslands to forests and wetlands.Forests and grasslands play a vital role in global carbon cycling and climate change mitigation, with significant potential for carbon sequestration and emission reduction.This conclusion is also supported by previous research [50][51][52].In this study, we confirmed that land use transitions can increase the SOCS in the WSP.Wetlands contain a large number of microorganisms that enhance the carbon sequestration capacity of these ecosystems.Therefore, it is imperative for policymakers to adapt current management practices and policies to preserve and restore degraded ecosystems in order to safeguard the stability and sustainability of SOCS in the study area.Additionally, this study analyzed the variations in the SOCS among different land use types and soil depths.The findings provide data and theoretical support for enhancing the SOCS in the WSP.Further analysis will be conducted to investigate the effects of factors such as fertilizers and irrigation on the temporal and spatial distribution of the SOCS at different soil depths.

Uncertainties and Limitations
The results of this study were affected by the variance in the historical image data, primarily due to the inconsistency in the acquisition times of the images for the different years.This can have an impact on the vegetation indices, such as the NDVI and EVI, and consequently on the SOCS prediction results.These vegetation indices are sensitive indicators of climate change and land use change.By employing vegetation indices' data with higher accuracy, the estimation of the SOCS can be significantly refined [53,54].Moreover, the temporal-spatial distribution of the SOCS in the study area was synergistically affected by several factors.In this study, we only explored the effects of land use change on the SOCS from qualitative and quantitative perspectives.Further research is still needed to address the synergistic effects of climate change and LUC on the SOCS in the study area.The results of this study showed that there were few data for the regions with higher elevations (Figure 1).In addition, the SOC and SBD data were obtained without considering the specific areas where the land use types have changed.Therefore, if more representative data can be obtained, we can more accurately explore the impact of LUC on the SOCS.We will further investigate future climate and land use under different scenarios to map the distribution and evolution of the SOCS in the study area.
To verify whether the results of the long-term SOCS are reliable, we compared several measured datasets from previous studies and the HWSD in different years with our simulation results (Figure 6).Strong consistency was found between the simulated data and measured data.Moreover, our results confirmed that the average SOCD in the topsoil in the WSP was 1.96 t m −2 in 2020, which is lower than in Chen's study [36].This difference can be attributed to our consideration of the effect of gravel content and multiple land use types on the SOCD.Our research takes into account more natural limitations and provides more accurate data support.

Conclusions
In this study, we explored the temporal-spatial changes in the LUC on the SOCS in the Western Sichuan Plateau with 30 m spatial resolution.SSA-RFR models and multisource datasets were utilized to estimate the distribution of the SOCS from 1990 to 2020.The vertical distribution of the SOCS in the study area decreased with increasing soil depth, and the SOCS increased from southeast to northwest.The wetlands had the highest SOCD.The grassland, as the main land use type in the study area, had the highest total SOCS (137.32 Tg C in 2020).Forests and grasslands were the main sources of the SOCS in the study area and had a significant impact on the increase in the total SOCS.More importantly, land use transitions were confirmed to increase the SOCS in the Western Sichuan Plateau under climate change.The results of this study confirmed that the SOCS in the study area increased significantly since 2010, with an increase of 27.88 Tg C over the total SOCS in 2010.Further, this study highlighted the importance of land management goals and practices in reaching the dual carbon goals early under the background of climate change.Because of the change, the SOCS is also likely to be affected by land management practices including fertilization and cultivation practices.It is crucial to understand the role of land management processes and to discern the effects of management practices to complement the SOCS change estimates provided in this study.

Figure 1 .
Figure 1.(a) The location of the study area.(b) The monthly mean precipitation and temperature data for the study area in 2020.(c) The land use data in 2020 and the spatial distribution of the SOC sampling sites.(d) The map of the soil types in the study area.

Figure 2 .
Figure 2. Schematic diagram of sampling area.(a) Distribution of samples in each sampling area.(b) Schematic diagram of soil profile.

Figure 3 .
Figure 3. Workflow chart of this study.

Figure 4 .
Figure 4. SOC content in different soil depths.

Figure 5 .
Figure 5.The linear relationships between the SBD and the SOC (a), pH, and SOC (b).

Figure 6 .
Figure 6.The SOCD of different land use types in the 0-30 cm interval.The white dots indicate the median and the upper and lower lines denote the maximum and minimum values of the SOCD.

Figure 7 .
Figure 7.The importance of indicators in the SSA-RFR models for estimating the (a) SBD and (b) SOC in the different soil layers.

Figure 8 .
Figure 8. Validation of simulated data using measured data.(a) Validation of SBD.(b) Validation of SOC.

Figure 9 .
Figure 9.The spatial distribution of the SOCD at different soil depths in the WSP obtained using the SSA-RFR models.

Figure 11 .
Figure 11.LUC in the study area from 1990 to 2020.We analyzed the time series changes in the SOCS from 1990 to 2020 (Table5) based on the land use changes.The forest and grassland were the main sources of the SOCS.The total SOCS exhibited a decreasing trend during 1990-2010, and the total SOCS increased from 2010 to 2020.The total SOCS increased by 18.03 Tg C from 1990 to 2020 and by 27.88 Tg C from 2010 to 2020.The grassland was the largest region that contributed to the SOCS in the study area during 1990-2020.The SOCS of the cropland was influenced by agricultural management practices and exhibited a decreasing trend.The SOCS of the built-up land with a vegetation region was affected by anthropogenic factors and exhibited a clear increasing trend.Land use change affects the activity and diversity of soil microorganisms and influences the stability and cycling processes of SOC in the cropland is the most affected by human activities, as shown by the fact that fertilization patterns during farming could accelerate the rate of decomposition of SOC.

Table 1 .
Relevant measured data of SOCS in China.

Table 2 .
Selected vegetation indices and calculation formulae.

Table 3 .
SSA-RFR model evaluation indicators for SBD and SOC prediction in different soil layers.

Table 4 .
Statistics of the SOCS in topsoil for the different land use types.

Table 5 .
SOCSs in the 0-30 cm soil layers for the different land use types.