The localities of Achupallas, Columbe, and Quimiag were selected because, collectively, they encompass the natural and transformed land-use systems evaluated in this study, including páramo, native forest, cropland, pasture, and exotic forest plantations. This selection allowed topsoil SOC stocks to be compared among land-use categories across three geographically distinct sectors of Chimborazo Province. Specifically, the localities were selected as follows: Achupallas (UTM coordinates X: 762669; Y: 9738743) in the south, Columbe (UTM coordinates X: 746154; Y: 9801254) in the central sector, and Quimiag (UTM coordinates X: 776275; Y: 9814701) in the northeastern part of Chimborazo Province (Figure 1). These localities predominantly feature Andosol soils of volcanic origin, characterized by low bulk density and high organic matter stabilization capacity (IUSS Working Group WRB, 2022), with textures ranging from loam to sandy loam and clay loam. This broadly comparable pedogenic context supported the evaluation of land-use-related differences in SOC stocks among the three localities.

In general terms, the main differentiating characteristics are as follows: Achupallas covers 101,600 ha, ranges from 2,000 to 4,440 m a.s.l., with a mean temperature of 12.5 °C, precipitation from 14 to 213 mm, and a relative humidity of 73%, and contains larger areas of endemic ecosystems such as forest and paramo.

Quimiag covers 13,949 ha, ranges from 2,400 to 5,319 m a.s.l., with a mean temperature of 16 °C, precipitation between 500 and 2,000 mm, and relative humidity of 79%. It includes the highest elevation (El Altar volcano) and the lowest population density of the three localities. It has two types of endemic ecosystems (forest and paramo) that are being moderately displaced by croplands and pastures.

Columbe covers 41,548 ha, ranges from 3,080 to 4,320 m a.s.l., with a mean temperature of 11 °C, precipitation from 14 to 1,750 mm, and a relative humidity of 73%. In contrast to the others, Columbe has the highest population density and only one natural ecosystem, the paramo (Figure 1).

In Ecuador, land uses in the high-Andean belt are clearly related to altitude and the environmental characteristics these ecosystems offer. The following systems were considered in this study (Figures 2, 3). Representative field photographs of these land-use systems, including natural ecosystems, agricultural areas, and burning practices associated with páramo conversion, are presented in Figure 2.

Natural paramos. The Ecuadorian paramo typically occurs between 3,000 and 4,800 m a.s.l., and is one of the most biodiverse ecosystems in the region, with species such as Calamagrostis intermedia and Chuquiraga jussieui, which are dominant and structurally important for the production of underground biomass, carbon capture, water regulation by their dense root system, and in the conservation of species adapted to extreme conditions, such as low temperatures and high humidity (Figures 3A,B).

The study followed an observational design, employing stratified sampling by land-use type and locality. The unequal number of samples among land-use categories reflected the actual field conditions in the high-Andean study areas, including the spatial distribution of land uses, access difficulties, topographic and logistical constraints, the time available for sampling, laboratory processing capacity, and available resources. Despite these constraints, the sampling design included the main land-use categories present across the three localities and allowed their SOC stocks to be compared within the scope of this study.

A 500 × 500 m grid was defined using ArcGIS 10.5. This grid was overlaid onto the map and aligned to the north for each locality, according to the stratification criteria (land use and location). Within each cell a 250 × 250 m sampling block was established. In each sampling block, three 20 × 20 m plots were set. Soil samples were collected at a depth of 0–30 cm, as this layer is particularly responsive to LULC changes and provides a standardized depth for comparison among land-use categories, thereby directly reflecting alterations driven by anthropogenic activities (Mahmood et al., 2024). Sampling was conducted using a soil auger, collecting five subsamples of 0.5 kg each from the four corners and the center of each plot. The subsamples were combined in a container, thoroughly homogenized, and used to obtain a composite soil sample weighing approximately 1 kg (Figure 4).

Each sample was labeled with UTM coordinates (X and Y), sampling date, locality, land-use type, and identification number, and was immediately transported to the laboratory for storage at room temperature and subsequent analysis. Soil sampling was carried out between October 2023 and May 2024, corresponding to the rainy season in Ecuador. Across the three localities, Achupallas, Columbe, and Quimiag, a total of 57 paramo sites, 27 natural forest sites, 30 pastures, 21 croplands, and 9 exotic forest plantations were sampled (Table 1).

