DOI: 10.5061/dryad.7m0cfxq7n
title: Warming modulates soil multifunctionality through assembly processes and co-occurrence patterns of arbuscular mycorrhizal fungal communities in drylands

description: Climate warming poses a threat to the functionality of global dryland
 ecosystems, while arbuscular mycorrhizal fungal (AMF) communities play a
 critical role in maintaining ecosystem stability and functioning. However,
 it remains unclear how the assembly mechanisms and co-occurrence network
 patterns of AMF communities respond to warming and regulate soil
 multifunctionality. We conducted a three-year warming experiment
 (+0.5－1.6°C) using open-top chambers (OTCs) in the Tengger Desert to
 investigate the assembly processes and co-occurrence network patterns of
 AMF in the rhizosphere soil of dominant shrubs, Artemisia ordosica and
 Caragana korshinskii. By measuring 15 soil functional indicators
 (including nutrient availability, biogeochemical cycles, microbial
 activity, and microbial productivity), we constructed a soil
 multifunctionality index and employed structural equation modeling (SEM)
 to uncover the mechanisms by which AMF communities regulate soil
 multifunctionality. Our results showed that warming increased the network
 complexity of AMF communities associated with C. korshinskii by 8%－145%,
 while decreasing that of the A. ordosica by 14%－80%. Stochastic processes
 dominated the assembly of AMF communities, with warming reducing the
 stochasticity of AMF communities in C. korshinskii while increasing it in
 A. ordosica. Notably, the two species showed distinct functional response
 pathways: For A. ordosica, warming-induced reduction in soil water content
 suppressed the mycorrhizal colonization rate (36%－79%) through a
 "stochastic assembly-network collapse" cascade effect,
 subsequently leading to decreased soil multifunctionality (54%－172%). In
 contrast, for C. korshinskii, the soil water reduction caused by warming
 maintained higher colonization rates (82%－198%) by enhancing AMF network
 complexity, thereby improving soil multifunctionality (43%－228%). We
 further proposed a conceptual framework that integrates niche theory into
 a mechanistic understanding of how stochastic processes and network
 complexity of AMF communities affect soil multifunctionality under climate
 warming. Overall, our study suggests that warming modulates soil
 multifunctionality by affecting the assembly processes and network
 complexity of AMF communities, exhibiting species-specific responses. This
 provides crucial theoretical support for understanding the mechanisms
 underlying functional evolution in drylands and for formulating climate
 adaptation strategies.
