This study was conducted at the Haibei Research Station of the Chinese Academy of Sciences (37°29′N, 101°12′E; 3,220 ma.s.l.), China. This region has a typical continental semi-humid plateau climate with long cold dry winters and short warm wet summers. The average annual temperature is −1.7°C and ranges from −15.2°C in January to 9.9°C in July. The average annual precipitation is 561 mm, and more than 80% of the annual precipitation is concentrated during the growing season. The growing season in this region is from May to August. The vegetation of this region is dominated by graminoids and forbs, such as Kobresia humilis Sergiev., Stipa aliena Keng., Elymus nutans Griseb., Poa pratensis L., Oxytropis ochrocephala Bunge., Medicago ruthenica (L.) Trautv., and Gentiana straminea Maxim.

The field experiment was established in July 2011 in a fenced 50 × 50 m² flat area that was previously used as a winter pasture to examine how the alpine meadow ecosystem responds to gradient warming. Twenty-five plots composed of five temperature treatments (ambient temperature and four warming treatments) with five repetitions were distributed into subplots divided into five rows and five columns based on a random block design (Supplementary Figure S1A). These four warming treatments were achieved with open-top chambers that were designed as four types of sample rounds from large to small, corresponding to the magnitude of warming from weak to strong (Supplementary Figure S1A). These OTCs were constructed with 1 mm thick fiberglass and were 40 cm high. The diameters of the bottom/top of the chambers were 2.05 m/1.60 m for W1, 1.75 m/1.30 m for W2, 1.45 m/1.00 m for W3, and 1.15 m/0.70 m for W4. During the experiment, all the open-top chambers were maintained at the plots throughout the year (Supplementary Figure S1). Decagon ECT sensors (Decagon Devices, Pullman, WA, United States) were used to record the daily air temperature at 20 cm high and 10 cm aboveground and soil temperature at a depth of 15 cm belowground. After monitoring, these four types of OTC treatments (W1–W4) successfully produced increased soil temperature gradients. The mean increase in soil temperature at the beginning of growing season (May) was 0.73°C, 1.60°C, 1.73°C, and 2.25°C in W1 to W4, respectively, while the mean increase in soil temperature at the end of growing season (August) was 0.52°C, 0.97°C, 1.49°C, and 1.57°C, respectively (Supplementary Table S1). Compared with the beginning of the growing season (May), the soil temperature at the end of the growing season (August) increased by 5.82°C (W0), 5.61°C (W1), 5.19°C (W2), 5.58°C (W3) and 5.14°C (W4; Supplementary Table S1).

Soil samples were collected at the beginning (on May 3) and end of the growing season (August 28) in 2019, respectively. Three soil cores (5.0 cm in diameter and 25 cm deep) were randomly collected in each plot and homogenized into one composite sample (2 seasons × 5 treatments × 5 repetitions = 50 soil samples). Each soil composite sample was then divided into two subsamples. One was used to extract DNA, and the other was used to determine the soil physicochemical and microbial properties. All the soil samples were stored at 4°C during the sampling process (within 48 h) and were transferred into −20°C storage once they had been transported to the laboratory. The plant community survey was only performed at the end of growing season because it was difficult to identify plant species when the plant community had just turned green at the beginning of growing season. The plant community in each plot was surveyed using a 50 × 50 cm² sample quadrat. All the shoots in each quadrat were cut at the ground surface and grouped by plant species. The shoots of each species were then weighed after drying at 80°C for 48 h, and aboveground biomass of plant community was calculated as sum of shoot biomass of each species.