Soil samples were collected from the 0–30 cm layer to standardize comparisons among land-use categories and to maintain consistency with the default 30 cm depth used in the IPCC Tier 1 approach for estimating SOC stocks in mineral soils. This depth is widely used in SOC studies because SOC stocks in the surface mineral layer can be particularly responsive to changes in land use and soil management (IPCC et al., 2019; Lim et al., 2018). Consequently, the reported SOC stocks refer exclusively to the 0–30 cm layer and do not represent total SOC stocks throughout the soil profile.

Soil bulk density (SBD) was determined from one undisturbed core sample collected from each plot (total of 144 samples) at a depth of 0–30 cm, corresponding to the same plots used for SOC sampling. Using a 100 cm³ stainless steel cylinder, each cylinder was labeled, sealed with plastic film to prevent moisture loss, placed in a resealable polyethylene bag (Ziploc), and transported to the laboratory in portable coolers (FAO, 2023). Upon arrival at the laboratory, each sample was removed from the cylinder and weighed to obtain the wet weight, then oven-dried at 105 °C for 24 h to determine the dry weight (DW). Soil bulk density was calculated as: SBD = D W   g V   c m 3

In addition, soil samples were air-dried for 3–6 days depending on their initial moisture, and sieved to 2 mm. The pH was determined on these samples using a soil-to-water ratio of 1:2.5 (w/v) (FAO, 2021).

The particle-size distribution was determined using the Bouyoucos hydrometer method on 100 g of oven-dry soil, employing a dispersing solution to ensure complete breakdown of soil aggregates and to quantify the proportions of sand, silt, and clay (Gee and Bauder, 2018; Motsara and Roy, 2008).