description: <strong><em>Study site and field experimental
 design</em></strong> Samples were collected
 at the Shapotou Desert Research and Experiment Station of the Chinese
 Academy of Sciences (SDRES, CAS; 37°32′ N, 105°02′ E), located on the
 southeastern edge of the Tengger Desert in North China. The mean annual
 temperature at an altitude of 1330 m over the past ten years is
 approximately 8.6°C, ranging from mean average temperatures of -6.9°C in
 winter (January) to 24.3°C in summer (July). The mean annual precipitation
 is approximately 186 mm, with 80% falling during the warmer growing season
 (May－September). Long-term observations from SDRES (1955–2016) show that
 the average annual temperature has increased by 0.5°C while the average
 annual precipitation has decreased by 4 mm over the last decade (Li et
 al., 2021; Li et al., 2023). The experimental design of
 this site has been described in detail in previous studies (Li et al.,
 2021; Li et al., 2023). In brief, the open-top chambers (OTCs) were
 deployed in the field in 2019, with their lower 50 cm buried underground.
 Seedlings were transplanted into the OTCs in May 2020, followed by the
 initiation of official monitoring. These experimental setups were
 specifically designed to model the temperature increase in the Tengger
 Desert over the previous decade and the current one. To accomplish this,
 we systematically adjusted the shape and size of the OTCs to simulate the
 targeted warming conditions. Ultimately, we selected four OTCs of
 different sizes from a series of existing devices and used them to
 simulate four treatments in the experiment. The four different sizes of
 OTCs are: a large octagonal column of 130 cm × 130 cm × 200 cm (+0.5°C
 warming); a small octagonal column of 100 cm × 100 cm × 200 cm (+0.9°C
 warming); a large quadrilateral of 130 cm × 130 cm × 200 cm (+1.2°C
 warming); a small quadrilateral of 100 cm × 100 cm × 150 cm (+1.6°C
 warming). In addition, we also set up control groups that were at least 2
 meters away from the OTCs to minimize bias. These control plots replicated
 the OTC subsurface conditions (buried 50 cm deep) and underwent identical
 pretreatment, installation, and recovery protocols, differing only in the
 absence of aboveground OTC structures. Air temperature, precipitation, and
 soil water at 0-5 cm depth in each OTC and the control group were
 continuously monitored throughout the experiment using an automated
 weather station (U30 HOBO data logger, Onset Computer Corporation, Bourne,
 MA, USA). Analysis of the data showed that the mean annual temperature of
 large, medium, and small OTCs increased by 0.5°C, 0.9°C, 1.2°C, and 1.6°C,
 respectively, compared to the control.
 <strong><em>Soil sampling and property
 measurement</em></strong> All field
 sampling for this study was conducted in mid-September 2023, three years
 after the implementation of the OTC. The experimental plants, <em>A.
 ordosica</em> and <em>C. korshinskii</em>, were
 cultivated as follows: seeds were collected from mature parent plants of
 the same species in the study region, screened for uniformity in size,
 shape, and thousand-seed weight, and germinated in a controlled nursery.
 Seedlings were bare-root transplanted into experimental plots in May
 2020 after reaching 7－8 cm in height. Monoculture treatments contained
 either <em>A. ordosica</em> (1 plant/m²) or <em>C.
 korshinskii</em> (2 plants/m²), planted in pre-cleared sandy
 substrate with no pre-existing vegetation. No other plant species were
 introduced into the OTCs during the experiment. In each
 OTC and control treatment, we collected the above-ground parts and
 below-ground root systems of individually planted <em>A.
 ordosica</em> and <em>C. korshinskii</em> plants. Four
 to five plants per treatment were randomly selected for collection, and
 their root systems were excavated to a depth of 50 cm. Rhizosphere soil
 collected from the 20－40 cm soil layer of these plants was thoroughly
 mixed to form a single representative sample. In total, 60 samples (5
 treatments × 2 plant species× 6 replicates) were collected. Specifically,
 during sampling, excess soil was shaken off from the roots by hand,
 leaving approximately 1 mm of soil still attached to the roots. The roots
 were then cut into approximately 5 cm segments using sterile tools, sealed
 in special 50 ml tubes, and immediately placed in liquid nitrogen for
 quick freezing. These roots were transported to the laboratory to isolate
 roots and rhizosphere soil using the methods outlined by Zhong et al.
 (2022). The rhizosphere soil suspension isolated from the roots was stored
 in a 50 ml Falcon tube and refrigerated at -80°C for subsequent DNA
 extraction. In addition, the other fresh soil samples were partly present
 at 4°C to measure microbial biomass and enzyme activity and partly
 air-dried for soil physicochemical analyses. Roots (diameter &lt; 1
 mm) used for analysis of AMF colonization rate were rinsed with distilled
 water in the laboratory and placed in centrifuge tubes.
 Both above-<strong> </strong>and below-ground biomass
 collected were oven-dried (60°C) to a constant weight.
 <strong><em>AMF colonization of
 roots </em></strong> Fine roots (root
 diameter &lt; 1 mm) were separated from the sieved soil, cut into 1 cm
 segments, and washed in a 10% (w/v) potassium hydroxide (KOH) solution for
 30 min at 90°C (Wu et al., 2023). After rinsing the roots with KOH, the
 samples were acidified in 2% HCl (v/v) for 30 min and then immersed in
 0.05% (w/v) trypan blue for 30 min. Finally, the percentage of root length
 colonized by AMF, i.e., the mycorrhizal colonization rate (%), was
 determined using the magnified intersection method (McGonigle et al.,
 1990). <strong><em>Amplicon sequencing and
 bioinformatics analysis </em></strong>
 Following the manufacturer's recommendations, DNA was
 extracted from 0.25 g of frozen soil for each soil sample using the
 PowerSoil® DNA Isolation Kit (MoBio Laboratories, Inc., Carlsbad, CA,
 USA). For Illumina MiSeq sequencing, 18S rRNA gene fragments were
 amplified by a two-step PCR procedure using GeoA-2 and AML2 primers (Lee
 et al., 2008) and NS31 and AMDGR primers (Simon et al., 1992; Sato et al.,
 2005), respectively. The first PCR reaction used a 25 µl reaction mixture
 with 2.5 µl DNA template, 1.25 µl GeoA-2 and AML2 primers (10 µM each), 25
 µl 2 × SYBR Premix Ex Taq (Takara, Japan), and 20 µl sterile water. The
 first-step PCR product was then diluted 10-fold to serve as the template
 DNA for the second PCR reaction. The primer pairs in the second PCR were
 supplemented with the 454 pyrosequencing adapters and the 8-bp barcode for
 multiplexing. PCR reaction conditions were as follows: 95℃ for 15 min;
 five cycles at 42℃ for 30 s, 72℃ for 90 s, 92℃ for 45 s; 35 (first
 reaction) / 20 (second reaction) cycles at 65℃ for 30 s, 72℃ for 90 s, 92℃
 for 45 s; 65℃ for 30 s and with a final extension at 72℃ for 10 mins. PCR
 products were separated on a 2% agarose gel by electrophoresis, purified
 using the AxyPrep DNA Gel Extraction Kit (Axygen Biosciences, Union City,
 CA, USA), and quantified using the Quantus™ Fluorometer (Promega, USA).
 The purified PCR products were pooled at equimolar concentrations and
 subjected to paired-end sequencing on the Illumina MiSeq PE300 sequencing
 platform (Illumina, San Diego, USA) by standard protocols, as conducted by
 Majorbio Bio-Pharm Technology Co., Ltd (Shanghai, China).
 Raw sequences were filtered using Quantitative Insights into
 Microbial Ecology (QIIME) 1.9.0-dev pipeline (http://www.qiime.org)
 (Caporaso et al., 2010). Briefly, sequences with a quality score below 30,
 less than 200 bp in length, containing ambiguous bases, and no valid
 primer or barcode were excluded from further analysis. The presence of
 probable chimeras was examined using UCHIME by the MaarjAM Glomeromycotina
 database as a reference, and all possible chimeric sequences were removed
 from the dataset. Non-chimeric sequences were clustered into different
 operational taxonomic units (OTUs) at a 97% identity threshold using
 USEARCH v.8.0 based on the UPARSE pipeline (Edgar, 2013). Representative
 sequences from each OTU were selected and searched against the National
 Center for Biotechnology Information (NCBI) nucleotide database and the
 MaarjAM 18S rRNA gene database. The number of sequences per sample was
 normalized to the minimum sample size using the SUB. SAMPLE command in
 MOTHUR to account for the effect of different sequence numbers on the
 analysis of AMF communities. All raw sequencing data have been deposited
 in the NCBI Sequence Read Archive (SRA) under study accession number
 PRJNA988244. <strong><em>Assessment of soil
 multifunctionality </em></strong> We
 evaluated soil multifunctionality using 15 variables related to nutrient
 availability, biogeochemical cycles, microbial activity, and microbial
 productivity (Hu et al., 2021; Qiu et al., 2021). Specifically, these
 variables included total organic carbon (TOC), total nitrogen (TN), total
 phosphorus (TP), nitrate nitrogen
 (NO<sub>3</sub><sup>-</sup> - N), ammonium
 nitrogen (NH<sub>4</sub><sup>+ </sup>- N),
 available nitrogen (AN), dissolved organic carbon (DOC), dissolved organic
 nitrogen (DON), β-1,4-glucosidase (βG), 1,4-β-Dcellobiohydrolase (CBH),
 β-1,4-N-acetylglucosaminidase (NAG), L-leucine aminopeptidase (LAP),
 alkaline phosphatase (ALP), and microbial biomass carbon (MBC), and
 nitrogen (MBN) contents. In addition to obtaining quantitative
 multifunctional indices for each warming treatment, we first standardized
 all individual soil functions using z-score transformations (Hu et al.,
 2021). These standardized z-scores were then averaged to determine the
 soil multifunctional indices under each warming treatment (Maestre et al.,
 2012). Extracellular soil enzyme activities such as β-1,4-glucosidase,
 cellobiose hydrolysis enzyme, β-1,4-N-acetylglucosaminidase, L-leucine
 aminopeptidase, and alkaline phosphatase were determined using fluorometry
 as described in previous studies (Bell et al., 2013). The pH, TOC, TN, AP,
 NO<sub>3</sub><sup>-</sup> - N,
 NH<sub>4</sub><sup>+ </sup>- N, and SWC were
 measured using the Qiu et al. (2021) and Xie et al. (2024) techniques. MBC
 and MBN were determined using the chloroform fumigation extraction method,
 and their conversion coefficients were 0.45 and 0.54, respectively
 (Brookes et al., 1985).
 <strong><em>Statistical
 analysis</em></strong> The AMF
 communities' α-diversity was assessed using the Quantitative Insights
 into Microbial Ecology (QIIME, version 2) platform (Zhu et al., 2024). We
 then averaged the standardized Shannon, Simpson, Chao1, and OTU richness
 indices to obtain an overall α-diversity index for the AMF communities of
 <em>A. </em><em>ordosica</em> and <em>C.
 korshinskii</em>. The β-diversity of the AMF communities was
 calculated based on the Bray-Curtis dissimilarity between samples using
 the ‘adonis’ function in the ‘vegan’ package. The overall differences in
 AMF communities under each warming treatment were determined using
 constrained principal coordinate analysis (CPCoA), analysis of similarity
 (ANOSIM), multiresponse permutation procedures (MRPPs), and permutation
 multivariate analysis of variance (PERMANOVA), all based on Bray-Curtis
 dissimilarity. One-way analysis of variance (ANOVA), accompanied by
 Duncan's post-hoc test, was used to compare the significant
 differences in soil single function, soil multifunctionality, plant
 attributes, AMF colonization, and α-diversity under different warming
 treatments. Co-occurrence networks for AMF communities
 under different treatments were constructed using the “WGCNA” package
 based on the Spearman correlation matrix. Only robust correlations with a
 Spearman correlation coefficient | ρ | &gt; 0.60 and a
 <em>p</em>-value &lt; 0.05 after correction by the
 Benjamini-Hochberg procedure were incorporated into the network analysis.
 The topological properties of the network were determined using the
 "igraph" package in R, which included the number of nodes and
 edges, average clustering coefficient, average path length, network
 diameter, graph density, network diameter, average degree, and modularity.
 Finally, the co-occurrence networks of AMF were visualized using Gephi
 (version 0.9.2) and Cytoscape (version 3.9.1) software.
 Using the "Picante" package in R (Kembel et al., 2010),
 the β-nearest taxon index (βNTI) and Bray–Curtis–based Raup–Crick index
 (RC<sub>bray</sub>) were computed to determine the relative
 importance of deterministic and stochastic processes in AMF community
 assembly. Briefly, ecological processes were divided into deterministic
 (homogeneous selection and variable selection) and stochastic (dispersal
 limitation, homogeneous dispersal, and drift) processes based on βΝΤΙ and
 RC<sub>bray</sub> values (Guo et al., 2022). When βNTI
 &lt; -2 and βNTI &gt; +2 represent community assembly dominated by
 homogeneous and variable selection, respectively (Stegen et al., 2015). In
 addition, the relative influences of homogenizing dispersal and dispersal
 limitation processes were quantified by |βNTI| &lt; 2,
 RC<sub>bray </sub>&lt; -0.95, and RC<sub>bray
 </sub>&gt; 0.95, respectively, and the influences of
 “non-dominated” processes were assessed by |βNTI| &lt; 2 and
 |RC<sub>bray</sub>| &lt; 0.95 (Stegen et al.,
 2015). We performed ordinary least squares linear
 regressions between AMF community α-diversity, β-diversity, network
 complexity, assembly processes, mycorrhizal colonization rate, and soil
 multifunctionality, and conducted Pearson correlation analysis between the
 α-diversity, β-diversity, network complexity, and single function. In
 addition, the linear regression correlations between AMF assembly
 processes, diversity, and network complexity were investigated. We then
 used partial least squares structural equation modeling (PLS-SEM) to
 assess the hypothesized direct and indirect relationships between warming
 degree, AMF diversity, network complexity, community assembly processes,
 mycorrhizal colonization rate, and soil multifunctionality. The first step
 in PLS-SEM analysis is to construct an a priori model based on known
 causal relationships between potential drivers of soil multifunctionality
 (Hu et al., 2021; Zhu et al., 2024). Before embarking on modeling, it is
 vital to conduct a bivariate correlation analysis between all the
 variables to ensure that the relationship between them is sufficiently
 linear. Finally, the model was evaluated based on fit indices such as
 Cronbach's alpha, composite reliability, and average variance
 extracted (AVE), with values greater than 0.7, 0.6, and 0.5, respectively
 (Zhu et al., 2024), indicating that the model performs well regarding
 internal consistency, reliability, and convergent validity. In addition,
 the Goodness of Fit (GoF) index was used to assess the model's
 overall predictive ability.
description: # Warming modulates soil multifunctionality through assembly processes and
 co-occurrence patterns of arbuscular mycorrhizal fungal communities in
 drylands Dataset DOI: [10.5061/dryad.7m0cfxq7n](10.5061/dryad.7m0cfxq7n)
 ## Description of the data and file structure This dataset contains data
 collected during a three-year field warming experiment (+0.5-1.6°C)
 conducted in the Tengger Desert, China, using Open-Top Chambers (OTCs).
 The study investigated the effects of warming on arbuscular mycorrhizal
 fungal (AMF) communities in the rhizosphere soil of two dominant shrub
 species (*Artemisia ordosica* and *Caragana korshinskii*), and their links
 to soil multifunctionality. Data includes: (1) AMF community composition
 data derived from high-throughput sequencing (e.g., OTU/ASV tables,
 taxonomy), (2) measurements of 15 individual soil functional indicators
 related to nutrient availability, biogeochemical cycles, microbial
 activity, and productivity, (3) calculated soil multifunctionality
 indices, (4) mycorrhizal colonization rates, (5) soil moisture content,
 and (6) data used for AMF co-occurrence network analysis and community
 assembly process calculations (e.g., null model results). ### Files and
 variables #### File: NT-otu_taxon_otu.full.csv **Description:** OTU Data
 Sheet for Caragana korshinskii（NT） ##### Variables * Variables: Column
 Variable Name Description Data Type Example Values Notes A domain
 Top-level taxonomic rank (always Eukaryota) Categorical d__Eukaryota
 Prefix d_\_ denotes domain B kingdom Taxonomic kingdom (always Fungi)
 Categorical k__Fungi Prefix k_\_ denotes kingdom C phylum Taxonomic phylum
 (always Glomeromycota) Categorical p__Glomeromycota Prefix p_\_ denotes
 phylum D class Taxonomic class within Glomeromycota Categorical
 c__Glomeromycetes, c__Archaeosporomycetes Prefix c_\_ denotes class E
 order Taxonomic order Categorical o__Glomerales, o__Diversisporales Prefix
 o_\_ denotes order F family Taxonomic family Categorical f__Glomeraceae,
 f__Claroideoglomeraceae Prefix f_\_ denotes family G genus Taxonomic genus
 Categorical g__Glomus, g__Claroideoglomus Prefix g_\_ denotes genus H
 species Taxonomic species or unclassified designation Categorical
 s__sp.*VTX00098, s__unclassified_g__Glomus Prefix s*\_ denotes species I
 otu OTU identifier (unique code for each fungal sequence variant)
 Categorical OTU215, OTU5 J–AI N_* Sample abundance columns (51 columns,
 e.g., N_0_5, N_1_5_1) Integer 0, 171, 680, 37253 Raw sequence counts per
 sample. 0 = no detection. AJ Total Sum of sequences across all samples for
 each OTU Integer 532, 2, 19120 Total abundance per OTU AK Prevalence
 Proportion of samples where OTU was detected Decimal (0–1) 0.2667, 0.0333,
 0.4667 Calculated as (samples with OTU) / (total samples) AL Percent
 Relative abundance of OTU across all sequences Decimal (0–1) 0.000476,
 1.79e-6, 0.0171 (Total sequences of OTU) / (all sequences) Key
 Abbreviations: OTU: Operational Taxonomic Unit (molecular proxy for a
 species/taxon). N_*: Sample IDs (e.g., N_0_5 = sample from specific
 treatment/replicate). unclassified: Incomplete taxonomic assignment (e.g.,
 s__unclassified_g__Glomus = species unclassified within genus *Glomus*).
 #### File: YH-otu_taxon_otu.full.csv **Description:** OTU Data Sheet for
 *Artemisia ordosica*（YH） ##### Variables * Variables: Column Variable Name
 Description Data Type Example Values Notes A domain Top-level taxonomic
 rank (always Eukaryota) Categorical d__Eukaryota Prefix d_\_ denotes
 domain B kingdom Taxonomic kingdom (always Fungi) Categorical k__Fungi
 Prefix k_\_ denotes kingdom C phylum Taxonomic phylum (always
 Glomeromycota) Categorical p__Glomeromycota Prefix p_\_ denotes phylum D
 class Taxonomic class within Glomeromycota Categorical c__Glomeromycetes,
 c__Archaeosporomycetes Prefix c_\_ denotes class E order Taxonomic order
 Categorical o__Glomerales, o__Diversisporales Prefix o_\_ denotes order F
 family Taxonomic family Categorical f__Glomeraceae,
 f__Claroideoglomeraceae Prefix f_\_ denotes family G genus Taxonomic genus
 Categorical g__Glomus, g__Claroideoglomus Prefix g_\_ denotes genus H
 species Taxonomic species or unclassified designation Categorical
 s__sp.*VTX00098, s__unclassified_g__Glomus Prefix s*\_ denotes species I
 otu OTU identifier (unique code for each fungal sequence variant)
 Categorical OTU215, OTU5 J–AI N_* Sample abundance columns (51 columns,
 e.g., N_0_5, N_1_5_1) Integer 0, 171, 680, 37253 Raw sequence counts per
 sample. 0 = no detection. AJ Total Sum of sequences across all samples for
 each OTU Integer 532, 2, 19120 Total abundance per OTU AK Prevalence
 Proportion of samples where OTU was detected Decimal (0–1) 0.2667, 0.0333,
 0.4667 Calculated as (samples with OTU) / (total samples) AL Percent
 Relative abundance of OTU across all sequences Decimal (0–1) 0.000476,
 1.79e-6, 0.0171 (Total sequences of OTU) / (all sequences) Key
 Abbreviations: OTU: Operational Taxonomic Unit (molecular proxy for a
 species/taxon). N_*: Sample IDs (e.g., N_0_5 = sample from specific
 treatment/replicate). unclassified: Incomplete taxonomic assignment (e.g.,
 s__unclassified_g__Glomus = species unclassified within genus *Glomus*).
 #### File: #### **DataYH-FE-2025-00579.csv** **DataNT-FE-2025-00579.csv**
 **Data_AMF_colonization-FE-2025-00579.csv** **DataYH1-FE-2025-00579.csv**
 **DataNT1-FE-2025-00579.csv**
 **DataYH_soil_and_plant__properties-FE-2025-00579.csv**
 **DataNT_soil_and_plant__properties-FE-2025-00579.csv**
 **Description:** Raw measurements of 15 soil function indicators, plant
 properties, and AMF colonization data for two plant species (*Artemisia
 ordosica*, code: YH; *Caragana korshinskii*, code: NT) across multiple
 treatments. ##### Variables * **Sheet: YH** * **Description**: Soil
 function data for *Artemisia ordosica* (YH) with 6 replicates per
 treatment (CK, T1, T2, T3, T4). * **Variables**: | Column | Variable |
 Unit | Description | | :----- | :-------- | :--------- |
 :--------------------------------------- | | A | Sample ID | - | e.g.,
 `CK-1` (Control, Replicate 1) | | B | S-β-GC | μg·g⁻¹·h⁻¹ | Soil
 β-glucosidase activity | | C | S-C1 | μg·g⁻¹·h⁻¹ | Soil cellobiohydrolase
 activity | | D | S-NAG | μg·g⁻¹·h⁻¹ | Soil N-acetyl-β-glucosaminidase
 activity | | E | S-LAP | ng·g⁻¹·h⁻¹ | Soil leucine aminopeptidase activity
 | | F | S-ACP | μg·g⁻¹·h⁻¹ | Soil acid phosphatase activity | | G | TOC |
 g/kg | Total organic carbon | | H | TN | g/kg | Total nitrogen | | I |
 NO₃⁻-N | mg/kg | Nitrate nitrogen | | J | NH₄⁺-N | mg/kg | Ammonium
 nitrogen | | K | AP | mg/kg | Available phosphorus | | L | AN | mg/kg |
 Available nitrogen | | M | DOC | mg/kg | Dissolved organic carbon | | N |
 DON | mg/kg | Dissolved organic nitrogen | | O | MBC | mg/kg | Microbial
 biomass carbon | | P | MBN | mg/kg | Microbial biomass nitrogen | *
 **Sheet: NT** * **Description**: Soil function data for *Caragana
 korshinskii* (NT) with 6 replicates per treatment (same structure as YH).
 * **Variables**: Identical to sheet "YH". * **Sheet: YH-1** *
 **Description**: Calculated z-scores and soil multifunctionality indices
 for YH. * **Variables**: | Column | Variable | Unit | Description | |
 :--------- | :--------------------------- | :------------- |
 :------------------------------------------------ | | B-P, R-AD | Raw data
 | Original units | Same as sheet "YH" (columns B-P) | | C, E,
 G... | z-score | - | Standardized value: `(sample - mean) / std. dev.` | |
 AF | averaged z-score | - | Sum of all 15 z-scores | | AG | soil
 multifunctional indices | - | Soil multifunctionality index: `AF / 15` | *
 **Sheet: NT-1** * **Description**: Calculated z-scores and soil
 multifunctionality indices for NT. * **Variables**: Identical to sheet
 "YH-1". * **Sheet: soil multifunctionality** * **Description**:
 Summary statistics (mean ± SE) and ANOVA results for soil/plant
 properties. * **Variables**: | Column | Variable | Unit | Description | |
 :----- | :------------------------ | :------------- |
 :------------------------------------------------------- | | B | Species |
 - | `Artemisia ordosica` (YH) or `Caragana korshinskii` (NT) | | C |
 Treatments | - | Control, T1, T2, T3, T4 | | D-R | Soil enzymes &
 properties | Original units | Same as sheets "YH/NT" (e.g., βG =
 S-β-GC) | | S | pH | - | Soil pH | | T | SWC | % | Soil water content | |
 U | RB | g | Root biomass | | V | AB | g | Above-ground biomass | | W |
 R/S | - | Root-to-shoot ratio | | X | AMF colonization | % | Arbuscular
 mycorrhizal fungi colonization | | Y | F-statistic (F) | - | ANOVA F-value
 | | Z | P-value (P) | - | ANOVA significance (e.g., `&lt;0.001`) | *
 **Sheet: Soil and plant properties** * **Description**: Sample-level plant
 biomass and soil properties. * **Variables**: | Column | Variable | Unit |
 Description | | :----- | :-------- | :--- |
 :------------------------------------- | | A | Sample ID | - | e.g.,
 `CK-1` | | C | AB\_YH | g | Above-ground biomass (YH) | | D | RB\_YH | g |
 Root biomass (YH) | | E | R/S\_YH | - | Root-to-shoot ratio (YH) | | H |
 AB\_NT | g | Above-ground biomass (NT) | | I | RB\_NT | g | Root biomass
 (NT) | | J | R/S\_NT | - | Root-to-shoot ratio (NT) | | O | pH\_YH | - |
 Soil pH (YH) | | Q | pH\_NT | - | Soil pH (NT) | | T | SWC\_NT | - | Soil
 water content (NT, fraction: 0–1) | | V | SWC\_YH | - | Soil water content
 (YH, fraction: 0–1) | * **Sheet: AMF colonization** * **Description**:
 Arbuscular mycorrhizal fungi colonization rates. * **Variables**: | Column
 | Variable | Unit | Description | | :----- | :-------- | :--- |
 :-------------------- | | A | Sample ID | - | e.g., `CK-1` | | B | YH | %
 | AMF colonization (YH) | | C | NT | % | AMF colonization (NT) | --- ####
 **Key Abbreviations** * **S-β-GC**: β-glucosidase * **S-C1**:
 Cellobiohydrolase * **S-NAG**: N-acetyl-β-glucosaminidase * **S-LAP**:
 Leucine aminopeptidase * **S-ACP**: Acid phosphatase * **TOC**: Total
 organic carbon * **TN**: Total nitrogen * **NO₃⁻-N**: Nitrate nitrogen *
 **NH₄⁺-N**: Ammonium nitrogen * **AP/AN**: Available phosphorus/nitrogen *
 **DOC/DON**: Dissolved organic carbon/nitrogen * **MBC/MBN**: Microbial
 biomass carbon/nitrogen * **AMF**: Arbuscular mycorrhizal fungi * **SWC**:
 Soil water content * **RB/AB**: Root biomass/Above-ground biomass *
 **R/S**: Root-to-shoot ratio * **SMF**: Soil multifunctionality index #
 File: # **network.R** Description: R code. ## Access information Other
 publicly accessible locations of the data: * NCBI Data was derived from
 the following sources: * All raw sequencing data have been deposited in
 the NCBI Sequence Read Archive (SRA) under study accession number
 PRJNA988244.