The soil moisture was measured gravimetrically after drying at 105°C for 12 h. The soil pH was measured with a pH electrode (Sartorius PB-10; Sartorius AG, Göttingen, Germany). SOC, dissolved organic carbon (DOC), and nitrogen were measured using an Elementar analysis system (CHNS; Elementar Analysensysteme GmbH, Langenselbold, Germany; Cabrera and Beare, 1993; Jones and Willett, 2006; Makarov et al., 2013). Both soil available nitrogen (as ammonia and nitrate) and available phosphorus were measured using FIAstar analyzer 5,000 (FOSS, Hillerød, Denmark). The rate of soil respiration for each warming treatment was measured 10 days per month with an LI-8100A Automated Soil CO₂ Flux System (LICOR, Inc., Lincoln, NE, United States) throughout the growing season (from May to August; Ren et al., 2021). SOC storage was calculated based on the SOC, bulk density, and soil sampling depth (Emde et al., 2021). Both microbial biomass carbon (MBC) and microbial biomass nitrogen (MBN) were measured using the chloroform fumigation direct extraction method (Brookes et al., 1985; Vance et al., 1987). The activities of four soil extracellular enzymes that decompose SOC, including the labile carbon decomposing enzymes amylase and sucrase and the recalcitrant carbon decomposing enzymes cellulase and peroxidase, were assayed as described by Xun et al. (2018), Yin et al. (2014), and Jiang et al. (2021).

Total DNA was extracted from 0.5 g of each soil sample using an E.Z.N.A.® Soil DNA Kit (Omega Bio-tek, Norcross, GA, United States) following the manufacturer’s instructions. The primer pair combination of ITS1F (5'-CTTGGTCATTTAGAGGAAGTAA-3')-ITS2R (5'-GCTGCGTTCTTCATCGATGC-3') was then used to determine the number of gene copies of fungal DNA (ITS), while the primer pair combination of 338F (5'-ACTCCTACGGGAGGCAGCAG-3')-806R (5'-GGACTACHVGGGTWTCTAAT-3') were used to determine the number of gene copies of bacteria (16S rRNA) using a QuantStudio™ 6 Flex Real-Time PCR System (Applied Biosystems, Waltham, MA, United States; Majorbio Bio-Pharm Technology Co., Ltd., Shanghai, China). The target fragment of the primer combination ITS1F-ITS2R was ~300 bp that targeted the ITS1 region (White et al., 1900), while the target fragment of the primer combination 338F-806R was ~468 bp in the V3-V4 region (Mori et al., 2014). The primers used for high-throughput sequencing of fungi (ITS) and bacteria (16S rRNA) were the same as the primers for qPCR. The amplification with primer pair ITS1F-ITS2R was conducted in a 20 μl reaction that consisted of 10 ng for the total DNA template, 0.8 μl of each primer (5 μM), 2 μl 10 × TransStart FastPfu buffer, 2 μl dNTPs (2.5 mM), and 0.2 μl TaKaRa rTaq DNA Polymerase. The thermal cycling conditions were as follows: 95°C for 3 min, 27 cycles of 95°C for 30 s, 55°C for 30 s 72°C for 45 s, and 72°C for 10 min. The amplification with primer pair 338F-806R was similar to those of the ITS1F-ITS2R amplifications with the exception of the use of 0.4 μl TransStart FastPfu DNA Polymerase. The PCR reactions were performed in triplicate. The sizes of both fungal and bacterial amplicons were confirmed separately using 2% agarose gel electrophoresis, purified using an AxyPrep DNA Gel Extraction Kit (Axygen Biosciences, Union City, CA, United States), and quantified using a Quantus™ Fluorometer (Promega, Madison, WI, United States). Purified amplicons were pooled in equimolar amounts and paired-end sequenced on an Illumina MiSeq PE300 platform (Illumina, San Diego, CA, United States) by Majorbio Bio-Pharm Technology Co., Ltd. (Shanghai, China) using a NEXTFLEX® Rapid DNA-Seq Kit (Bioo Scientific Corporation, Austin, TX, United States) for library construction. All the raw sequencing reads were quality-filtered by FASTP (ver. 0.20.0; Chen et al., 2018) and merged by FLASH (ver. 1.2.7; Magoc and Salzberg, 2011). Operational taxonomic units (OTUs) with a 97% similarity cutoff (Stackebrandt and Goebel, 1994; Edgar, 2013) were clustered using UPARSE (ver. 7.1; Edgar, 2013). The taxonomy of each representative OTU sequence was analyzed using RDP Classifier (ver. 2.2; Wang et al., 2007) against the ITS database (UNITE ver. 7,¹ for fungal OTUs) and the 16S rRNA database (Silva ver. 138,² for bacterial OTUs) using a confidence threshold of 0.7. The sequence raw datasets in this study were deposited in the NCBI Sequence Read Archive (SRA) with accession numbers: PRJNA719774 (ITS-BioProject) and PRJNA719776 (16S rRNA-BioProject).