Soil organic matter (SOM) was determined by loss-on-ignition on ground material (sieved at 425 and 212 µm). Five grams from the 212–425 µm fraction were weighed into pre-tared 30 mL crucibles, dried at 105 °C for 24 h, weighed ( W 105 ℃ , then ashed at 450 °C for 2 h and weighed again ( W 450 ℃ (Martínez et al., 2017; Nelson and Sommers, 1996; Hoogsteen et al., 2018). Prior to weighing, samples were placed in a desiccator to reach room temperature and avoid moisture uptake; the mass loss corresponds to SOM. SOM  % = W 105 ℃   g ‐ W 450 ℃   g   W 105 ℃   g   *   100

The loss-on-ignition (LOI) method was selected due to its low operating cost for a large number of analyses, making it appropriate for large-scale soil evaluations. The LOI method has been widely used to estimate soil organic matter (SOM) and infer SOC dynamics and losses across different soil types and ecosystems (Chambers et al., 2024; Hoogsteen et al., 2018; Jensen et al., 2018; Kupka and Gruba, 2022). Previous studies comparing the LOI and DUMAS (dry combustion) methods have reported no significant differences in their results (Ng et al., 2025; Cargua Catagña et al., 2017). However, the Van Bemmelen conversion factor (1.724) is commonly used to convert SOM to SOC, although this factor may vary depending on soil mineralogy (Pribyl, 2010; Allison, 2016; Minasny et al., 2020). Since most soils in the central highlands are of volcanic origin (Andosols), the LOI values should be interpreted as estimates, not absolute measurements. SOC  % = S O M   *     1 1 , 724

Statistical analyses were conducted in the R environment (v. 4.3.0) (R Core Team, 2023), using the tidyverse and ggplot2 packages for data processing and visualization (Wickham et al., 2019; Wickham, 2016). SOC stock differences were evaluated using linear models (one-way ANOVA) fitted with the lm function, with land use treated as a fixed categorical factor.

Given the observational nature of the study and the spatial distribution of samples, analyses were performed on a priori defined data subsets to assess SOC variability at different scales: (i) between natural and managed ecosystems, (ii) among land-use categories, and (iii) within each locality (Achupallas, Columbe and Quimiag), using independent models.

This approach corresponds to a stratified (multi-scale) analysis and allows the characterization of SOC variability across analytical scales without fitting a single global model.

Model assumptions were evaluated using residual diagnostics (Q–Q plots and residuals vs. fitted values). When necessary, data transformations were applied to improve model fit. Multiple comparisons were performed using estimated marginal means with Tukey adjustment (α = 0.05), implemented in the emmeans package (Lenth, 2023). Because design-based sampling weights or inclusion probabilities were not available, unweighted models were used. The fitted linear-model framework allowed comparisons among categories with unequal sample sizes, while comparisons involving the least represented categories were interpreted cautiously.

Residual spatial autocorrelation was evaluated separately for each locality using Moran’s I and a row-standardized spatial weights matrix based on the four nearest neighbors derived from the UTM coordinates. Significance was assessed using 9,999 permutations, and sensitivity was examined using three to six nearest neighbors. Residual semivariograms were also examined. Because significant residual spatial dependence was detected only in Columbe, its locality-specific model was refitted using generalized least squares (GLS) with an exponential spatial correlation structure based on Euclidean distances among sampling locations. For this spatial model, SOC stock values were square-root transformed, and pairwise comparisons were adjusted using the Holm method. Moran’s I was subsequently recalculated on the normalized residuals of the adjusted model. The original linear models were retained for Achupallas and Quimiag. These analyses were performed using the spdep, gstat, and nlme packages.

The study areas exhibit distinct geological characteristics; Columbe and Quimiag possess relatively uniform andesitic and tuffaceous compositions, whereas Achupallas displays greater lithological complexity, including metamorphic and sedimentary formations (Table 2).

The pH showed differences by use and locality, ranging between 5.39 ± 0.54 and 6.83 ± 0.07. In general, natural ecosystems exhibited lower pH values compared to intervened systems (p < 0.05).

The SBD exhibited an inverse relationship with Organic Matter (OM) content across all study areas (Table 2). Natural ecosystems, which recorded the highest OM concentrations (reaching up to 45.39 ± 14.23 in Achupallas forests), consistently showed the lowest SBD (0.29 ± 0.06–0.75 ± 0.10 g cm⁻³). In contrast, intervened areas were characterized by higher SBD (0.57–1.31 g cm⁻³) and generally lower OM levels.

Soil texture differed clearly among localities: sandy loam predominates in Columbe and Achupallas, whereas clay loam prevails in Quimiag. Coarser textures in Columbe and Achupallas may facilitate aeration but limit water retention, contrasting with the finer texture of Quimiag that favours water retention.

SOC stock variation was examined at three complementary scales of interpretation: (i) ecosystem level (natural vs. intervened systems), (ii) land-use category, and (iii) land use within each locality.

The SOC stock exhibited pronounced differences between natural ecosystems and intervened systems across the entire study area. On average, natural ecosystems (paramo and natural forests) displayed significantly higher SOC stocks (264.1 ± 101.3 Mg C ha⁻¹) than intervened systems (176.0 ± 82.55 Mg C ha⁻¹; p < 0.001) (Figure 5A), indicating a substantial depletion of soil carbon associated with land-use change.

Statistical analysis demonstrated significant overall differences in SOC stocks between natural and anthropogenic systems. While natural forests and paramo ecosystems exhibited the highest SOC stocks with no significant differences between them, land-use conversion to productive systems resulted in a significant reduction of the carbon pool (p < 0.001). Pastures showed the most marked decline in SOC stocks and presented significantly lower values than natural ecosystems, whereas croplands and exotic forest plantations showed intermediate values and did not differ significantly from either natural ecosystems or pastures (Figure 5B).

SOC stocks varied among localities, reflecting the influence of spatial context and land-use type (Figure 6). Residual spatial autocorrelation was not significant in Achupallas (Moran’s I = −0.106, p = 0.446) or Quimiag (I = 0.021, p = 0.600). In Columbe, the initial model showed weak but significant positive residual spatial dependence (I = 0.209, p = 0.008). After fitting the spatial GLS model, residual autocorrelation was no longer significant (I = 0.031, p = 0.549). Sensitivity analyses using three to six nearest neighbors and the examination of residual semivariograms supported these results. Accordingly, the original linear models were retained for Achupallas and Quimiag, whereas the spatial GLS model was used for Columbe. In Achupallas, the highest SOC stocks were observed in exotic forest plantations (277.3 ± 15.80 Mg C ha⁻¹), followed by paramo (243.1 ± 84.34 Mg C ha⁻¹) and natural forests (220.5 ± 24.60 Mg C ha⁻¹), with no significant differences among these systems. In contrast, grasslands exhibited the lowest mean SOC stock (182.6 ± 31.90 Mg C ha⁻¹), indicating a lower SOC storage under this land use. Croplands were not present at this site. In Quimiag, natural ecosystems accumulated the highest SOC stocks, with particularly elevated values in paramo (396.5 ± 138.2 Mg C ha⁻¹) and natural forests (340.4 ± 102.7 Mg C ha⁻¹), with no significant differences between them. Páramo soils exhibited significantly higher SOC stocks than croplands (253.9 ± 15.0 Mg C ha⁻¹) and grasslands (249.3 ± 65.2 Mg C ha⁻¹), whereas natural forests showed intermediate values and did not differ significantly from either páramo or the anthropogenic land uses (Figure 6). Finally, in Columbe, where paramo constitutes the dominant natural ecosystem, this system exhibited a markedly higher SOC stock (248.0 ± 92.8 Mg C ha⁻¹) than intervened land uses. Pastures (116.4 ± 48.9 Mg C ha⁻¹), croplands (110.8 ± 33.3 Mg C ha⁻¹), and particularly exotic forest plantations (56.9 ± 14.2 Mg C ha⁻¹) showed pronounced reductions in soil carbon storage according to the spatial GLS model (all Holm-adjusted p-values <0.001), providing clear evidence of substantial SOC losses associated with paramo conversion.

The differences observed in SOC stocks between natural ecosystems and intervened systems resulted in substantial contrasts in potential CO₂-equivalent emissions associated with land-use change. When comparing the mean SOC values obtained for natural ecosystems (paramo and natural forests) and intervened systems (grasslands, croplands, and exotic forest plantations), an average SOC loss of 88.1 Mg C ha⁻¹ was estimated following land-use conversion (Table 3). This loss corresponds to a potential CO₂-equivalent emission estimate of approximately 323 Mg CO₂-eq ha⁻¹.

In addition, the sensitivity analysis showed that using conversion factors of 1.9 and 2.0 reduced the estimated mean SOC loss from 88.1 Mg C ha⁻¹–80.0 Mg C ha⁻¹ (−9.3%) and 76.0 Mg C ha⁻¹ (−13.8%), respectively. The corresponding estimates of potential CO₂-equivalent emissions decreased from 323.1 Mg CO₂-eq ha⁻¹ to 293.2 and 278.5 Mg CO₂-eq ha⁻¹, respectively. Despite these reductions, the overall pattern was still evident, as natural ecosystems continued to have higher SOC stocks than transformed land-use systems.

To incorporate the temporal dimension of land-use change, multitemporal land-cover assessments between 2000 and 2020 were used to estimate an annualized soil carbon loss. By combining field-measured SOC differences and the transformed areas in each locality derived from the LULC analysis (Damián-Carrión et al., 2026), an average annualized SOC loss estimate of 4.4 Mg C ha⁻¹ yr⁻¹ was estimated. This value corresponds to an annualized potential CO₂-equivalent emission estimate of approximately 16.15 Mg CO₂-eq ha⁻¹ yr⁻¹ per hectare converted during the 20-year period.

At the local scale, annual SOC losses and estimated potential CO₂-equivalent emissions varied considerably among locations. This heterogeneity reflects both the specific magnitude of carbon depletion and the extent of ecosystem transformation identified across the different study locations (Table 4). In Achupallas, where the conversion of páramo and forest was most extensive, SOC losses reached 70.0 Gg C year⁻¹, equivalent to 256.7 Gg CO₂-eq year⁻¹. In Columbe, the transformation of the moorland generated estimated emissions of 193.0 Gg CO₂-eq year⁻¹ (52.6 Gg C year⁻¹), while in Quimiag, the associated emissions were considerably lower, reaching only 5.28 Gg CO₂-eq year⁻¹ (1.44 Gg C year⁻¹), consistent with the smaller area affected during the 20-year study period.

At the regional scale, paramo ecosystems in Chimborazo Province cover approximately 196,053 ha, as estimated by the Ecuadorian governmental environmental authority for 2023, collectively storing around 51,299 Gg C as SOC. Under a theoretical scenario involving the complete conversion of these natural ecosystems to anthropogenic land uses (grasslands, croplands, and exotic forest plantations), potential emissions associated with land-use change would reach approximately 9,405 Gg CO₂-eq yr⁻¹, representing a substantial fraction of Ecuador’s annual national emissions. These estimates correspond to potential emissions derived from SOC stock losses and do not represent directly measured atmospheric emissions.

In summary, the observed reductions in SOC stocks across high-Andean ecosystems represent potential CO₂-equivalent emissions. These values varied according to land-use type, as well as the intensity and spatial extent of ecosystem conversion in each study area.

Converting natural ecosystems to agricultural land significantly alters key soil parameters. Paramo and forest soils had more acidic pH than intervened systems, as a consequence of the alkalinizing effect of ash from pre-cultivation burning (Arunrat et al., 2024; Elakiya et al., 2023; McLaughlin, 2009; Ratier Backes et al., 2021) (Figure 2D) and other inputs applied by farmers (e.g., CaCO₃, CaO₂) (Enesi et al., 2023; Olego et al., 2021; Torres Ramos et al., 2014; Wenyika et al., 2025). Increases in pH affect the soil’s capacity to store C by modifying stabilization mechanisms of organic matter and favouring microbial activity (Malik et al., 2018). In general, the transformation of acidic soils to near-neutral pH reduces the solubility of Fe, Al, and Mn and the reactivity of their hydroxides, and therefore the stability of organo-mineral complexes that protect soil organic matter (Hargrove, 1986; Wen et al., 2025), while increasing the bioavailability of nutrients such as P and N. These soil-level changes favour the mineralization of organic matter by microorganisms and the consequent loss of SOC stock (Wang and Kuzyakov, 2024).

Similarly, the higher compaction observed in intervened soils (SBD up to 1.3 g cm⁻³) relative to paramo and native forest soils (SBD: 0.29–0.75 g cm⁻³) indicates loss of porosity and infiltration as a consequence of SOC loss in the Andean region (FAO y ITPS, 2015; Patiño et al., 2021) and specifically in Ecuadorian paramos (Beltran-Davalos et al., 2022; Comas et al., 2017; Paez-Bimos et al., 2022). This structural collapse is mediated by a robust inverse coupling between SBD and OM, wherein the depletion of the organic matrix facilitates the reorganization and closer packing of mineral particles. Given the common pedogenic origin of these sites, these differences are primarily attributed to land-use intensification and altitude-driven climatic conditions rather than inherent textural variations (Angst et al., 2023; Jungkunst et al., 2022).

Soil texture is a key edaphic factor influencing the physical protection and stability of SOC (Six et al., 2002). Although soil texture does not change with land-use conversion, it regulates organic matter aggregate formation and the adsorption of organo-mineral particles. This physical restriction decreases microbial mineralization and decomposition rates, thereby increasing SOC residence time and enhancing SOC stability. In this study, the prevalence of loam textures may contribute to SOC stabilization through physical protection mechanisms (Yao et al., 2023). However, variations in SOC from a lithological perspective can be explained by differences in parent material. In Quimiag, soils have developed over pyroclastic materials rich in reactive minerals and exhibiting clay loam textures, which promote more effective SOC stabilization. In contrast, Achupallas is characterized by soils derived from agglomerates, siliceous lavas, and metamorphic materials, with sandy loam textures that provide lower mineral protection for SOC. Finally, Columbe originates from a combination of andesitic lavas, tuffs, and shales, resulting predominantly in sandy loam textures and intermediate conditions, where land use plays a decisive role in SOC conservation or loss (Torn et al., 1997; Yang et al., 2020).

Previous studies have highlighted the high SOC content of high-Andean soils (Beltran-Davalos et al., 2025) (Table 5), considering them among the most important C reservoirs on the planet (Beltran-Davalos et al., 2022; Fernandez Pérez et al., 2019). Volcanic-ash-derived soils characteristic of the region show a notable capacity to accumulate and store organic carbon, contributing substantially to the global balance of soil carbon (García et al., 2019; Tonneijck et al., 2010). A comparison of our results with the values presented in Table 5 shows that the mean SOC stock obtained in this study (264.1 Mg C ha⁻¹ at 0–30 cm) was higher than those reported for the páramos of Peru and Colombia, as well as one of the values recorded in Bolivia. The second value for Bolivia (260 Mg C ha⁻¹) was similar to ours, which was also close to the upper limit of the range documented for Ecuador (125–275 Mg C ha⁻¹). Therefore, the result obtained is among the highest values recorded for the surface layers of the high-Andean soils included in Table 5, highlighting the importance of the páramos and native forests of the Chimborazo province as carbon reservoirs. However, the accumulated stock of SOC tends to increase with depth; this allows us to contextualize the value of 530 Mg C ha⁻¹ reported for northern Ecuador in profiles of 150–200 cm depth, compared to the first 30 cm considered in this study. SOC stock in paramo soils can be twelve times higher than in forest plantations (Tonneijck et al., 2010) such as pine plantations (Pinus patula, Pinus radiata; 12.5–14.8 Mg C ha⁻¹) (Quiroz Dahik et al., 2021) or eucalyptus plantations (Eucalyptus sp.; 18.2 Mg C ha⁻¹) (Carrión-Paladines et al., 2022), and even higher than in tropical ecosystems such as lowland tropical rainforest (100–200 Mg C ha⁻¹) and tropical montane forest (∼74 Mg C ha⁻¹) (Suarez, 2012). Our results show that SOC stock in Andean paramo soils exhibits high spatial variability (Zhu et al., 2019) exceeding 50% among localities (Figure 6). Despite this spatial heterogeneity, a consistent pattern emerges in which natural ecosystems maintain substantially higher SOC stocks than anthropogenic land uses, as observed in Figure 5. These results support our first hypothesis that natural ecosystems (paramo and native forests) store significantly higher SOC stocks than anthropogenic land uses such as croplands, pastures, and exotic forest plantations. The high SOC stocks in paramo ecosystems are primarily driven by the combination of low temperatures and high soil moisture, which inhibit microbial decomposition rates (Alavi-Murillo et al., 2022). These environmental constraints, coupled with the high recalcitrance of Andean vegetation, promote a greater accumulation of organic matter in paramo soils than in forest areas.

Furthermore, the results obtained show that it is a fragile reservoir and that substantial SOC losses may occur when natural soil is transformed into agricultural soil (Table 4). This pattern suggests that land-use conversion and intensification reduce the mechanisms that protect SOC in volcanic soils, particularly through disturbance of soil structure and loss of organic matter protection (Staß et al., 2025; Kumar et al., 2024). Therefore, the observed differences among land-use categories emphasize the sensitivity of these SOC-rich soils to anthropogenic disturbance.

A comparative analysis of SOC stocks (Table 4) reveals a consistent pattern of lower stocks in agricultural and livestock systems and exotic forest plantations than in natural Andean ecosystems. Despite methodological variations in sampling depth among the cited studies, our results show an estimated mean SOC difference of 88.1 Mg C ha⁻¹, which exceeds the previously documented range of 12.9–57 Mg C ha⁻¹. The greater difference is mainly related to the high SOC stocks of the natural reference ecosystems, paramos and native forests (264.1 Mg C ha⁻¹), rather than to exceptionally low values in the transformed land uses (176.0 Mg C ha⁻¹). In addition, our analysis integrates two natural ecosystems and three transformed land-use categories, whereas previous studies generally evaluated specific land-use conversions. Therefore, the 88.1 Mg C ha⁻¹ represents a mean difference among current land-use categories rather than a temporal loss rate measured directly. The magnitude of this estimated SOC difference not only reflects a decline in soil quality but also provides the basis for estimating the associated potential CO₂-equivalent emissions. These estimates do not represent directly measured CO₂ fluxes.

The use of a fixed Van Bemmelen factor introduces uncertainty into SOC stock estimates and potential CO₂-equivalent estimates. The sensitivity analysis showed that conversion factors of 1.9 and 2.0 reduced these estimates by 9.3% and 13.8%, respectively, compared with the baseline factor of 1.724. This suggests that the baseline estimates may be moderately overestimated if the actual conversion factor for these soils is closer to 1.9 or 2.0. However, Armijos-Arcos et al. (2025) reported a factor of 1.729 for Ecuadorian páramo soils, which is very close to that used in this study; therefore, the actual magnitude of the bias cannot be determined without an independent measurement of SOC.

In addition, the reported estimates should be interpreted considering that high-Andean soils may store substantial amounts of SOC below the depth evaluated in this study (0–30 cm). Although this standardized depth allowed consistent comparisons among land-use categories and the assessment of a soil layer that is particularly responsive to changes in land use and management, SOC stocks in deep organic-rich profiles and peat-forming environments, where substantial carbon accumulation may occur below 30 cm, were not assessed (Tonneijck et al., 2010; Comas et al., 2017). Therefore, the SOC stocks reported here refer exclusively to the 0–30 cm layer and should not be extrapolated to the entire soil profile.

The potential CO₂-equivalent emissions estimated in this study were derived from SOC stock differences associated with land-use change, following the stock-difference approach recommended by the IPCC guidelines for the AFOLU sector. This relationship supports our second hypothesis that lower SOC stocks in transformed land uses, relative to natural reference ecosystems, correspond to substantial estimated SOC losses and associated potential CO₂-equivalent emissions. Under this framework, decreases in SOC stocks can be expressed as a potential net carbon release to the atmosphere, an approach that is widely applied in regional and global assessments when direct measurements of soil CO₂ fluxes are not available (IPCC et al., 2019). However, as these estimates rely on a “space-for-time” substitution rather than continuous longitudinal monitoring, they represent potential rather than observed fluxes; consequently, our results should be interpreted as an approximation of long-term landscape transformation impacts rather than a record of transient soil carbon dynamics. Accordingly, the observed differences in SOC stocks support the interpretation that land-use conversion contributes to SOC losses, although the observational design does not allow direct causal attribution.

The average SOC loss of 88.1 Mg C ha⁻¹ (equivalent to 323 Mg CO₂-eq ha^-1) observed following the conversion of paramo and natural forests to anthropogenic land uses falls within the upper range reported by global land-use change studies, which document reductions of 20%–50% of the original SOC stocks, particularly in tropical and mountainous regions characterized by high initial soil carbon contents (Beillouin et al., 2023; Don et al., 2011; Guo and Gifford, 2002). In the context of high-Andean volcanic soils, which are characterized by substantial carbon accumulation under natural conditions (Alavi-Murillo et al., 2022; Staß et al., 2025), SOC losses can be particularly pronounced following extensive disturbances affecting soil depths of 50–200 cm (Fu et al., 2025), which is consistent with the magnitude of the potential emissions estimated in this study.

For comparison, previous studies have reported SOC loss rates for other páramo and Andean forest ecosystems. For instance, Thompson et al. (2021) reported SOC loss rates in intermediate paramo ecosystems of 0.045% C yr⁻¹, corresponding to <0.1 Mg C ha⁻¹ yr⁻¹, whereas Calderón-Loor et al. (2020) estimated losses ranging between 1.5 and 2.5 Mg C ha⁻¹ yr⁻¹ in abandoned forests and paramos. These values provide useful context for the annualized estimate obtained in this study; however, direct quantitative comparison should be interpreted with caution because our estimate was derived from SOC stock differences rather than from repeated temporal measurements. From a national perspective, the magnitude of SOC stock losses estimated here represents approximately 26% of the 36,088 Gg CO₂ emissions reported for Ecuador’s production system by Buenaño et al. (2023), and nearly half of Ecuador’s net greenhouse gas emissions in 2018 (16,280 Gg CO₂-eq), which were largely associated with land-use change and forestry activities (Meurer and Soria, 2024). The high emission rates estimated in this study, which exceed the values reported for lower altitude ecosystems, are probably due to the inherent vulnerability of high-altitude SOC stocks (Chen et al., 2024; Kuhry et al., 2022). While the natural accumulation of microbial residues in these cold environments creates organo-mineral complexes (stable micro- and macro-aggregates) that physically protect the carbon from mineralization, SOC may remain vulnerable to further losses (Hemingway et al., 2019; Luo et al., 2024). Intensive management can trigger greenhouse gas (GHG) emissions (Kopittke et al., 2024; Tobiloba et al., 2025) and the collapse of soil aggregates (Liang et al., 2025). This pattern may be associated with soil structural disruption; while the natural accumulation of microbial residues creates stable organo-mineral complexes that physically protect carbon (Han et al., 2025), intensive management disrupts these aggregates and exposes previously shielded organic matter to rapid mineralization. This process is further exacerbated by nutrient imbalances that constrain microbial carbon use efficiency, leading to higher CO₂ emissions despite potential carbon inputs (Wiesmeier et al., 2019).

The annualized SOC loss estimate for the analyzed 20-year period can be interpreted in light of evidence indicating that most SOC losses occur during the first decades following land-use conversion, before ecosystems reach a new quasi-equilibrium state after several decades or even more than a century (Sanderman et al., 2017; Smith, 2008). Building on this, our results reflect findings in high-Andean volcanic soils, where the rapid depletion of labile C fractions during the first three decades confirms this phase as a critical period of accelerated mineralization (Staß et al., 2025). Consequently, the annualized estimate should be interpreted cautiously in relation to rates reported for tropical and montane ecosystems undergoing early stages of conversion and is consistent with the high vulnerability of SOC in high-Andean ecosystems subjected to anthropogenic pressures (Poeplau et al., 2011; Zhao et al., 2025).

The regional-scale projections presented in this study are not a retrospective account of historical degradation; rather, they represent theoretical scenarios of potential emissions derived from the total or partial mobilization of current SOC stocks. While our baseline is informed by multitemporal assessments, these projections should not be interpreted as actual current emissions. Complete conversion of paramo ecosystems is unlikely due to legal, environmental, and socioeconomic constraints; however, these scenarios allow the quantification of the carbon storage ecosystem service provided by these systems and highlight their relevance for nature-based climate change mitigation strategies (Griscom et al., 2017; Rumpel and Chabbi, 2021).

Estimates of potential emissions are subject to uncertainties associated with soil spatial heterogeneity, temporal variability in SOC decomposition rates, and errors inherent to satellite-based land-use classification (Brasika et al., 2025). Previous studies indicate that classification errors using these approaches typically range between 5% and 15%, or are expressed as up to 95% classification confidence, introducing a comparable level of uncertainty in regional-scale emission estimates (García et al., 2019; Olofsson et al., 2014). Despite these limitations, the magnitude of SOC losses estimated in this study is consistent with ranges reported for land-use change in high-Andean ecosystems. Our results therefore support the view that the degradation and conversion of paramo and montane forest ecosystems may substantially reduce their capacity to function as long-term carbon reservoirs, potentially contributing to increased CO₂ emissions at regional scales (Carrillo-Rojas et al., 2019; Staß et al., 2025).

Land-use change was consistently associated with lower soil organic carbon reserves in the high Andean ecosystems of the central Ecuador highlands, confirming that paramos and native forests function as important SOC reservoirs. Their conversion to anthropogenic land uses, such as croplands, pastures, and exotic forest plantations, was associated with substantial SOC depletion, highlighting the high sensitivity of these ecosystems to disturbances and the crucial importance of their conservation for maintaining soil carbon storage.

SOC losses can be expressed as potential CO₂-equivalent emissions, indicating that the degradation of high-Andean ecosystems may reduce their capacity to retain soil carbon. While uncertainties remain regarding spatial soil heterogeneity and land-use classification, these results emphasize the important role of paramo and montane forest ecosystems as long-term SOC reservoirs and their strategic importance for climate change mitigation and sustainable soil management in the tropical Andes. Therefore, the protection and sustainable management of these ecosystems with high SOC content should be considered a priority within regional land-use planning strategies, particularly for climate change mitigation. Future research should incorporate repeated long-term measurements, deeper soil layers, and spatially explicit assessments to better quantify SOC dynamics and reduce uncertainty in potential CO₂-equivalent emission estimates.