the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
Asymmetric patterns of soil carbon mineralization during thaw slump recovery: Stable fast-pool size but rebounding decomposition rate
Abstract. In areas of thermokarst landslides on the Qinghai-Tibet Plateau, the decline in soil CO2 emissions over time is accompanied by a decrease in the proportion of active carbon and an increase in the proportion of inert carbon, suggesting that the soil organic carbon pool capacity and its mineralization rate may have changed. The lack of quantitative analysis regarding how actual active carbon pool capacity and mineralization rates regulate CO2 emissions has led to some uncertainty in CO2 emission estimates. In light of this, this study utilized a publicly available dataset of a 23-year recovery sequence from thermokarst landslides in the Qilian Mountains. By introducing a dual-pool index model, the cumulative CO2 release process was decoupled into fast-pool capacity (C0,fast) and fast-pool decomposition rate constant (kfast), and their driving factors were identified using hierarchical variance decomposition (LMG) and partial least squares structural equation modeling (PLS-SEM). The results indicate that over the 23-year recovery period, C0,fast did not change significantly, while the turnover rate kfast rebounded significantly. Furthermore, C0,fast was primarily driven by total microbial biomass (r = 0.74, p = 0.003), whereas kfast was influenced by the chemical quality of organic carbon (r = 0.60, p = 0.023). Both LMG and PLS-SEM confirmed that the chemical quality of organic carbon is the primary driver of variation in kfast (contributing approximately 30 %), while microbial indicators contributed less than 15 %. Furthermore, as the number of recovery years increased, δ¹³C exhibited a significant negative trend (β = −0.761, p < 0.001), and both this trend and the subsequent recovery in pH were significantly correlated with kfast (r = −0.53 and 0.51, respectively). Soil carbon mineralization during the recovery from thermokarst landslides exhibits an asymmetric pattern characterized by “stable pool capacity and a rebound in mineralization rate.” This finding provides new kinetic evidence for assessing carbon dynamics following permafrost disturbance.
- Preprint
(1124 KB) - Metadata XML
-
Supplement
(766 KB) - BibTeX
- EndNote
Status: final response (author comments only)
-
RC1: 'Comment on egusphere-2026-3656', Anonymous Referee #1, 04 Aug 2026
-
AC1: 'Reply on RC1', Hanhan Li, 11 Aug 2026
Responses to General Comments:
We sincerely thank the reviewers for their thorough and constructive review of this manuscript. We are very pleased that the reviewers recognized the high relevance of the topic, “Assessing the CO2 turnover time in thermal karst soils,” to this journal. We fully accept the criticisms regarding the manuscript’s current shortcomings in scientific novelty and quality, and we regard them as important guidance for refining this study. We believe these revisions will significantly enhance the academic contribution of this manuscript. We are committed to implementing these changes, and below we provide a point-by-point response to your valuable comments.
Responses to specific comments:
We would like to thank the reviewers for their critical and constructive comments, which have helped us clarify the interpretation of our model results and enhance the rigor of the paper. Below, we respond to each point in turn.
- Clarification Regarding the Term “C0,fast” and the Cumulative CO2 Release Paradox
Thank you for your comments; we apologize for any confusion this may have caused. First, we agree that more precise terminology is needed for C0,fast. C0,fast is a parameter derived from fitting the cumulative CO2 curve during the cultivation period; it represents the potential mineralizable capacity of the fast pool under controlled cultivation conditions, rather than a direct measurement of the actual fresh carbon pool in the field. We will revise the statement to read: “C0,fast did not undergo a significant change.” Second, regarding the apparent paradox of “why a seemingly stable fast carbon pool capacity is accompanied by a decrease in cumulative CO₂ release,” this should be understood in terms of the dynamic equilibrium of C0,fast, which reflects the depletion of readily degradable carbon early on and the subsequent continuous replenishment of carbon sources. The key point here is that these compensatory inputs are primarily dominated by chemically complex, recalcitrant substrates. As these substrates are gradually degraded by microorganisms, they release soluble metabolites or microbial residues with relatively fast turnover rates. In the dual-pool model fitting, these products are kinetically classified as fast-pool components; however, functionally, they differ from the readily available polysaccharides commonly found in soil prior to disturbance. Therefore, although C0,fast is maintained, the average activity and mineralization efficiency per unit of carbon have decreased.
- Missing Discussion of Slow-Catalyst Parameters
Thank you for your comments and suggestions. We will address the lack of significant differences in the slow-catalyst parameters (C0,slow , and kslow) in our revised manuscript.
- Explanation of the Apparent Contradiction with the Results of Mu et al. (2024)
We thank the reviewer for pointing out the potential contradiction between our findings and those of Mu et al. (2024), and we apologize for any confusion this may have caused. We would like to clarify here that the two studies are not contradictory but rather represent complementary interpretations of the same recovery sequence from different perspectives. Mu et al. (2024) used physical fractionation and static chemical analysis methods, such as 13C NMR, to reveal the compositional evolution of the soil carbon pool during the recovery of a thermal karst landslide: the content of particulate organic carbon (POC) decreased with the duration of recovery, while inert carbon fractions (such as mineral-associated organic carbon, MAOC, and alkyl carbon) gradually accumulated. Our study does not dispute this trend in the transformation of the static carbon pool. Instead, by performing kinetic decoupling of the cumulative CO2 release curve using a two-compartment exponential model, we provide two key insights from the perspective of carbon turnover dynamics. (i) The absence of significant changes in C0,fast indicates that, although the original readily degradable carbon was largely lost after the landslide as revealed by kinetic fitting, the fast pool did not disappear as a result. Instead, the composition of the fast pool underwent a qualitative shift during the recovery process, with the lost readily degradable carbon being replaced by new components. While these new components are indeed more “stubborn” than the original active carbon, they still exhibit a significantly faster turnover under cultivation conditions compared to the inert carbon in the slow pool; therefore, the model classifies them as the “fast pool.” (ii) The significant rebound in kfast reveals the recovery dynamics of the fast pool’s turnover rate. In the early post-landslide period (years 1-12), kfast remained suppressed; by year 23, kfast had recovered significantly, reaching or even slightly exceeding the control level. This does not imply that the decomposition rate has returned to its pre-landslide state, but rather indicates that microbial mineralization is gradually recovering, even though the substrate is no longer the original high-quality organic carbon but rather a substitute active carbon with more refractory chemical properties. Mu et al. (2024) addressed the questions of “what the carbon pool is and how it changes,” while we used C0,fast and kfast to further address the questions of “whether the carbon pool’s turnover rate has recovered and whether its capacity has changed accordingly.” We will integrate these perspectives into the Discussion section of the revised manuscript.
Responses to Technical comments:
Thank you for your technical suggestions. We have clarified and revised the text in response to your questions and recommendations, and these improvements have effectively strengthened the reliability of our statistical inferences. Below, we address each point individually.
- Clarification Regarding Missing Samples and Total Sample Size
We appreciate the reviewer’s comments and apologize for any confusion this may have caused. The original in vitro experimental dataset (Mu & Mu, 2025) included a total of 24 soil samples. However, one of these samples, 1-HT2-B2 (Year 12, subsoil, replicate 2), had no cumulative CO2 release records in the publicly available dataset; this sample may have been lost during the incubation process or was not measured. Consequently, this sample could not be included in any model fitting from the outset. This represents a limitation in the availability of the original research data itself, rather than a subjective decision made during our analysis. All analyses were based on the 23 samples with complete cumulative mineralization curves. Furthermore, among the 23 samples with valid data, one sample (1-HT3-A1, 23-year topsoil) had fewer than four valid observation points (data missing after day 69), making it impossible to reliably estimate the two-tank model (which requires at least four observation points for its four parameters). Therefore, this sample was fitted using only a single-tank model, with both the fast-tank and slow-tank parameters set to missing. Ultimately, 22 samples were successfully fitted to the two-tank model and used for the analysis of all fast- and slow-tank parameters (C0,fast, kfast, C0,slow, kslow). The single-pool sample was retained solely for the analysis of total pool parameters (C0,total and ktotal), as these parameters are available in both the single-pool and dual-pool models. This was not explicitly stated in the original manuscript due to an oversight on our part; we will clarify this in the Methods section of the revised manuscript.
- Regarding the Inclusion of Single-Pool Model Parameters in Analysis of Variance and the Handling of Unbalanced Designs
We appreciate the reviewer’s comments and apologize for any confusion this may have caused. The mixed use of single-pool and dual-pool models may lead to design imbalance. We clarify this as follows: In the main analysis of the fast/slow pool parameters (C0,fast, kfast, C0,slow, kslow) in the manuscript, only 22 dual-pool samples were used. This is because, among the 23 valid samples, one sample (1-HT3-A1) had fewer than four valid observations (fewer than the number of parameters in the dual-pool model) and could only be fitted using a single-pool model; consequently, neither its fast-pool nor slow-pool parameters could be obtained. Therefore, the analysis of the fast- and slow-pool parameters was based on the 22 samples for which the dual-pool model was successfully fitted, and there was no mixing of single-pool and dual-pool parameters in the calculations. Second, the main analysis of the total reservoir parameters (C0,total and ktotal) utilized all 23 samples, as both single-reservoir and dual-reservoir models can output total reservoir parameters. Furthermore, to verify whether the inclusion of the single-pool sample mentioned above would affect the conclusions, we conducted an additional sensitivity analysis: we excluded the data related to that single-pool sample from the total-pool parameter analysis, retained only the 22 dual-pool samples, and re-ran relevant analyses such as ANOVA, LMG hierarchical decomposition, and PLS SEM. The results showed that the inclusion or exclusion of the single-pool sample in the analysis of C0,total and ktotal did not affect any of the core conclusions of this study (SI Table S15).
- Data Distribution and ANOVA Assumptions
Thank you for your comments and suggestions. We fully agree that ANOVA assumptions should be tested when data are non-normal, and we will include the following clarification in the revised manuscript: Among the 22 valid dual-pool samples, kfast values for 3 samples reached the upper limit of the dual-pool model (1.0 d-1), accounting for 13.6% of the total sample size. Based on the reviewer’s suggestion, we performed Shapiro-Wilk normality tests on the kinetic parameters. The results showed that all parameters deviated significantly from a normal distribution (p < 0.05), suggesting that the normality assumption of ANOVA may not be fully satisfied for these parameters. We subsequently conducted a Kruskal-Wallis nonparametric test, the results of which showed a significant cross-over with those of the ANOVA: kfast was significant in the ANOVA (F₃,₁₄ = 5.11, p = 0.013), but did not reach the significance threshold in the Kruskal–Wallis test (H = 6.28, p = 0.099); whereas C0,fast was not significant in the ANOVA (p = 0.056) but was significant in the Kruskal-Wallis test (H = 8.29, p = 0.041); The other parameters (C0,total, C0,slow, kslow) were not significant in either test. Although the p for kfast in the nonparametric test was slightly higher than 0.05, it remained at a relatively low level, and Tukey’s HSD post hoc comparisons under the ANOVA framework confirmed that kfast for the 23-year treatment was significantly higher than that of the other treatments (23 years vs. 1 year: p = 0.006; 23 years vs. 12 years: p = 0.029; 23 years vs. control: p = 0.014). C0,fast was significant in the nonparametric test, further confirming that a treatment effect is indeed present in the data and is not driven by the assumption of a single method.
- Selection of Upper and Lower Limits in Model Fitting and Their Impact on ANOVA
We appreciate the reviewers’ comments and apologize for any confusion caused regarding the basis for our selection of upper and lower limits and the comprehensiveness of our sensitivity analysis. We have provided the following clarification and revisions: The decomposition rate constants for activated carbon pools reported in the literature are generally low. Schädel et al. (2014), based on a three-pool model inversion of circumpolar Arctic permafrost, reported that the turnover time for the fast pool at 5°C was 0.21-0.48 years, corresponding to kfast ≈ 0.005-0.013 d-1; Paul et al. (2006) also indicated, based on measurements of temperate agricultural soils, that the k values for the activated carbon pool range from 0.007 to 0.03 d-1. In this study, the upper limit of kfast was set at 1.0 d-1, which is higher than the values reported in the aforementioned literature, primarily based on the following three considerations: First, differences in sample properties. The low kfast values reported in the literature mostly originate from old permafrost that has undergone long-term decomposition or highly stable agricultural soils, where the easily decomposable components have been largely depleted. In contrast, the samples in this study were collected from fresh permafrost exposed by a thermal karst landslide on the Qinghai-Tibet Plateau; the organic matter had not undergone sufficient decomposition. Chemical analysis revealed high O-alkyl C content and a low C/N ratio (Mu et al., 2024), indicating a high proportion of active carbon and a decomposition rate that may be significantly faster than that of mature soils reported in the literature. Reserving sufficient parameter space for such highly active samples is a necessary prerequisite for avoiding artificial suppression of data variability. Second, avoid statistical bias caused by artificial truncation. In nonlinear fitting, overly tight upper limits can cause parameter estimates to be forcibly truncated, distorting differences between groups. In this study, the measured median kfast was approximately 0.3 d-1, with a maximum of approximately 0.8 d-1; only 13% of the samples reached the upper limit of 1.0 d-1, and even after relaxing the upper limit, the kfast values for these samples did not increase indefinitely, indicating that 1.0 d-1 is sufficient to accommodate the natural variability of the data. If the upper limit were lowered to near the mean value reported in the literature (e.g., 0.1 d-1), the kfast values for most samples would be truncated at the same threshold, and the analysis of variance would reflect model constraints rather than biological differences. Third, sensitivity analysis indicates that the upper limit does not alter the core conclusions. We systematically varied the upper limit of kfast within the range of 0.5-5.0 d-1 and repeated all fittings and analyses of variance. The results show that, under all upper limits, the main effect of recovery duration on kfast remained significant (p = 0.003-0.028). This indicates that the core conclusions are insensitive to the choice of upper limit, and the setting of 1.0 d-1 is robust and verifiable. Finally, we will include in the revised manuscript the information requested by the reviewers regarding the variation in effect sizes and the p-values for the depth effect.
Schädel C, Schuur EAG, Bracho R, Elberling B, Knoblauch C, Lee H, Luo Y, Shaver GR, Turetsky MR (2014) Circumpolar assessment of permafrost C quality and its vulnerability over time using long-term incubation data. Glob Change Biol, 20: 641-652. https://doi.org/10.1111/gcb.12417
Paul EA, Morris SJ, Conant RT, Plante AF (2006) Does the Acid Hydrolysis–Incubation Method Measure Meaningful Soil Organic Carbon Pools?.Proceedings - Soil Science Society of America, 70:1023-1035. https://doi.org/10.2136/sssaj2005.0103.
- Uncertainty in Parameter Estimates
Thank you for your comments and suggestions. We have added a section on parameter uncertainty via Bootstrap analysis to the revised manuscript: The current results of the parameter uncertainty Bootstrap analysis (n = 5000) show that, after accounting for the residual errors from the dual-library model fit for each sample, the median ANOVA p-value for kfast based on recovery years is 0.090, with a 95% confidence interval of [0.001, 0.679]. Across all Bootstrap iterations, approximately 36.6% yielded significant results with p < 0.05. Although this proportion does not exceed half, it is far higher than the expected significance rate of 5% in the absence of an effect, indicating that the year effect maintains a low p-value level under most parameter perturbations. It should be noted that the Bootstrap method only conveys uncertainty in the curve-fitting process, whereas the core conclusions of this study were cross-validated by multiple statistical methods: mixed models (p = 0.003), Tukey HSD post hoc comparisons (p < 0.05), a large effect size (Cohen’s d > 1.4), sensitivity analysis of parameter upper bounds (p < 0.05 for all upper bounds), and the PLS-SEM permutation test (p = 0.005) all consistently support the significant effect of recovery duration on kfast. Therefore, the Bootstrap results can be interpreted, to some extent, as reinforcing the robustness of our conclusions.
- Nested Design
We thank the reviewers for pointing out that the original one-way ANOVA did not account for the nested structure of topsoil and subsoil samples from the same sampling point. This was indeed an oversight in our initial manuscript. In subsequent revisions, we plan to address this issue by following these steps: We identified 12 independent plots based on sample numbers as random-effects factors and reanalyzed all kinetic parameters using a linear mixed model (LMM). The fixed effects were years of recovery, soil layers, and their interaction, while the random intercept was the plot. For the key parameter kfast, the mixed-model results showed that the main effect of years of recovery remained significant (C(Years) [T.23 years]: coefficient = 0.589, z = 2.95, p = 0.003), consistent with the conclusions of the ordinary ANOVA. The depth effect for other parameters (e.g., C0,fast) was also more significant in the mixed model (p = 0.001). The results of the adjusted mixed model will be presented in the Results section of the main text and in the appendix tables; the identification of the nested structure and the process of constructing the mixed model will also be described in detail in the Methods section.
Respond to data analysis and interpretation:
Thank you for your detailed comments on the analysis and interpretation of the data; they are crucial in helping us further clarify the logic behind our data analysis, refine the boundaries of our inferences, and improve the interpretation of our model. Regarding the issues you raised, our responses are as follows:
- Regarding the inference linking minor δ¹³C differences to plant carbon input
Thank you for your comments and suggestions. We acknowledge that we did not directly measure the δ13C of plant tissues, nor did we perform source apportionment. However, the observed trend of δ13C shifting toward negative values along the recovery sequence is consistent with theoretical expectations for rhizosphere sediments. Werth & Kuzyakov (2010) provided a comprehensive review of ¹³C fractionation at the root–microbial–soil interface, noting that CO2 from C3 plant root respiration is approximately -2.1 ± 2.2‰ poorer in ¹³C than root tissue, and that microbial processing further alters the 13C signature of SOM (microbial biomass in C3 soils is enriched by +1.2 ± 2.6‰ relative to SOM). These fractionation patterns indicate that even minor inputs of fresh plant-derived carbon can leave a detectable isotopic fingerprint in SOM. Although the magnitude of the changes in this study is small and does not constitute conclusive proof, it provides supportive, indirect evidence for vegetation restoration replenishing the fast pool. We will adjust the wording in the revised manuscript to explicitly frame the δ13C trend as an “observation that raises a hypothesis” rather than conclusive evidence.
Werth M, Kuzyakov Y (2010) 13C fractionation at the root–microorganisms–soil interface: A review and outlook for partitioning studies. Soil Biology & Biochemistry, 42:1372-1384. https://doi.org/10.1016/j.soilbio.2010.04.009.
- Interpretation of the pH increase from 7 to 8 and improvements in the abiotic environment
Thank you for your comments and suggestions. We agree that microbial decomposition typically reaches its optimum under near-neutral pH conditions. However, in the specific context of this thermokarst sequence, the key lies in the trajectory of pH recovery rather than its absolute value. In the early stages of the landslide, the exposure of deep, acidic parent material caused a sharp drop in surface soil pH (from ~7.5 to ~6.8), which constituted significant acid stress for microorganisms. The subsequent rebound to ~7.8–8.0 essentially represented relief and recovery from this state of stress. Although a pH of 8.0 may be slightly above the optimal range for certain acidophilic microbial groups, compared to the acidified state following the landslide, this alkalization trend likely optimized the overall physicochemical environment by enhancing the availability of nutrients such as phosphorus. Furthermore, our PLS-SEM model also confirmed a significant positive effect of pH on kfast (β = 0.515, p = 0.041). This suggests that, within the specific study interval, the process of recovery toward neutral values was a key factor in promoting the restoration of decomposition function. We will clarify this important distinction between “pH recovery” and “optimal pH” in the revised manuscript.
- Regarding the Specification of the PLS-SEM Model Structure
Thank you for your comments; we apologize for any confusion this may have caused. Regarding the reviewer’s concern about the path between chemical and microbial biomass, we clarify as follows: In the baseline PLS-SEM model, the path between organic carbon chemical quality and microbial biomass was specified from the outset (Fig. 4) and was not omitted. The standardized coefficient for this path is 0.920, but it did not reach statistical significance (Bootstrap p = 0.094), indicating that, under the data conditions of this study, the direct effect of organic carbon chemical quality on total microbial biomass is not yet significant; however, this does not affect the significant direct effect of organic carbon chemical quality on kfast. Regarding concerns about “whether other parameters were included,” we specifically constructed an extended PLS-SEM model in addition to the main model, incorporating soil organic carbon stable isotope (δ13C) and pH as latent variables (SI Table S10). The results showed that both δ13C and pH had significant path effects on kfast (path coefficients of -0.572 and 0.515, respectively; p < 0.05). Environmental variables such as soil moisture and texture were not included in the final structural equation model to avoid overparameterization, as they exhibited weak correlations with kfast (|r| < 0.36) and showed no significant differences across treatments. We will discuss this point more explicitly in the revised manuscript.
Responses to Minor comments:
We thank the reviewers for raising these critical questions, as their comments have helped significantly improve the accuracy and clarity of the manuscript.
- Regarding Terminology, Novelty of Data, and Definition of Basic Concepts
Thank you for your comments and suggestions. Regarding the use of the terms “recovery” and “evolution”: We fully agree with the reviewer’s observation: even 23 years after the landslide, the physical landscape of the landslide area remains in a disturbed state, and true landscape restoration often takes hundreds of years. In the original manuscript, we used the term “recovery” specifically to refer to the process of rebuilding biogeochemical functions, rather than to denote a return to the original geological morphology prior to disturbance. To maintain consistency with the terminology used in Mu et al. (2024) and to avoid implying that the ecosystem has fully recovered, we will uniformly replace the term “recovery” used to describe the time series with “post-landslide succession” or adopt the term “evolution” from the original literature in the revised manuscript. Regarding vegetation conditions, we acknowledge that the initial manuscript lacked detailed descriptions of plant height or composition; in subsequent revisions, we will draw upon the vegetation (e.g., Carex, Kobresia) reported by Mu et al. (2024) to support the context of biological colonization.
- Regarding the reproducibility of Pearson correlation coefficients:
We fully agree with the reviewer’s comments. The Pearson correlation analysis in the original manuscript did, to a large extent, replicate the correlations between soil properties reported by Mu et al. (2024), and these findings did not provide additional incremental information for this study. In subsequent revisions, we will retain only the Pearson correlation results related to the newly fitted kinetic parameters (C0,fast, kfast, C0,slow, kslow), focusing on presenting the relationships between these parameters and the chemical quality of organic carbon and microbial indicators.
- Regarding the definition of “activated carbon”
We thank the reviewer for pointing out the importance of defining this concept. In the context of this paper, “active carbon” does not refer to a specific physical fraction (such as particulate organic carbon, POC) or a particular single chemical substance, but rather constitutes an operational definition. We will clearly define it in the Introduction: active carbon refers to that portion of the carbon pool (C0,fast) representing the potential for mineralization, as inferred from a two-compartment kinetic model, whose core characteristic is a high turnover rate (kfast). From a physicochemical perspective, this pool may be heterogeneous, comprising both easily degradable substrates (e.g., O-alkyl carbon) and relatively complex microbial-derived components (e.g., alkyl carbon) that can be enzymatically degraded within specific cultivation time scales. Through this kinetics-based definition of “activity,” this study aims to capture the dynamic characteristics of soil carbon cycling during the post-landslide succession process, which serves as a valuable complement to the static chemical analysis by Mu et al. (2024). We will clearly articulate this distinction in the revised manuscript.
Responses to Line comments:
We thank the reviewer for such a thorough review of our manuscript; these comments are crucial for improving the quality of our work. Below are our responses to each point:
- Regarding the term “thermal melting”:
We agree with the reviewer’s point. “Thermal melting” may be ambiguous. We will uniformly revise it throughout the text to “thaw” or “thawing,” following the terminology used in Mu et al. (2024).
- Regarding the definition of “active carbon pool capacity”:
We thank the reviewer for pointing out that the definition of “the capacity of the active carbon pool” in the Introduction was unclear. In this study, this term refers to the size of the active carbon pool that can be mineralized (i.e., the proportion or total amount of mineralizable carbon), rather than its turnover rate. In the revised manuscript, when “active carbon pool capacity” is first mentioned in the Introduction, we will directly provide its physical meaning and model expression, and clearly distinguish it from the “decomposition rate constant,” so that readers can accurately understand these two core concepts before proceeding to the Methods section.
- Regarding the mixed use of the term “Plot”:
We thank the reviewer for pointing out this ambiguity. In the revised manuscript, we will follow your suggestion and change “plot” to “site” on Line 129 to eliminate the ambiguity.
- Regarding the placement of interpretive content and supporting references:
We agree with the reviewer’s comments. We will move some speculative or highly interpretive statements to the Discussion section and add corresponding references to support our arguments.
- Regarding the term “land-use”:
The term “land-use” on line 290 is indeed inappropriate for describing naturally occurring thaw slump sequences. We plan to revise it to “thaw slump age gradient,” as used by Mu et al. (2024), to accurately reflect the research context.
- Regarding the relevance of references to thaw slump recovery:
We thank the reviewer for pointing out the discrepancies in the references. We will carefully examine the context of the relevant references, adjust the wording, and ensure that previous research is accurately presented.
- Regarding the interpretation of the PLS-SEM model:
We agree with the reviewer’s precise interpretation. The PLS-SEM model does indeed primarily reveal a strong association between the A/O ratio and kfast. We will revise the relevant statements (such as sentences in the “Results” or “Discussion” sections) to avoid implying that the A/O ratio is the sole driving factor. Instead, we will emphasize that “within the framework of this model, the A/O ratio is the most significant variable explaining the variation in kfast” and acknowledge the limitation that the model does not include all potential factors.
- Speculation regarding the enhancement of alkyl carbon signals by microbial residues:
We agree with the reviewer’s comment that the statement on line 405 of the original manuscript was indeed too definitive and went beyond the scope of direct evidence provided by our current data and the cited literature. Based on the theoretical framework of the microbial carbon pump (MCP) proposed by Liang et al. (2017), residues produced by microbial turnover may be selectively preserved through chemical structure or mineral adsorption. In the revised manuscript, we will amend the relevant sentence to a more cautious, speculative statement, such as: “In the late stages of recovery, the continuous production and selective preservation of microbial residues may contribute to the observed enhancement of alkyl carbon signals.”
Liang C, Schimel JP, Jastrow JD (2017) The importance of anabolism in microbial control over soil carbon storage. Nature Microbiology, 2:17105. https://doi.org/10.1038/nmicrobiol.2017.105.
- Literature support regarding pH changes caused by mineral weathering:
We thank the reviewer for pointing out this critical omission in the literature citations. Upon careful verification, we found that of the two references cited in the original manuscript, Olefeldt et al. (2016) primarily mapped the distribution of thermokarst landscapes in the circumpolar Arctic and estimated carbon stocks, while Walker et al. (2010) focused on exploring the application of chronosequences in ecological succession studies; neither study focused on the geochemical mechanisms of pH recovery following thermokarst landslides. To correct this error and strengthen the foundation of our argument, we have cited literature that is directly based on thermokarst landslides on the Tibetan Plateau and has observed the mechanisms of pH changes (Men et al. 2026; Wang et al. 2023) and have rewritten the relevant paragraphs. Specifically: thermal melting and landslides expose deep-seated acidic parent material, causing a post-disturbance decline in surface soil pH; as the recovery process progresses, the pH gradually rebounds, which is most likely attributable to the combined effects of multiple processes, including the exposure of alkaline parent material, the loss of organic matter, and the migration of soluble salts—accumulated in the low-lying terrain of the landslide area—toward the surface soil through evaporation.
Wang L, Liu G, Ma P, Cheng Z, Wang Y, Li Y, Wu X (2023) Effects of thaw slump on soil bacterial communities on the Qinghai-Tibet Plateau. Catena, 232: 107342. https://doi.org/10.1016/j.catena.2023.107342.
Men X, Fu Z, Zeng L, Wang L, Du W, Gao S, Zhang G, Jiang G, Wu Q (2026) Permafrost collapse alters alpine ecosystem development via microtopographic modification on the Qinghai-Tibet Plateau. Catena, 262, 109692. https://doi.org/10.1016/j.catena.2025.109692
- Discussion on microbial activity at pH 7 and pH 8:
We agree with the reviewer’s reference to the core conclusions of Rousk et al. (2010) and Sinsabaugh et al. (2008). We will revise the text to emphasize that the key point of this study is the “recovery process” of pH from the acidified state (~6.8) following the landslide toward alkalinity, rather than asserting that pH = 8 is superior to pH = 7.
- Regarding the ambiguous phrasing “exhibits a positive effect”:
We thank the reviewer for pointing out the ambiguity of this phrasing. The original statement, “organic carbon chemical mass exhibits a positive effect,” did indeed fail to specify the object of that effect. According to the core findings of this study, organic carbon chemical mass—particularly the chemical reorganization process characterized by alkyl carbon accumulation—plays a key driving role in the recovery of kfast. Therefore, during the subsequent revision process, we rewrote the Conclusions section, replacing the ambiguous “positive effect” with a clear causal statement to clarify that organic carbon mass influences the mineralization process by driving the recovery of kfast.
- Regarding the axis labels and color legend in Figure 1:
We appreciate the reviewer’s suggestion. In the revised manuscript, we will increase the font size of the axis labels and scale labels in Figure 1 and ensure that the image resolution meets the requirements. Additionally, we will clearly explain the meaning of the color coding in the figure legend.
- Regarding the error bars and mean points in Figure 2:
We will add an explanatory note to the figure legend of Figure 2 to ensure that readers can distinguish between the error bars (SE) and the mean points.
- Regarding the update to Note S1:
We thank the reviewer for carefully identifying this discrepancy in the data. Upon verification of Table S14, we confirm that soil depth and its interaction do indeed have a significant effect. The original text in Note S1 stated, “whereas the effects of soil depth and their interaction were not significant (Table S14).” We will correct this to read: “whereas the effects of soil depth and their interaction were significant (Table S14).” We thank the reviewer for helping us correct this error.
Citation: https://doi.org/10.5194/egusphere-2026-3656-AC1
-
AC1: 'Reply on RC1', Hanhan Li, 11 Aug 2026
-
RC2: 'Comment on egusphere-2026-3656', Anonymous Referee #2, 24 Sep 2026
Dear authors and editors,
Please find my notes on this manuscript below.
General Comments:
The researchers attempt to tackle an important unknown in modeling Arctic carbon pools: the impact of thermokarst activity on carbon cycling and turnover, and whether the processes following such disturbances are synchronous and/or monotonic. They explored this via dual-pool index models where they decoupled carbon capacity and carbon decomposition/turnover rates. They concluded that while the fast carbon pool (“capacity”) remains relatively stable (no significant changes along the thermokarst recovery sequence), the carbon decomposition rate does not, rebounding significantly across the thermokarst gradient. Further, they found that “fast carbon” pools were mainly impacted by microbial biomass whereas the decomposition rate was mainly influenced by the chemical quality of the organic carbon. They suggest adding a substrate quality adjustment to future modeling, and highlight the important dynamic between initial carbon supply and composition, microbial activity, vegetation carbon inputs in later recovery stages, pH, and the phenology of microbial community composition and abundance throughout the recovery period. This topic is highly relevant to this journal, and highly applicable and important to improving Arctic carbon modeling; I greatly enjoyed reading this manuscript and think these findings can significantly improve Arctic carbon modeling efforts. I thought overall the text was well written and easy to follow, however, there are some key distinctions that need to be made to improve clarity and there are also some areas that need to be strengthened and amended before publication (discussed below).
Specific comments:
-As this manuscript builds off of a prior piece and uses the same data, but found contradictory conclusions (the 2024 paper finds thaw slumps led to a loss of labile C and a build up of recalcitrant C, whereas this study finds the opposite), the authors need to address these opposing findings and explore what may have led to such different conclusions. This is the most important piece missing from this manuscript.
-A brief characterization of “active carbon”, and “fast” vs “slow” pools would improve accessibility for readers not as familiar with these distinctions and the implications for turnover. As these terms are largely undefined in text, the specific meaning for each of these is vague throughout the manuscript.
-Plots for sites/locations along the thermokarst gradient and when referencing soil cores need to be better distinguished.
-It would have been nice to see a brief description/explanation of the microbial variables assessed (PLFA, F/B ratio, etc.) when stated in the methods, and perhaps a brief explanation of why these microbial parameters are important/useful indicators.
-How is C0 (capacity of fast carbon) characterized throughout the incubation period? Are differences in initial soil carbon vs. soil carbon produced by microbial activity further along in the incubation distinguished (i.e., is new carbon from mineralization distinguished from the existing initial bioavailable carbon)?
-How does vegetation come into play? The researchers mention that the control is an undisturbed meadow, but did any of the recovery stages have colonizing vegetation? Was the potential effect of vegetation colonization (or lack thereof) considered when comparing with the control or other recovery stages? Vegetation as a source of new carbon is mentioned in the discussion, but was vegetation colonization observed in some or all recovery stages? What kind of vegetation? Certain vegetation can not only enhance soil carbon but also aid in N-fixation, thereby improving soil quality for microbial turnover. I think this should be more pointedly addressed as being excluded from the proposed model, despite it also being an important mechanism for carbon cycling along thermokarst gradients. Additionally, the different types of vegetation (if present) should be briefly mentioned as the specific types and species may have important implications for carbon cycling.
-Was increased surface water runoff or transport observed at these thaw slump recovery sites? Was potential leaching or lateral transport of carbon considered?
-How do the analyses account for the use of shallow vs deeper soil from the cores? This should be expanded upon a bit more, especially since authors note that were significant differences in carbon pools in topsoil compared to subsoil in lines 268 – 270.
-Authors should provide further insight as to why the Cslow and kslow did not significantly differ along the gradient – what was the soil composition of the Cslow pool? Were there differences in soil moisture, temperature, microbial activity, etc., between the topsoil and subsoil, and the Cfast and Cslow pools? It makes sense that bioavailable carbon be most impacted by thermokarst disturbances, yet it is interesting that in 23 years the more recalcitrant carbon pools and decomposition rates remained stable / largely unchanged, though research indicates that microbial communities can shift to take advantage of recalcitrant carbon.
- Why would kfast significantly increase with recovery duration, while the main mechanism of decomposition here (microbial activity) appears to decline with recovery period and also potentially negatively impacts kfast (Fig 4)? In Line 280, authors note that the positive correlation between kfast and the A/O ratio contradicts previous findings, but perhaps unique microbial communities that are more efficient at decomposing less labile carbon sources can be found along thermokarst gradients. This, and other reasons the decomposition may accelerate despite rising A/O ratios, is discussed in the discussion section, and I realize the microbial community specifics are out of the scope of the manuscript in its current form, but the addition of this – even if only briefly – could greatly strengthen this argument.
-In the statistical testing, were models assessed for meeting assumptions (parametric vs nonparametric)? Not many ecological datasets are cleanly parametric and often require non-parametric or mixed-effects approaches; if this dataset was not normally distributed or didn’t have normality of homogeneity, then it could make these analyses unreliable. Further, did the authors consider testing for autocorrelation (due to repeated measurements of the same cores over a period of time), or conducting timeseries-specific analyses that take the autocorrelation of repeated measurements from the same core into account?
-How did authors reconcile the majority of the samples being analyzed with dual-pool models and the one core with a single-pool model?
-Line 500: Given the authors note the importance of pH, why was pH not included in this substrate quality adjustment term?
Technical comments:
-Fig 1: Include a legend for the datapoints / line colors here.
-Fig 1: “Land-use” isn’t the appropriate word for this gradient – I recommend simply using the timeframe range or associating it with early, mid, and late-stage post-thermokarst.
-Fig 2: There is quite high variability among the 3 data points for C0,fast topsoil during recovery year 12 – what are some potential reasons for this?
-Fig 4: Why only these parameters in the PLS-SEM figure, rather than including all those showing high impacts in Fig 3b?
-Line 338 – 339: Authors note soil moisture decreased with restoration duration – was this topsoil, subsoil, or an average of both? Were potential interaction effects of soil moisture on microbial community at different stages of restoration also assessed?
-Line 405 and thereafter: Much of this section raises good, probable points but seems to be largely speculation as it does not include citations. Citations should be added here where applicable, and if no evidence exists, language should be adjusted to make clear that these statements are speculative.
-Lines 469 – 471: Quantifying these changes in pH within the text here would be helpful. What was the pH directly after the disturbance? What was the rate of change across years? Was the change linear, and synchronous among stages / top or subsoil?
-In lines 498 – 503 authors suggest adding a substrate quality adjustment function for fast pool rate constants based on A/O, but within this research they indicated that the general trend of higher A/O correlating with lower decomposition rates did not hold, potentially influenced by vegetation carbon inputs and microbial residues. How would this adjustment function take this possible deviation into account?
Line comments:
-Line 149: There seems to be an errant period here where there should be a comma.
-Lines 167 – 168: The “fast” in “kfast” is not written as a subscript.
-Lines 182 - 187: Tone shifts to present tense when rest of text is in past tense.
-Line 251: The semicolon should be a period.
-Line 259: The semicolon should be a period.
-Line 306: The semicolon should be a period.
-Line 355: The CO,fast are not written as subscript.
-Line 358: Kfast not written as subscript.
-Line 486-487: Seems like the sentence is incomplete or needs to be re-worded/joined with the next sentence.
-Table S14: Indicates a significant (or borderline) effect of soil depth (k upper bound = 0.5) though is said not to be significant in text.
Citation: https://doi.org/10.5194/egusphere-2026-3656-RC2 -
AC2: 'Reply on RC2', Hanhan Li, 29 Sep 2026
General Comments:
The researchers attempt to tackle an important unknown in modeling Arctic carbon pools: the impact of thermokarst activity on carbon cycling and turnover, and whether the processes following such disturbances are synchronous and/or monotonic. They explored this via dual-pool index models where they decoupled carbon capacity and carbon decomposition/turnover rates. They concluded that while the fast carbon pool (“capacity”) remains relatively stable (no significant changes along the thermokarst recovery sequence), the carbon decomposition rate does not, rebounding significantly across the thermokarst gradient. Further, they found that “fast carbon” pools were mainly impacted by microbial biomass whereas the decomposition rate was mainly influenced by the chemical quality of the organic carbon. They suggest adding a substrate quality adjustment to future modeling, and highlight the important dynamic between initial carbon supply and composition, microbial activity, vegetation carbon inputs in later recovery stages, pH, and the phenology of microbial community composition and abundance throughout the recovery period. This topic is highly relevant to this journal, and highly applicable and important to improving Arctic carbon modeling; I greatly enjoyed reading this manuscript and think these findings can significantly improve Arctic carbon modeling efforts. I thought overall the text was well written and easy to follow, however, there are some key distinctions that need to be made to improve clarity and there are also some areas that need to be strengthened and amended before publication (discussed below).
Response: Thank you very much for your positive feedback and constructive suggestions regarding this study, as well as for your recognition of the manuscript’s quality and value. Based on your overall comments, we will revise the manuscript and further clarify the scope of application for relevant terminology, mechanisms, and conclusions throughout the text. Please see our subsequent reply for specific details.
Specific comments:
-As this manuscript builds off of a prior piece and uses the same data, but found contradictory conclusions (the 2024 paper finds thaw slumps led to a loss of labile C and a build up of recalcitrant C, whereas this study finds the opposite), the authors need to address these opposing findings and explore what may have led to such different conclusions. This is the most important piece missing from this manuscript.
Response: We thank the reviewer for pointing out a potential contradiction between our findings and those of Mu et al. (2024), and we apologize for any confusion this may have caused. First, we would like to clarify that the two studies are not contradictory but rather offer complementary interpretations of the same recovery sequence from different perspectives. Mu et al. (2024) used physical fractionation and static chemical analysis methods, such as ¹³C NMR, to reveal the compositional evolution of the soil carbon pool during the recovery of a thermal karst landslide: the content of particulate organic carbon (POC) decreased with the duration of recovery, while inert carbon fractions (such as mineral-associated organic carbon, MAOC, and alkyl carbon) gradually accumulated. Our study does not dispute this evolution of carbon pool composition. Instead, by performing kinetic decoupling of the cumulative CO₂ release curve using a two-compartment exponential model, we provide two key insights from the perspective of carbon turnover dynamics. (i) Changes in C0,fast indicate that, although the original readily degradable carbon was lost in large quantities after the landslide as revealed by kinetic fitting, the fast pool did not disappear as a result. Instead, the composition of the fast pool underwent a qualitative shift during the recovery process, with the lost readily degradable carbon being replaced by new components. While these new components are indeed more “stubborn” than the original active carbon, they still exhibit a significantly faster turnover under cultivation conditions compared to the inert carbon in the slow pool; therefore, the model classifies them as the “fast pool.” (ii) The significant rebound in kfastreveals the recovery dynamics of the fast pool’s turnover rate. In the early post-landslide period (years 1–12), kfast remained suppressed; by year 23, kfast had recovered significantly and reached or even slightly exceeded the control level. This does not imply that the decomposition rate has returned to its pre-landslide state, but rather indicates that microbial mineralization is gradually recovering, even though the substrate is no longer the original high-quality organic carbon but rather a substitute active carbon with more persistent chemical properties. Mu et al. (2024) addressed the questions of “what the carbon pool is and how it changes,” while we used C0,fast and kfast to further address the questions of “whether the carbon pool’s turnover rate has recovered and whether its capacity has changed accordingly.” We will integrate these perspectives in the discussion section of the revised manuscript and cite Mu et al. (2024) to demonstrate the further findings of this study building upon the work of Mu et al. (2024).
-A brief characterization of “active carbon”, and “fast” vs “slow” pools would improve accessibility for readers not as familiar with these distinctions and the implications for turnover. As these terms are largely undefined in text, the specific meaning for each of these is vague throughout the manuscript.
Response: Thank you for bringing this issue to our attention. In the current draft, the terms “active carbon,” “fast pool,” and “slow pool” are primarily used as working terms without clear definitions upon their first appearance. This can indeed cause confusion for readers unfamiliar with the two-pool dynamics model and may also hinder their understanding of the concept of carbon pool turnover. We will provide clear definitions for these terms in the revised version. Specifically, “active carbon” refers broadly to the portion of soil organic carbon that turns over relatively quickly, is readily utilized by microorganisms, and is sensitive to disturbances and environmental changes. Meanwhile, “fast pool” and “slow pool” are kinetic carbon pools classified based on the two-pool index model: the fast pool has a large decomposition rate constant (kfast), mineralizes rapidly under cultivation conditions, and reaches its potential mineralizable capacity (C0,fast) relatively early; the slow pool, on the other hand, has a smaller decomposition rate constant (kslow) and undergoes continuous mineralization at a lower rate throughout the entire incubation period. It is important to emphasize that this distinction between fast and slow pools is an apparent kinetic classification, reflecting the model’s fit to the mineralization curve during the incubation period; it does not directly correspond to carbon fractions that can be physically separated by chemical methods.
-Plots for sites/locations along the thermokarst gradient and when referencing soil cores need to be better distinguished.
Response: Thank you for pointing out this issue. In the current manuscript, the distinctions between study sites, plots, soil cores, and soil horizons along the thermoklinic landslide gradient are not sufficiently clear. According to the source literature, Mu et al. (2024), in subsequent revisions, we plan to clearly distinguish between “site” (study area/geographic location), “plot” (a study site at a specific recovery stage), “quadrat” (repeated quadrats within a study site), “core/profile” (soil core or profile), and “composite sample” (analytical sample mixed by depth). This will ensure a clear distinction between different locations, recovery stages, replicate quadrats, and soil layers along the thermokarst gradient.
-It would have been nice to see a brief description/explanation of the microbial variables assessed (PLFA, F/B ratio, etc.) when stated in the methods, and perhaps a brief explanation of why these microbial parameters are important/useful indicators.
Response: Thank you for your suggestion. We agree that the first mention of PLFA, the F/B ratio, and other terms in the Methods section of the original manuscript did not include brief definitions or their ecological significance, which hinders reader understanding. In the revised manuscript, we will add the following at the point of first mention: Microbial activity is a key driver of soil carbon mineralization. Specifically, Total PLFA reflects living microbial biomass; the F/B ratio can be used to characterize changes in the relative abundance and community structure of fungi and bacteria, suggesting adjustments in microbial resource utilization strategies; and BG activity reflects the enzymatic degradation potential of microorganisms toward active carbon components. Therefore, we will include Total PLFA, the F/B ratio, and BG as microbial variables in the PLS-SEM model to characterize the potential impact of microorganisms on soil carbon mineralization at three levels: biomass, community structure, and enzymatic function.
-How is C0 (capacity of fast carbon) characterized throughout the incubation period? Are differences in initial soil carbon vs. soil carbon produced by microbial activity further along in the incubation distinguished (i.e., is new carbon from mineralization distinguished from the existing initial bioavailable carbon)?
Response: Thank you for raising this important question. C0,fas is a model parameter obtained by fitting the 248-day cumulative CO2-C release curve using a two-pool index model. The specific fitting equation is: C(t)= C0,fast (1-e−kfast(t))+C0,slow(1-e-kslow(t)). C0,fas represents the total potential amount of carbon that the fast pool can contribute to mineralization throughout the entire incubation period, while C0,slow represents the total potential amount of carbon that the slow pool can contribute to mineralization; C0,total=C0,fast+C0,slow. Therefore, the unit of C0,fast is mg CO₂-C g⁻¹ SOC, and it describes the apparent capacity of the fast mineralization fraction under given incubation conditions. This parameter is obtained through a nonlinear least-squares fit; it is not independently measured, nor does it represent the actual, directly distinguishable initial amount of bioavailable carbon present in the soil at the start of the incubation.
It should be noted that the current model and code cannot distinguish between the bioavailable carbon already present in the initial soil and the carbon generated, recycled, or re-entered into the fast pool by microbial activity during the incubation. The model assumes that the sizes of the two carbon pools remain constant during the incubation period, with CO₂ release controlled solely by a first-order mineralization rate constant; it does not explicitly simulate microbial growth, death, detritus formation, extracellular enzyme dynamics, stimulation effects, or the process of microbial synthesis of new carbon. Therefore, the cumulative CO₂-C release measured during the incubation represents a net mineralization signal, reflecting the combined contribution of the initial organic carbon and the carbon that may have been turned over and utilized by microorganisms during the incubation period; however, it is not possible to quantitatively distinguish which portion originates from the initial active carbon and which from new carbon produced by microbial activity.
In the revised manuscript, we will add a paragraph in the Methods section providing a clear definition of C0,fast and explaining the model assumptions. We will also clarify in the Limitations section that future studies may need to incorporate isotope tracing or microbial residue markers to further elucidate the sources and turnover processes of fast pool carbon.
-How does vegetation come into play? The researchers mention that the control is an undisturbed meadow, but did any of the recovery stages have colonizing vegetation? Was the potential effect of vegetation colonization (or lack thereof) considered when comparing with the control or other recovery stages? Vegetation as a source of new carbon is mentioned in the discussion, but was vegetation colonization observed in some or all recovery stages? What kind of vegetation? Certain vegetation can not only enhance soil carbon but also aid in N-fixation, thereby improving soil quality for microbial turnover. I think this should be more pointedly addressed as being excluded from the proposed model, despite it also being an important mechanism for carbon cycling along thermokarst gradients. Additionally, the different types of vegetation (if present) should be briefly mentioned as the specific types and species may have important implications for carbon cycling.
Response: Thank you for raising this very important issue. Vegetation colonization and its species composition are among the key mechanisms of carbon cycling in the recovery sequence of thermokarst landslides. The current draft does not sufficiently address this topic; it is indeed necessary to more clearly explain the role of vegetation in the recovery process and to highlight the limitation that the current model does not explicitly account for vegetation colonization. According to the source literature, Mu et al. (2024), the vegetation type in the study area is alpine meadow, dominated primarily by Carex and Kobresia.
In the current analysis, the two-pool index model decouples the 248-day cumulative CO₂-C release curve into a fast pool and a slow pool; the resulting C0,fast and kfast are apparent parameters of net mineralization kinetics. These parameters may implicitly reflect the combined supply of microbial-available carbon during the cultivation period, including bioavailable carbon already present in the initial soil, residual organic matter from decomposition, and carbon potentially generated by microbial turnover and detritus recycling during cultivation; however, the model does not explicitly simulate vegetation inputs, root exudates, stimulation effects, or microbial carbon formation. Therefore, vegetation colonization—as a source of new carbon and a mechanism for improving the quality of substrates for microbial turnover—is not explicitly represented in the current model. The vegetation carbon input mentioned in our discussion is primarily based on δ¹³C results, which indirectly suggest that fresh plant carbon, particularly root exudates, may enter the fast pool during the late restoration phase; however, this evidence cannot distinguish specific vegetation types nor assess the potential promoting effect of nitrogen-fixing plants on nitrogen supply and microbial turnover.
We will address the issue of vegetation colonization more explicitly in the revised manuscript. Given that the natural meadows in the study area are dominated by two non-nitrogen-fixing vegetation types—Carex and Kobresia— we will briefly mention these two vegetation types in the discussion, noting that they may influence soil carbon turnover through root exudates, litter input, and nutrient cycling. However, since the current two-tank model does not explicitly simulate vegetation colonization or vegetation carbon input, the model parameters reflect the net mineralization resulting from the combined effects of vegetation input, microbial residue turnover, and environmental factors, rather than the independent effect of vegetation colonization. Second, the “Limitations” section explicitly acknowledges that the current study did not obtain specific species composition data for vegetation colonization during the recovery phase; consequently, it is not possible to assess the potential impacts of different vegetation types on rapid carbon pool capacity and decomposition rates. Future research could further integrate vegetation community surveys and plant functional traits to more directly elucidate the role of vegetation colonization in carbon cycling during the recovery of thermokarst landslides.
-Was increased surface water runoff or transport observed at these thaw slump recovery sites? Was potential leaching or lateral transport of carbon considered?
Response: Thank you for raising this question. We apologize, but it is important to clarify that neither the dataset used in this study nor the current analysis includes observational data on surface runoff, lateral carbon migration, or leaching fluxes from the thermokarst collapse recovery sites. Therefore, we are unable to quantitatively assess the contributions of these processes to rapid carbon pool capacity and decomposition rates, and the current two-tank model does not incorporate hydrological migration or leaching processes. However, in this environmental context, potential carbon leaching or lateral migration may be an important factor influencing carbon mineralization. We plan to briefly address potential carbon leaching or lateral migration when discussing the carbon balance of the fast pool in the revised manuscript; for example, in Section 4.1, we will note that following a thermal collapse, some readily degradable carbon may be lost through surface runoff and lateral migration, resulting in an initial loss of active carbon; however, during the recovery process, new inputs of alternative active carbon—such as root exudates and microbial debris—continually replenish the fast pool, preventing a significant decline in C0,fast. At the same time, leaching may selectively remove some readily available substrates, leaving the remaining fast-pool carbon with relatively more persistent chemical properties. This may further promote functional adaptation of the microbial community and increase mineralization activity per unit biomass, thereby supporting the significant rebound in kfast in 2023. We explicitly state that this is an inference.
-How do the analyses account for the use of shallow vs deeper soil from the cores? This should be expanded upon a bit more, especially since authors note that were significant differences in carbon pools in topsoil compared to subsoil in lines 268 – 270.
Response: We agree that the treatment of topsoil and subsoil needs to be explained more clearly in the manuscript. In this study, topsoil (0–15 cm) and subsoil (15–30 cm) were collected simultaneously during each restoration phase; both were taken from different depths of the same soil core, constituting stratified sampling. Since the topsoil and subsoil inherently exhibit significant differences in organic carbon content, microbial biomass, moisture, and oxygen conditions, we treated “Depth” as a fixed factor in the analysis of variance (ANOVA) and included the “Years × Depth” interaction term to examine the effects of restoration duration, soil depth, and their interaction on kinetic parameters. When the interaction was significant, it indicated that the effect of restoration duration varied with depth; therefore, we conducted separate Tukey HSD post hoc tests at each depth level; When the interaction was not significant, it indicated that the effect of restoration duration was consistent across different depths; therefore, we conducted post hoc tests uniformly for significant main effects. This corresponds to the significant differences in carbon pools between topsoil and subsoil that you mentioned; for example, C0,fast was significantly higher in topsoil than in subsoil. In the graphical presentation, Figure 2 displays cumulative mineralization curves in panels organized by depth, while Figure 3 distinguishes topsoil and subsoil using different colors and shapes.
However, in the subsequent driver analysis (LMG and PLS-SEM), due to the limited total sample size (n = 23), we followed the SEM modeling approach described in Mu et al. (2024) and pooled the topsoil and subsoil samples for analysis to maintain statistical power. This means that depth differences were explicitly tested and reported in the ANOVA but were not treated as independent stratification variables in the driver analysis pathways. We acknowledge this as a limitation and will clarify this point more explicitly in the revised manuscript. We will subsequently detail how “Depth” was handled in model fitting in the Methods section and explain the impact of the topsoil–subsoil distinction on our conclusions in the Discussion. Thank you again for raising this important point.
-Authors should provide further insight as to why the Cslow and kslow did not significantly differ along the gradient – what was the soil composition of the Cslow pool? Were there differences in soil moisture, temperature, microbial activity, etc., between the topsoil and subsoil, and the Cfast and Cslow pools? It makes sense that bioavailable carbon be most impacted by thermokarst disturbances, yet it is interesting that in 23 years the more recalcitrant carbon pools and decomposition rates remained stable / largely unchanged, though research indicates that microbial communities can shift to take advantage of recalcitrant carbon.
Response: Thank you for raising this issue. We agree that the revised manuscript needs to further explain why C0,slow and kslow do not show significant changes along the recovery gradient. First, it should be noted that since Mu et al. (2024) have already reported differences between the topsoil and subsoil, the revised manuscript will include additional discussion on the potential differences between C0,slow and C0,fast based on the findings of this study. Second, we need to further elucidate the relationships among C0,slow, C0,fast, kslow, and kfast. Currently, based on the results from the dual-pool dynamic model, we believe that the turnover mechanisms of soil organic carbon after a landslide exhibit distinct pool-specific differentiation. For the fast pool, although a large amount of readily degradable carbon was lost in the early stages of the landslide, C0,fast exhibited significant changes with both the duration of recovery and soil depth, indicating that the fast pool did not disappear but rather underwent a qualitative replacement: newly introduced substitute active carbon, such as root exudates and microbial residues, replenished the lost high-quality substrate; Although these new components are chemically more refractory than the original active carbon, they still exhibit relatively rapid turnover under cultured conditions and are therefore classified as a fast pool by the model. Meanwhile, kfast was suppressed in the early post-landslide period but had significantly rebounded to control levels—and even slightly exceeded them—by year 23. This is not a simple restoration of the original decomposition function but is more likely the result of microbial mineralization capacity readapting to the new substrates.
In contrast, the dynamics of the slow pool were markedly different. C0,slow showed no significant changes with either recovery duration or depth, suggesting that the slow pool, as a highly inert background pool, maintained a relatively stable total mass over the decades-long recovery period. This may stem from a dynamic equilibrium between the transfer of relatively more active components from the original slow pool to the fast pool and the simultaneous replenishment of new, more refractory carbon. More importantly, kslowvaries significantly only with soil depth, not with restoration duration. This suggests that the decomposition of inert carbon in the slow pool is controlled to a greater extent by abiotic factors—such as temperature, moisture, oxygen diffusion, and the degree of mineral protection in deep soil layers—while the recovery of microbial mineralization capacity has a relatively limited impact on it. Of course, this does not mean that microorganisms are completely unable to utilize slow-release carbon; rather, it likely reflects the weak mineralization signals from the slow-release pool during the 248-day incubation, high estimation uncertainty, and the current indicators’ insufficient ability to capture the degradation function of the slow-release pool. Therefore, the stability of C0,slow and kslow along the recovery gradient should not be interpreted as the slow-release pool being completely unaffected by disturbances, but rather as the slow-release pool having a slow turnover, a lagged response, and being dominated by abiotic protection. We will include the above explanation in the revised manuscript.
- Why would kfast significantly increase with recovery duration, while the main mechanism of decomposition here (microbial activity) appears to decline with recovery period and also potentially negatively impacts kfast (Fig 4)? In Line 280, authors note that the positive correlation between kfast and the A/O ratio contradicts previous findings, but perhaps unique microbial communities that are more efficient at decomposing less labile carbon sources can be found along thermokarst gradients. This, and other reasons the decomposition may accelerate despite rising A/O ratios, is discussed in the discussion section, and I realize the microbial community specifics are out of the scope of the manuscript in its current form, but the addition of this – even if only briefly – could greatly strengthen this argument.
Response: Thank you for raising this very important question. The apparent contradiction you pointed out does indeed require a clearer explanation: kfast increases significantly with restoration duration, yet total microbial biomass decreases over the same period; furthermore, in the PLS-SEM analysis, the path from microbial biomass to kfast is not significant and is even negative in direction. This may suggest that total microbial biomass is not the primary factor controlling changes in kfast. PLFA reflects the size and coarse composition of the microbial community; it is not equivalent to metabolic activity, enzyme synthesis capacity, carbon utilization efficiency, or the ability to degrade specific substrates per unit of biomass. In the later stages of recovery, microbial biomass may be lower, but the mineralization activity per unit of biomass may be higher, or the community structure may have undergone functional adaptation, thereby maintaining or even increasing the turnover rate of the fast pool despite the decline in biomass. Our results also indirectly support this: total PLFA showed no significant correlation with kfast; the combined independent contribution of microbial indicators in LMG was less than 15%, whereas the chemical quality of organic carbon—particularly PC1 and the A/O ratio—was the primary driver of kfast; in PLS-SEM analysis, the direct effect of recovery duration on kfast was significant, whereas the microbial biomass pathway was not significant.
We will further explain the seemingly counterintuitive result that the A/O ratio is positively correlated with kfast in the discussion section. Alkyl carbon measured by NMR is not a chemically homogeneous substance; it includes both hard-to-degrade components of plant origin, such as long-chain lipids and cutin, as well as relatively easily utilizable components of microbial origin, such as branched-chain fatty acids, peptidoglycans, and some protein side chains. The continued accumulation of microbial residues during the late stages of restoration can simultaneously increase both the alkyl carbon signal and the A/O ratio; however, a significant portion of these residues still consists of active components that can be mineralized relatively quickly. At the same time, the influx of fresh plant carbon resulting from vegetation restoration—particularly root exudates—may enhance microbial metabolic activity through a stimulatory effect, thereby accelerating the decomposition of pre-existing organic carbon. Therefore, an increase in kfast does not necessarily require an increase in microbial biomass; rather, it more likely reflects the combined effects of changes in substrate chemical quality, functional adaptation of the microbial community, and increased activity per unit of biomass.
Thank you for proposing the unique mechanism whereby microbial communities may more efficiently degrade less active carbon sources; we fully agree that it helps explain the observations in this study. Although the current PLFA data provide only information on microbial biomass and rough taxonomic groups—and cannot directly identify changes at the species level or in functional genes—we will briefly address this possibility in the Discussion section of the revised manuscript, noting that during the later stages of recovery, the microbial community may undergo structural adjustments, shifting toward groups better suited to utilizing more recalcitrant substrates, such as certain fungi, actinomycetes, or specific degradative bacteria, or by adjusting carbon utilization efficiency to allocate more carbon toward respiration rather than biomass synthesis. We will also explicitly state that this explanation constitutes a mechanistic inference based on the existing literature, which cannot be directly verified by the current data and requires future validation using techniques such as metagenomics, metatranscriptomics, or stable isotope probing. We believe that even a brief addition of this paragraph will significantly strengthen the completeness of the argument regarding the kfast rebound mechanism and address the issue of microbial functional adaptation that you raised. Thank you again for raising this important point.
-In the statistical testing, were models assessed for meeting assumptions (parametric vs nonparametric)? Not many ecological datasets are cleanly parametric and often require non-parametric or mixed-effects approaches; if this dataset was not normally distributed or didn’t have normality of homogeneity, then it could make these analyses unreliable. Further, did the authors consider testing for autocorrelation (due to repeated measurements of the same cores over a period of time), or conducting timeseries-specific analyses that take the autocorrelation of repeated measurements from the same core into account?
Response: Thank you for raising this important statistical issue. We fully agree that ecological datasets often do not meet the assumptions of parametric tests; therefore, it is essential to test for normality and homogeneity of variance beforehand and to use nonparametric or mixed-effects methods when necessary. At the same time, since the same soil core or culture sample is measured repeatedly at multiple time points, autocorrelation must also be taken into account. We have conducted analytical validation for both of these issues and will report on them in detail in a subsequent revised manuscript.
First, for each dynamic parameter, we conducted hypothesis tests, including the Shapiro-Wilk normality test for residuals and the Levene test for homogeneity of variances, and then selected parametric or nonparametric methods based on the results. For C0,fast, the residuals satisfied the assumptions of normality (Shapiro-Wilk, p = 0.791) and homogeneity of variances (Levene, p = 0.622); therefore, a parametric two-way ANOVA with Tukey’s HSD post hoc test was used. The results showed that the main effect of recovery duration was marginally non-significant (p = 0.056), while the main effect of depth was significant (p = 0.001), with topsoil values significantly higher than those of subsoil. For kfast, the residuals also satisfied the assumptions of normality (p = 0.418) and homogeneity of variances (p = 0.269); therefore, parametric ANOVA and Tukey’s HSD were applied. The results showed that the main effect of restoration duration was significant (p = 0.013), with 23 years significantly higher than 1 year and 12 years, while the effect of depth was not significant. For C0,slow, the residuals satisfied the assumptions of normality (p = 0.110) and homogeneity of variance (p = 0.387); parametric ANOVA was used, and the results showed that neither the recovery period nor the depth had a significant effect. For C0,totaland kslow, the residuals did not follow a normal distribution (Shapiro-Wilk p=0.047 and p=0.021); therefore, we switched to the nonparametric SRH test and Mann-Whitney-Holm post hoc comparisons. The results showed that neither the recovery time nor the depth effect for C0,total was significant, while for kslow, only the depth effect was marginally significant (p=0.065).
Second, to address the autocorrelation that may arise from repeated measurements of the same soil core or the same cultivation sample at multiple time points, we further re-fitted the cumulative mineralization curves for each sample using an AR(1) iterative GLS mixed model, explicitly incorporating first-order autocorrelation into the error structure, and then conducted hypothesis tests and ANOVA again based on the re-estimated kinetic parameters. The results showed that the core conclusions remained robust after accounting for autocorrelation and were even more clearly established. For kfast, the recovery duration effect was significant under the AR(1) model (p = 0.002), which was stronger than the OLS result (p = 0.013); the value at 23 years remained significantly higher than those at 1 and 12 years. For C0,fast, the recovery duration effect became significant under the AR(1) model (p = 0.015), while the depth effect remained significant (p = 0.001). For C0,total and C0,slow, the recovery period and depth effects remain insignificant under the AR(1) model (p > 0.05). For kslow, the depth effect is significant under the nonparametric test for AR(1) (p = 0.017), while the recovery period effect remains insignificant. We further conducted a whitening residual diagnosis; the proportion of significant Ljung-Box tests for the original residuals was 47.8%, which decreased to 13.0% after whitening, indicating that the AR(1) model effectively removed most of the autocorrelation. In addition, we performed a sensitivity analysis on parameter upper bounds and found that the recovery-period effect for kfast remained significant under all upper bound settings, suggesting that the core conclusions do not depend on arbitrarily set parameter boundaries.
The current core conclusions remain consistent across parametric and nonparametric methods as well as after AR(1) correction. We will include the aforementioned validation results in a subsequent revised draft. Thank you once again for raising this important point, which has helped us significantly improve the rigor and reliability of our statistical analysis.
-How did authors reconcile the majority of the samples being analyzed with dual-pool models and the one core with a single-pool model?
Response: Thank you for your comments, and we apologize for any confusion this may have caused. Regarding the possibility that mixing single-pool and dual-pool model parameters could lead to design imbalance, we would like to clarify as follows: In the main analysis of the fast/slow pool parameters (C0,fast, kfast, C0,slow, kslow) presented in the manuscript, only 22 dual-pool samples were used. This is because, among the 23 valid samples, one sample (1 HT3 A1) had fewer than four valid observations (fewer than the number of parameters in the dual-pool model) and could only be fitted using a single-pool model; consequently, neither the fast-pool nor the slow-pool parameters could be obtained for this sample. Therefore, the analysis of the fast/slow-pool parameters was based on the 22 samples for which the dual-pool model was successfully fitted, and there was no mixing of single-pool and dual-pool parameters in the calculations. Second, the main analysis of the total reservoir parameters (C0,total and ktotal) utilized all 23 samples, as both single-reservoir and dual-reservoir models can output total reservoir parameters. Furthermore, to verify whether the aforementioned single-pool sample would affect the conclusions, we conducted an additional sensitivity analysis: we excluded the data related to that single-pool sample from the total-pool parameter analysis, retained only the 22 dual-pool samples, and re-ran relevant analyses such as ANOVA, LMG hierarchical decomposition, and PLS SEM. The results showed that the inclusion or exclusion of the single-library sample in the analysis of C0,total and ktotal did not affect any of the core conclusions of this study (SI Table S15).
-Line 500: Given the authors note the importance of pH, why was pH not included in this substrate quality adjustment term?
Response: Thank you for pointing this out. We agree that pH is an important factor in this study, and the extended PLS-SEM analysis also shows that pH has a significant positive effect on kfast (β = 0.515, p = 0.041). pH is not included in the currently proposed substrate mass adjustment term kfast=kbase×(A/O)γ, primarily because we wish to distinguish between the “chemical mass effect of the substrate” and “non-biological environmental regulatory effects.” The A/O ratio represents the chemical quality of organic carbon itself, whereas pH, temperature, moisture, and other factors are environmental conditions that modulate degradation rates by influencing microbial activity, enzyme function, and substrate accessibility—rather than being components of substrate quality. Therefore, we treat pH as an independent environmental control variable rather than part of the substrate quality adjustment term, which is also consistent with the independent direct path of pH to kfast in the extended PLS-SEM analysis. In other words, kbase already implicitly accounts for the fundamental control effects of physical environmental factors such as temperature and moisture, while pH can act as another independent environmental moderator, operating in parallel with the substrate quality effect.
However, we fully agree that pH may indirectly modulate the effect of substrate quality on kfast by altering substrate availability or microbial community composition. Therefore, in the revised manuscript, we will add the following clarification in the Discussion section: the current parameterization scheme is intended to avoid confounding environmental effects with substrate quality effects; in the future, pH could be further incorporated into the model as an independent modifying factor—for example, by writing kfast=kbase×(A/O)γ×f(pH)—or by including pH in the environmental control term of kbase. However, we will explicitly state that this extension requires data from more sites and longer recovery sequences to validate.
Technical comments:
-Fig 1: Include a legend for the datapoints / line colors here.
Response: Thank you for pointing this out. We agree that Figure 1 currently lacks a clear explanation of the colors used for the data points and the fitted lines, which may leave readers uncertain about what the different colors represent. In the revised manuscript, we will add a legend to Figure 1 to clearly explain that the different colors represent the three replicate samples under each treatment, the scatter points indicate the observed cumulative CO₂-C release values for each sample, and the solid lines represent the optimal model-fitted curves for that sample—that is, the two-compartment model or the one-compartment model. Additionally, we will include a note in the figure caption explaining the correspondence between the colors, scatter plots, and fitted lines to ensure readers can accurately understand the meaning of each element in the figure. Thank you again for your suggestion.
-Fig 1: “Land-use” isn’t the appropriate word for this gradient – I recommend simply using the timeframe range or associating it with early, mid, and late-stage post-thermokarst.
Response: Thank you for pointing this out. We agree that “land-use” is not an appropriate term to describe the recovery gradient of thaw slumps, as this study does not address changes in land-use types but rather the recovery time series following the occurrence of thaw slumps. In the revised manuscript, we will uniformly replace “land-use durations” with “time since thaw slump initiation” or “recovery stage,” and clearly specify their corresponding categories as Control, 1-year (active/early stage), 12-year (middle stage), and 23-year (stable/late stage). Relevant wording in the figure captions and main text will also be updated accordingly to prevent readers from misinterpreting the recovery stages as changes in land use. Thank you again for your suggestion.
-Fig 2: There is quite high variability among the 3 data points for C0,fast topsoil during recovery year 12 – what are some potential reasons for this?
Response: Thank you for pointing out this phenomenon. The possible causes are primarily as follows. First, in the 12th year, the site was in the mid-transition phase of recovery from a thermokarst collapse. Vegetation recovery, microbial community reorganization, and substrate supply may not yet have reached a stable state. As a result, there was high spatial heterogeneity in microtopography, soil moisture, pH, organic carbon input, and microbial biomass within the study plots, leading to significant differences in rapid storage capacity among the three replicate plots. Second, C0,fast is a fitting parameter in the two-pool model; its estimation depends on the shape of the cumulative mineralization curve during the early stages of cultivation and involves a certain trade-off with kfast. If the mineralization rates or curve shapes differ among the three replicate samples during the early stages of cultivation, the model’s allocation of fast pool capacity may vary significantly. Third, the small sample size and high variability among the three replicate sites increase the standard error, which may also be one of the reasons why the main effect of recovery years on C0,fast was not significant (p = 0.056).
In the revised manuscript, we will examine the original cumulative mineralization curves and fitting parameters for these three replicate samples to determine whether there are any individual outliers or issues with the identifiability of the model parameters. If necessary, we will provide the original curves for the replicate samples, the fitting results, and a sensitivity analysis after removing outliers in the supplementary materials. Finally, it should be noted that although C0,fast exhibited considerable variability in the 12-year topsoil data, the core conclusion of this study—that kfast rebounds significantly with recovery duration and is significantly higher in the 23-year data compared to the 1-year and 12-year data—remains robust after parametric tests, nonparametric tests, and AR(1) correction. Therefore, this variability does not affect the main conclusions. Thank you again for raising this important issue.
-Fig 4: Why only these parameters in the PLS-SEM figure, rather than including all those showing high impacts in Fig 3b?
Response: Thank you for raising this question. It is important to note that LMG and PLS-SEM serve different purposes. LMG is a variance decomposition method used to quantify the relative importance of each variable among multiple candidate variables for kfast; therefore, it allows for the inclusion of as many candidate indicators as possible; In contrast, PLS-SEM is a theory-driven path model used to test the causal pathways of recovery years, organic carbon chemical quality, and microbial attributes on kfast. Consequently, variable selection must be based on a priori assumptions, conceptual consistency, and model parsimony; it is not possible to simply include all high-contribution variables from the LMG analysis in the path diagram. First, the latent variables in the current Fig. 4 were constructed based on the following considerations. “Recovery Time” is an ordinal code for the duration of recovery, representing the time required for thermal slump recovery. “Organic Chemistry” represents the chemical quality of organic carbon and consists of O-alkyl C, Alkyl C, and the A/O ratio, as these three indicators directly reflect the composition and degree of decomposition of carbon functional groups as determined by NMR, consistent with the theoretical concept of substrate quality. Microbial Biomass represents microbial biomass, community structure, and enzyme activity, and is composed of Total PLFA, the F/B ratio, and BG, as PLFA reflects microbial biomass and taxonomic structure, while BG may reflect enzyme activity related to readily degradable carbon. kfast uses the model-fitted values directly. This path structure corresponds to our previously proposed hypothesis: restoration duration directly influences the chemical composition of organic carbon and microbial properties, both of which further influence kfast, while restoration duration also has a direct effect on kfast.
Second, to verify whether the PLS-SEM in Fig. 4 omitted any important variables, we conducted further validation. Collinearity diagnostics for the original variables indicated that these variables could not simultaneously be included in a stable path model. We calculated the VIF for all candidate raw variables: TN was 2089.65, SOC was 1970.51, Alkyl/O-alkyl was 1897.52, Alkyl C was 1421.65, C/N was 615.94, NAG was 165.70, BG was 142.30, O-alkyl C was 131.80, F/B ratio was 7.69, and Total PLFA was 6.92, while DOC was the only variable with a VIF less than 5 (VIF = 2.42). This indicates severe multicollinearity among SOC, TN, C/N, O-alkyl C, Alkyl C, A/O ratio, NAG, and BG; if all were included in the PLS-SEM model, the path coefficients would be extremely unstable or even impossible to interpret. Therefore, in Fig. 4, we adopted a latent variable strategy: “Recovery Time” represents the number of years required for recovery; “Organic Chemistry” consists of O-alkyl C, Alkyl C, and the A/O ratio, representing the chemical quality of organic carbon; “Microbial Biomass” consists of Total PLFA, the F/B ratio, and BG, representing microbial biomass, community structure, and enzyme activity; and “kfast” uses the model-fitted value directly. The maximum Inner VIF value for this model is 2.277, and there is no severe multicollinearity among all latent variables, indicating that the SEM structure is valid. The three indicators under “Organic Chemistry” in the Outer model have relatively high VIF values because they are highly correlated with each other; however, PCA has effectively reduced the dimensionality by extracting the first principal component. The VIF values for the microbial biomass indicators are all low (Total PLFA = 1.31, BG = 1.22, F/B = 1.18), indicating no multicollinearity issues.
We subsequently conducted three robustness tests. First, we added SOC, TN, C/N, DOC, and NAG as additional control variables to the original regression model. The F-test for ΔR² of the regression model with SOC, TN, C/N, DOC, and NAG added was not significant compared to the original regression model (F(5,3) = 1.7511, p = 0.3419), indicating that these additional variables did not significantly improve the model and that the main model did not exhibit obvious omitted variable bias. It should be noted that the degrees of freedom for residuals in this test were only 3, limiting its statistical power; therefore, we did not treat this as the sole evidence but instead made a combined judgment based on the subsequent extended PLS-SEM and the DOC-specific sensitivity analysis.
Second, we constructed an extended PLS-SEM model. Building upon the original variables—Recovery Time, Organic Chemistry, Microbial Biomass, and kfast—we added two latent variables: Soil C, N status (including SOC, TN, CN, DOC) and Enzyme Activity (including BG and NAG) to examine whether the omission of important variables has an independent effect on kfast. The maximum Inner VIF value in the extended model was 2.572, which is below the commonly used threshold of 5, indicating that there is no severe multicollinearity among the latent variables. The results show that the path coefficients of Recovery Time, Organic Chemistry, and Microbial Biomass on kfast are β = +0.697 (p = 0.0053), β = −0.463 (p = 0.0185), and β = −0.340 (p = 0.1280), respectively; both the direction of their effects and their significance levels are consistent with those of the main model. However, none of the newly added paths were significant; among them, the effects of Soil C, N status had an effect on kfast of β = −0.380 (p = 0.0953), and Enzyme Activity had an effect of β = −0.168 (p = 0.4539). This indicates that even when important omitted variables are included in the model, the direct effect of recovery time on kfast remains stable, and the newly added variables do not alter the core conclusions.
Third, considering that DOC makes a relatively high independent contribution in the LMG and is the only candidate variable among the original variables with a VIF less than 5 (VIF = 2.42), we further conducted a specific sensitivity analysis with DOC as an independent latent variable. The maximum Inner VIF of this model was 2.575, which is below the multicollinearity threshold. The results showed that the effect path of DOC status on kfast was not significant (β = +0.515, p = 0.1547), and it did not alter the effects or significance levels of Recovery Time, Organic Chemistry, and Microbial Biomass on kfast. This indicates that although DOC is correlated with kfast to some extent, it does not exhibit a significant independent direct effect in the path model and cannot replace the explanatory role of recovery years or organic carbon chemical quality in explaining kfast. A comprehensive analysis of the three robustness tests reveals that, although variables such as SOC, TN, C/N, DOC, and NAG show some independent contribution in the LMG variance decomposition, they are highly correlated with the organic carbon chemical quality and microbial variables in the current PLS-SEM model. Furthermore, the original variables have high VIF values and cannot be simultaneously included in a stable PLS-SEM model. After incorporating them into the extended model, they neither significantly improved explanatory power nor altered the current core pathways. Therefore, the variable selection in the current PLS-SEM is reasonable, and there is no apparent omission bias. We will provide further justification for this selection in the revised manuscript and include the complete results of the original variable VIF diagnostics, control variable regression, the extended PLS-SEM, a sensitivity analysis specifically for DOC, and Inner/Outer VIF calculations in the supplementary materials. Thank you once again for raising this important point.
-Line 338 – 339: Authors note soil moisture decreased with restoration duration – was this topsoil, subsoil, or an average of both? Were potential interaction effects of soil moisture on microbial community at different stages of restoration also assessed?
Response: Thank you for raising this question. Regarding the statement in lines 338–339 that “soil moisture decreases with restoration time,” it should be clarified that this result was observed separately in the topsoil and subsoil, and is not an average of the two. According to the source literature, Mu et al. (2024), compared to the control, soil moisture in the topsoil decreased by approximately 29% and in the subsoil by approximately 25% after 23 years of restoration, with both depths exhibiting a decreasing trend over the restoration period. We will explicitly note this in the revised manuscript to prevent readers from mistakenly assuming that the result represents an average across depths. Regarding the potential interaction between soil moisture and microbial communities at different restoration stages, we must candidly state that we have not specifically assessed this interaction in the current analysis. Our supplementary analysis shows that soil moisture is not significantly correlated with kfast (Table S11); therefore, moisture was not included as an independent latent variable or interaction term in the main PLS-SEM analysis. The structural equation model by Mu et al. (2024) also indicates that moisture primarily influences CO₂ release indirectly by affecting organic carbon mass, without directly testing the interaction between moisture and the microbial community. Therefore, we did not include it in the main model. Finally, we would like to thank you once again for raising this important point.
-Line 405 and thereafter: Much of this section raises good, probable points but seems to be largely speculation as it does not include citations. Citations should be added here where applicable, and if no evidence exists, language should be adjusted to make clear that these statements are speculative.
Response: Thank you for bringing this issue to our attention. We fully agree with your assessment: while the discussion section beginning on line 405 proposes many reasonable and plausible mechanisms, some of the statements do indeed lack direct citations or are phrased too definitively without sufficient supporting evidence, thereby taking on a speculative nature. In the revised manuscript, we will review this section paragraph by paragraph, adding citations where appropriate and adjusting the wording of statements lacking direct evidence to clearly indicate their speculative nature, thereby avoiding any misinterpretation by readers that these mechanisms have been fully substantiated by the data from this study or by existing literature. For example, we will cite Liang et al. (2017) and Kallenbach et al. (2016) to illustrate the contribution of microbial residues to alkyl carbon signals, and supplementing Cotrufo et al. (2013) and Sokol et al. (2019) to support the mechanisms linking microbial residues to the formation of mineral-bound organic matter.
Kallenbach, C. M., Frey, S. D., Grandy, A. S.: Direct evidence for microbial-derived soil organic matter formation and its ecophysiological controls. Nature communications., 7, 13630 , https://doi.org/10.1038/ncomms13630, 2016.
Liang, C., Schimel, J. P., Jastrow, J. D.: The importance of anabolism in microbial control over soil carbon storage. Nature microbiology., 2, 1-6, https://doi.org/10.1038/nmicrobiol.2017.105, 2017.
Cotrufo, M. F., Wallenstein, M. D., Boot, C. M., Denef, K., Paul, E.: The Microbial Efficiency-Matrix Stabilization (MEMS) framework integrates plant litter decomposition with soil organic matter stabilization: do labile plant inputs form stable soil organic matter? Global Change Biology., 19, 988-995, https://doi.org/10.1111/gcb.12113, 2013.
Sokol, N. W., Sanderman, J., Bradford, M. A.: Pathways of mineral-associated soil organic matter formation: Integrating the role of plant carbon source, chemistry, and point of entry. Global Change Biology., 25, 12-24, https://doi.org/10.1111/gcb.14482, 2019.
-Lines 469 – 471: Quantifying these changes in pH within the text here would be helpful. What was the pH directly after the disturbance? What was the rate of change across years? Was the change linear, and synchronous among stages / top or subsoil?
Response: Thank you for this suggestion. We agree that simply describing the trends in pH changes is insufficient, and we will supplement the revised manuscript with quantitative pH results based on the existing dataset. Specifically, this will include: the mean pH values, standard errors, and sample sizes for the topsoil and subsoil at each restoration stage (control, 1 year, 12 years, and 23 years); the magnitude of the pH decline relative to the control one year after disturbance; the pH changes and average annual rates of change for the periods 1–12 years, 12–23 years, and 1–23 years; and linear and quadratic regressions for topsoil and subsoil, respectively, to test whether pH changes linearly with restoration duration. The current analysis shows that pH does not change linearly with recovery duration but rather exhibits a quadratic trend, initially decreasing and then rebounding after disturbance. The recovery duration effect is significant (p = 0.031), and the depth effect is significant (p = 0.027); the direction of change in topsoil and subsoil is generally synchronized, but there are depth-dependent differences in absolute levels and the magnitude of change. We will add these specific values and statistical results near lines 469–471 in the main text and provide the complete pH data table and regression results in the Supplementary Materials. Thank you again for raising this important point.
-In lines 498 – 503 authors suggest adding a substrate quality adjustment function for fast pool rate constants based on A/O, but within this research they indicated that the general trend of higher A/O correlating with lower decomposition rates did not hold, potentially influenced by vegetation carbon inputs and microbial residues. How would this adjustment function take this possible deviation into account?
Response: Thank you for raising this critical question. The contradiction you point out is indeed significant: this study observed a positive correlation between the A/O ratio and kfast, contrary to the conventional assumption that “a higher A/O ratio leads to slower decomposition.” This is likely related to vegetation carbon input and the accumulation of microbial residues during the late recovery phase. Therefore, the substrate mass adjustment function we propose should not be understood as a formula with a fixed exponent that can be directly calibrated and applied; rather, it should be viewed as a conceptual framework in which the substrate mass sensitivity index γ itself can vary depending on the stage of recovery, carbon source, and the contribution of microbial residues.
Specifically, the adjustment function we propose is kfast=kbase×(A/O)γ, where kbase is the base rate controlled by physical factors such as temperature and moisture, and γ is the substrate mass sensitivity index. Under the traditional assumption, γ is typically negative, meaning that an increase in A/O reduces the decomposition rate. However, this study indicates that γ cannot simply be set as a constant or a fixed negative value. In the early stages of recovery, an increase in A/O may primarily reflect the accumulation of inert alkyl carbon of plant origin, at which point γ may still be negative; whereas in the later stages of recovery, an increase in A/O may be primarily driven by microbial residues, a significant portion of whose alkyl carbon can still be rapidly utilized under cultivation conditions, at which point the effective γ may approach zero or even be positive. Therefore, this adjustment function needs to allow γ to vary with the stage of recovery, or introduce a “alkyl carbon availability” correction factor ε, representing the proportion of rapidly mineralizable alkyl carbon in the A/O ratio; for example, it can be expressed as kfast=kbase×(A/O×ε)γ, where ε increases in the later stages of recovery as microbial residues accumulate. In this way, the function can accommodate both the traditional scenario where an increased A/O ratio inhibits decomposition and the scenario observed in this study, where an increased A/O ratio is accompanied by a rise in kfast.
It should be clarified that the data from this study are insufficient to directly calibrate γ or ε, nor can they distinguish the relative contributions of plant-derived inert alkyl carbon and microbial-derived active alkyl carbon in the A/O ratio. Therefore, in the revised manuscript, we will explicitly state that this adjustment function is a conceptual proposal for future Earth system models, rather than a parameterization scheme validated in this study; its core significance lies in highlighting that models should not assume a fixed relationship between substrate mass and decomposition rate, but should instead account for the modulatory effects of the recovery stage, vegetation carbon input, and microbial residues on substrate mass. Future research will need to combine microbial residue markers such as amino sugars, isotope tracing, and more intensive restoration sequences to quantitatively distinguish alkyl carbon from different sources and further constrain the range of variation for γ or ε. We will adjust the relevant wording in the revised manuscript to clarify this point and avoid any misunderstanding by readers that the function can be directly applied or has already been calibrated.
Line comments:
-Line 149: There seems to be an errant period here where there should be a comma.
Response: Thank you for pointing out this punctuation error. In the revised version, we have replaced the period with a comma in that instance to conform to the article’s stylistic guidelines.
-Lines 167 – 168: The “fast” in “kfast” is not written as a subscript.
Response: Thank you for pointing this out. In the revised draft, we will uniformly check the formatting of “kfast” throughout the text to ensure that subscripts are used in all instances, including the main text, figure captions, tables, and supplementary materials.
-Lines 182 - 187: Tone shifts to present tense when rest of text is in past tense.
Response: Thank you for pointing out the tense issue in the article. We will standardize the tense to the past tense in the revised version.
-Line 251: The semicolon should be a period.
Response: Thank you for pointing out the punctuation issue. We will review the sentence structure to ensure clarity.
-Line 259: The semicolon should be a period.
Response: Thank you for pointing out the punctuation issue. We will review the sentence structure and, at the same time, check the grammar and logical relationships between adjacent sentences.
-Line 306: The semicolon should be a period.
Response: Thank you for pointing out the punctuation issue. We will review the sentence structure to ensure clarity.
-Line 355: The CO,fast are not written as subscript.
Response: Thank you for pointing this out. We will correct it to C0,fast in the next revision and ensure the subscript format is correct.
-Line 358: Kfast not written as subscript.
Response: Thank you for pointing this out. We will correct it to “kfast” in the next revision and ensure the subscript format is correct.
-Line 486-487: Seems like the sentence is incomplete or needs to be re-worded/joined with the next sentence.
Response: Thank you for pointing this out. This sentence does indeed contain repetitions and is incomplete. We will combine and rewrite it in the revised draft.
-Table S14: Indicates a significant (or borderline) effect of soil depth (k upper bound = 0.5) though is said not to be significant in text.
Response: Thank you for pointing out this inconsistency. We have reviewed Table S14 and the relevant statements in the main text and confirm that the wording was indeed inaccurate. Table S14 shows that when the upper limit of kfast is set to 0.5 d⁻¹, the soil depth effect is marginally significant (Depth p = 0.045); however, under the baseline upper limit of 1.0 d⁻¹ and the more lenient upper limits of 2.0 and 5.0 d⁻¹, the depth effect was not significant (Depth p-values of 0.218, 0.413, and 0.445, respectively). In contrast, the recovery years effect was significant under all upper limit settings (Years p = 0.028, 0.003, 0.003, and 0.008, respectively). Therefore, the general statement in the main text that “depth is not significant” is inaccurate and should be qualified as “the depth effect is not significant under the baseline upper limit of 1.0 d⁻¹ and the more lenient upper limits.”
We further examined the results of the parameter upper-limit sensitivity analysis and found that when the kfast upper limit was set at 0.5 d⁻¹, approximately 21.7% of the kfast estimates for the samples reached this limit, indicating that this relatively strict upper limit imposed a constraint on some samples and may have distorted the estimates of the depth effect; whereas under the baseline upper limit of 1.0 d⁻¹, the proportion reaching the limit was 17.4%, 13.0% under limits of 2.0 and 5.0 d⁻¹, and 0% when no constraint was applied. Therefore, the marginal significance of the depth effect at 0.5 d⁻¹ is more likely an artifact caused by the upper limit constraint rather than a robust ecological signal. Based on this, we will make the following revisions in the revised manuscript: First, in the main text, we will no longer state in general terms that the depth effect is not significant, but will explicitly state that the depth effect is not significant under the baseline upper limit of 1.0 d⁻¹, but is marginally significant under the stricter upper limit of 0.5 d⁻¹, and that the depth effect is sensitive to parameter upper limits. Second, in the footnote for Table S14, we will supplement the text with the proportion of samples reaching the upper limit under each threshold, note that the depth effect at 0.5 d⁻¹ should be interpreted with caution, and point out that the effect of recovery duration on kfast is robustly significant under all thresholds, whereas the depth effect is unstable and does not affect the core conclusion of this study regarding recovery duration driving the rebound of kfast.
Citation: https://doi.org/10.5194/egusphere-2026-3656-AC2
-
AC2: 'Reply on RC2', Hanhan Li, 29 Sep 2026
Viewed
| HTML | XML | Total | Supplement | BibTeX | EndNote | |
|---|---|---|---|---|---|---|
| 293 | 114 | 42 | 449 | 40 | 30 | 29 |
- HTML: 293
- PDF: 114
- XML: 42
- Total: 449
- Supplement: 40
- BibTeX: 30
- EndNote: 29
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1
Dear editor and authors,
General Comments
Thank you for this interesting manuscript. The manuscript concerns a re-analysis of a dataset recently published in JGR Biogeosciences (https://doi.org/10.1029/2024JG008162). The main novel contribution of this current manuscript concerns the fitting of single and two pool model to the CO2 production during an incubation from Mu et al 2025. The evaluation of turnover time of CO2 produce by thermokarst soils is highly relevant for this journal. However, the manuscript currently falls short on its scientific novelty and quality. Below, I explain this further.
Specific comments
The main finding of the current manuscript describes that the total size of fast carbon does not change along the thaw slump gradient, while its turnover time (kfast) does. The authors state that C0 describes the capacity of the fast/slow carbon. In this decay model, C0 describes the initial decomposition rates during the incubation. Throughout the manuscript the authors use C0 to infer the size (or ‘capacity’) of the carbon pool and state that the size of the carbon does not change, while the actual fresh carbon pool has not been measured. More importantly, the soil along the slum age gradient produce less CO2 over a 189 day incubation than the control. How can the same size of fresh carbon pool produce more CO2? This needs to be discussed. Furthermore, any discussion on the lack of difference in parameters for the slow is currently missing but needs to be added.
Most importantly, the authors do not discuss the results and interpretation of the original publication, which state “We found thaw slumps led to a large loss of labile C, and the resistant C accumulated gradually along the thaw sequence”. In this current manuscript the authors claim the opposite, as discussed above, without discussing the difference in interpretation of the same data. This crucial difference needs to be addressed in order to substantiate their conclusions.
Technical comments
The main results from this manuscript rely on the analysis of the key parameters (C0,fast; kfast; C0,slow; kslow; C0,total) extracted from non-linear model fitting. In its current description, there are multiple biases in the performed ANOVA on those key parameters, that make interpretation of the statistical results unreliable and will need to be further improved.
Furthermore, I am also concerned about the analysis and interpretation of the data. For example;
As it stands, I do not see fit for publication of this manuscript for Biogeosciences. Below I state more line comments that need addressed for this manuscript.
Minor comments
The authors speak of recovery from the thaw slump event throughout the manuscript. However, the photos in the original manuscript show clear physical disturbance even 23 years after the occurrence of the thaw slump event, so it is unclear if recovery is ongoing or achieved at these sites. Furthermore, it is unclear what type of recovery is then meant. For example, no description of the plant height or composition is given. The original publication rather speaks of a evolution.
The presented pearson correlations highly overlap with the original manuscript, thereby do not provide new data and should be removed.
The authors should explain in the introduction what is meant by active carbon. Which physical/chemical parameters determine whether carbon is active?
Line comments
290 please change “land-use” into an adequate term to describe the thaw slump age gradient
405 “In the late stages of recovery, the continuous production and selective preservation of microbial residues can enhance alkyl carbon signals;” This is speculation and not shown in the reference. Please remove or adjust sentence accordingly.
Figure 1. Labels on axes are hard to read. Colours are not explained.
Figure 2. Error bars and mean points need to be explained
Notes S1. “whereas the effects of soil depth and their interaction were not significant (Table S4).” Please update this text to reflect the P-values shown in Table S14, where it shows a significant effect.