The statistical analyses were primarily conducted using R 4.0.5.³ A linear mixed-effects (LME) model (“lme” model) was applied to analyze the individual effects of warming and season on soil, plant, and microbial variables using the “nlme” and “nlme” functions in R (Pinheiro et al., 2019), in which both warming and season were treated as fixed effects and block as a random effect. Significant differences between treatments were tested using the post-hoc multiple comparisons Tukey method at a 95% significance level (p ≤ 0.05). A rarefaction curve of microbial species richness was estimated using the “rarecurve” function of the “vegan” package in R (Oksanen et al., 2019). A taxonomy Circlize of microbial community composition at the phylum level was summarized using the “tax_circlize” function of the “amplicon” package in R.⁴ Dissimilarities of the plant or microbial community were depicted with principal coordinate analysis (PCoA) ordination based on the Bray–Curtis dissimilarity index using the “cmdscale” function of the “stats” package in R (Oksanen et al., 2007). A permutational multivariate analysis of variance (PERMANOVA) was used to determine the effects of warming or season on both plant and microbial communities using the “adonis” function of the “vegan” package in R (Oksanen et al., 2019). Linear discriminant analysis (LDA) effect size (LEfSe ver. 1.0)⁵ was performed to assess significant differences of the fungal and bacterial taxa from Kingdom; Domain, not Kingdom, was used as the basic classification level for bacteria to genus levels (Segata et al., 2011).

The natural logarithm ratio was used to interpret the rate of change of all the variables measured during the growing season. For example, the rate of change of SOC during the growing season was calculated as ΔSOCcontent%=elnSOCAugSOCMay−1×100, where SOC_Aug and SOC_May represent the content of SOC in August and May, respectively (Emde et al., 2021). Considering that the score values of the PCoA axis of the plant or microbial communities cannot be characterized by the natural logarithm owing to a negative value, we applied the ratio of γ-diversity to α-diversity to characterize the rate of change of the composition of plant or microbial community during the growing season using the “diversity” function of the “vegan” package in R (Oksanen et al., 2019).

Linear regression analysis was performed to assess the relationships between the content of SOC or the rates of change of SOC with all plant, soil, and microbial variables using the linear mixed-effects regression with block as a random effect. A random forest test was performed to predict response variables that are significantly related to the content of SOC (or the rate of change of SOC) using the “rfPermute” function of the “rfPermute” package in R (Archer, 2013). A structural equation model (SEM) was also constructed using SPSS AMOS 21.0 (IBM, Inc., Armonk, NY, United States).

Warming had no significant effect on all the soil variables measured, but the season had significant effects on other soil physical and chemical properties except for soil ammonium nitrogen and nitrate nitrogen (Figure 1). The soil moisture and available phosphorus in May were lower than those in August, while the contents of SOC, DOC, dissolved organic nitrogen, and soil carbon storage in May were higher than those in August (Figure 1). Moreover, warming significantly reduced the aboveground biomass and species richness of the plant communities but did not affect the species composition of these communities (Supplementary Figure S2).

Except for W1, both the MBC and MBN of May were higher than those of August (Figures 2A,B). There was no significant difference in the ratio of MBC/MBN and MBC/SOC among the treatments (Figures 2C,D). The rates of soil respiration in May were lower than those in August under all the warming treatments (Figure 2E). The activity of amylase in May was lower than that in August, and the activity of sucrase in May was higher than that in August under the W0 treatment (Figures 2F,G). The activity of cellulase in both W0 and W3 in May was lower than that in August, while there were no significant differences in the activity of peroxidase between May and August under all the warming treatments (Figures 2H,I).

Warming had no significant effect on the fungal (ITS) gene copy number (Figure 2J), bacterial gene copy number (Figure 2K), and the ratio of fungal/bacterial gene copy number (Figure 2L), while seasons significantly affected the gene copy number of fungi (Figure 2J). The gene copy number of fungi in May under all the warming treatments (except for W4) was higher than that in August (Figure 2J). Rarefaction curves showed that the high-throughput sequencing effectively identified the species richness of both fungal and bacterial communities (Supplementary Figure S3). A total of 307 genera of 36 classes of fungi were detected in 13 phyla, and 407 genera of 81 classes of bacteria were detected in 36 phyla (Supplementary Table S2). Ascomycetes were the dominant fungal taxon, comprising 65.27% of the total number of reads, while Proteobacteria was the dominant taxon of bacteria, comprising 30.90% of the total number of bacterial reads (Supplementary Figure S4).

Warming had no significant effect on the species richness of fungal and bacterial communities (Figures 3A,B). The season had a significant effect on the richness of fungal species but had no significant effect on the richness of bacterial species (Figures 3A,B). The richness of fungal species was higher in May than that in August (except for W4) under all the warming treatments, while the richness of bacterial species in May was significantly lower than that in August under the W0 treatment (Figures 3A,B). In addition, warming had no significant effect on the taxonomic composition of fungal and bacterial communities, while season had a significant effect on them (Figures 3C,D). The relative abundances of 21 fungal OTUs and 310 bacterial OTUs were higher in May than in August, while the relative abundances of only four bacterial OTUs were lower in May than in August (Supplementary Figure S5).

Soil organic carbon positively correlated with the DOC, MBC, MBN, and sucrase activity but negatively correlated with the soil temperature, aboveground biomass, rate of soil respiration, and cellulase activity (Table 1). The random forest analysis confirmed that MBC, DOC, MBC/MBN, and microbial biomass nitrogen were response variables that highly correlated with changes in the SOC (Figure 4A). The structural equation models also showed that the season, but not warming, mediated changes in the SOC through MBC (Figure 5A).

The rate of change of soil and microbial properties were not significantly different among the warming treatments (Supplementary Figure S6). The rates of change of the SOC were positively related to the rates of seasonal change of DOC, MBC, and MBC/MBN (Table 2). Although warming had no significant effect on the rates of change of SOC, the SOC decreased significantly during the growing season only under the warming treatments (W2, W3, and W4; Supplementary Figure S6). The random forest analysis also confirmed that the rates of change of MBC, MBC/MBN, and DOC were variables that were associated with the change of SOC during the growing season (Figure 4B). The structural equation model that was conducted to measure the causal effect of warming on the change of SOC had also shown that warming failed to alter the change of SOC by the seasonal changes in soil temperature, but the change of MBC had a significant effect on it (Figure 5B).

Regardless of season, warming had no significant effect on the SOC (Figure 1A). This result does not support hypothesis I (H1) and is consistent with the findings of other studies (Yue et al., 2015; Zhang et al., 2016; Jia et al., 2017; Guan et al., 2018; Jia et al., 2019; Wang et al., 2019; Chen et al., 2021) that indicated that warming did not dominate the change of SOC. Unlike the warming treatments, the content of SOC in August was significantly lower than that in May (Figure 1A), and the storage of SOC of August decreased by an average of 18.53 million grams of carbon per hectare (Figure 1I), indicating that alpine meadow soil is a carbon source rather than a carbon pool during the growing season. This result supports hypothesis 2 (H2) and is similar to that of another study in a grassland ecosystem (Yang et al., 2021), suggesting that the season, rather than warming, dominated the change in SOC (Gavazov et al., 2017). Since a study has confirmed that the soil carbon flow is higher than that in the non-growing season, we hypothesize that there is an “interannual cycle” process. In the “interannual cycle” process, the content of SOC changes from a “loss state” of the growing season to an “increase state” of the non-growing season and finally returns to the original level at the beginning of the growing season in the following year. In this study, warming reduced the aboveground biomass of plants, which in turn reduced the return of plant litter to the soil during the non-growing season. Therefore, warming did not cause the loss of SOC under the current experimental conditions, but it was expected that the alpine meadow SOC could be lost in future long-term climate warming.

Temperature is a key limiting factor for high-altitude ecosystems (Duan et al., 2019). In this study, the negative correlation between the soil temperature and SOC further confirmed that changes in soil temperature induced by warming or seasons could be used to explain the change in SOC (Supplementary Table S1). In particular, the maximum increase in temperature produced by the warming treatment in May and August was 2.25°C and 1.57°C, respectively. This is well below the increase in temperature caused by the season (5.14–5.82°C), which also indicated that the insignificant effect of warming on SOC could be because the increase in soil temperature caused by warming did not reach the temperature threshold of organic carbon change (Feng et al., 2022). Structural equation modeling also confirmed this (Figure 5). However, since the seasonal changes of SOC are likely to be a cumulative process, whether this inference is suitable for this study still merits further exploration.

The rate of change of SOC during the growing season was not significantly related to the rate of change of soil temperature (Table 2). This was inconsistent with the third hypothesis (H3) but confirms the conclusions of other studies (Classen et al., 2015; Yue et al., 2015), which found that the inter-seasonal temperature difference caused by the warming treatment did not impact microbiota, and thus, failed to alter the seasonal dynamic changes of SOC. However, this result is inconsistent with the results of some studies (Xu et al., 2018; Tian et al., 2021). Tian et al. (2021) and Xu et al. (2018) found that elevated soil temperatures indirectly affected the seasonal dynamics of the content of SOC by altering soil nutrition, physical protection, and microbial biomass. Plant biomass (Classen et al., 2015), soil respiration (Tian et al., 2021), and microbiota (Classen et al., 2015; Xu et al., 2018, 2021) are often used to explain the seasonal changes of SOC. In this study, the aboveground biomass gradually decreased with increasing temperature (Supplementary Figure S2A). Since the aboveground biomass of all the treatments was close to zero at the beginning of the growing season, the reduction in aboveground biomass caused by warming also illustrates the decrease in the seasonal difference of aboveground biomass under warming treatments. This also indicates that the seasonal difference of aboveground biomass had decreased as the temperature increased. Thus, we inferred that owing to the small input of plant aboveground carbon, the seasonal difference of SOC should be reduced as the temperature increases. However, this inference was inconsistent with the results of this study (Figure 1A), indicating that it was unreasonable to only use the plant aboveground biomass to predict the changes in SOC. Considering that belowground biomass and root exudates also play an important role in driving the change of SOC (Jansson and Hofmockel, 2020), further research needs to comprehensively consider belowground biomass and root exudates to explore the effect of warming on the seasonal changes of SOC. In addition to plant biomass, soil respiration is highly sensitive to temperature, and an increase in the temperature will activate soil respiration and accelerate its rate (García-Palacios et al., 2021). Since warming led to a gradual decrease in the seasonal temperature difference in this study (Supplementary Table S1), the seasonal difference of SOC will decrease as the temperature increases. However, we had not observed a decrease in the seasonal difference of SOC across a temperature gradient (Supplementary Figure S6), indicating that the seasonal change of SOC in this study cannot also be explained by the seasonal change of soil respiration. The insignificant correlation between the rate of change of soil respiration and that of the SOC effectively supported this (Table 2). Moreover, the decomposition of soil microbiota was another important reason for the seasonal change in SOC (Chen et al., 2015; García-Palacios et al., 2021; Tian et al., 2021). On the one hand, changes in the soil microbial community composition will affect the SOC by changing the activity of soil extracellular enzymes (Tian et al., 2021). In this study, warming did not have a significant effect on the seasonal differences of fungal and bacterial communities and the activities of soil extracellular enzymes. In addition, there was no significant correlation between the rate of change of the soil microbial community and the rate of change of the soil enzyme activity and that of SOC (Table 2). These results indicate that changes in soil microbial communities and soil enzyme activities in this study cannot explain the insignificant effect of warming on the seasonal differences of SOC. Alternatively, MBC, as an internal component of SOC, directly or indirectly affects the change in content of SOC (Table 2). In this study, the rate of change of MBC was significantly related to the rate of change of SOC. In addition, the SEM models also showed that warming failed to mediate seasonal changes in SOC by altering the changes in MBC (Figures 4B, 5B), indicating that the insignificant seasonal changes in soil temperature on microbial biomass caused by warming in this study could be an important reason for the insignificant effect of warming on the seasonal changes in SOC.

Although MBC comprises a small proportion of SOC (Figure 2D), both the random forest analysis and structural equation model all showed that MBC significantly correlated with the change in SOC, indicating a high degree of coupling of MBC to the change in SOC (Figures 4A, 5A). This is consistent with the findings of another study (Wang et al., 2021). This occurs because MBC, as an internal component of soil carbon (Mcdaniel et al., 2014; Ogle, 2018; Huang et al., 2021), has been regarded to be the liable carbon fraction of SOC that is closely related to the availability of soil nutrients (Guan et al., 2018) and indirectly participates in the key process of the soil carbon cycle (Ma et al., 2020). In addition, our results had shown that neither season nor warming had a significant effect on the ratio of MBC/MBN. This was also confirmed by another study on an alpine meadow ecosystem on the Qinghai–Tibet Plateau, in which the MBC/MBN remained stable under the warming treatment (Wang et al., 2014). However, this result contradicts another study in a subtropical ecosystem that confirmed that warming produces positive effects on MBC/MBN by increasing the decomposition of soil organic matter (Gao et al., 2017). This inconsistency is primarily determined by the climate characteristics of the two types of ecosystems. In particular, the characteristic cold climates of alpine meadows do not facilitate the decomposition of organic matter, but the opposite is true for subtropical climates. In addition to the decomposition of organic matter, the microbial community dominated by fungi (or bacteria) should not be confirmed to be associated with a high (or low) ratio of MBC/MBN (Gao et al., 2021). However, the ratio of l g ITS/l g 16S rRNA gene copy numbers was always approximately equal to 1 under all the treatments, suggesting a co-dominance of fungi and bacteria. This also explained the stability of MBC/MBN ratio under every treatment in this study. However, what was anomalous was that although the MBC/MBN ratio changed owing to warming or season, it was shown to be another variable that was closely associated with the rate of change of SOC. However, both the random forest analysis and the structural equation model also showed that the rate of change of the content of SOC was determined by the rate of change of MBC rather than MBN (Figure 5). Therefore, we inferred that the contribution of MBC/MBN to the change in SOC is just another manifestation of the contribution of MBC and is not really driven by the ratio of MBC/MBN, indicating that the microbiological properties of soil that were related to the carbon cycle were “carbon-related” rather than “nitrogen-related.”

In addition, warming had no significant effect on the content of DOC regardless of season (Figure 1D). This result contradicts those of other studies (Six et al., 2004; Li et al., 2016; Tian et al., 2021). One of these studies found that warming can increase the content of DOC by increasing the dissolved organic compounds produced by plants or microbiota (Bai et al., 2020). However, other studies have found that warming can reduce the organic carbon content dissolved in the soil by promoting the formation of aggregates to solidify the DOC (Six et al., 2004; Li et al., 2016; Tian et al., 2021). In addition, the rate of change of DOC can mediate the rate of change of the SOC by changing the rate of change of MBC during the growing season (Figure 5B). These results are consistent with those of other studies (Zhang et al., 2015; Xiong et al., 2016; Bhople et al., 2019; Deltedesco et al., 2020) that further confirm that the DOC participates in the process of regulation of soil carbon cycle by affecting MBC-derived carbon as described in the carbon cycle model by Fan et al. (2021). However, since the content of DOC makes only minor contributions to the variation in the content of SOC (Chen et al., 2015), its impact on SOC stability was likely to be achieved only through its impact on MBC. Integrating this analysis, MBC, the ratio of MBC/MBN, and DOC are all variables that are significantly associated with the change in the content of SOC, but the effects of both MBC/MBN and DOC are essentially mediated by MBC. There was a high degree of coupling between MBC and soil organic carbon. In addition, MBC is more sensitive to temperature than SOC (Walker et al., 2018). Thus, strengthening the research on the temperature sensitivity of MBC can improve our ability to predict changes in SOC under climate warming.