MATERIALS AND METHODS
2.1
Site selection and soil collection
Field sampling was conducted in
Guizhou Province and the Guangxi Autonomous Region of southwestern
China. The spatial range was 104°82′–108°37′E, 22°42′–27°53′N (Fig.
1). The soil type was calcareous lithosols (limestone soil) according to
the Soil Taxonomy of China and the FAO/UNESCO system of soil
classification. Samples were collected across six counties from typical
karst regions with plantation forest (PF), shrubland (SH), and reference
cropland (CR). Plantation forest and shrubland comprised artificial and
natural vegetation restoration, respectively, and cultivation was
abandoned there in 2002–2003. Shuicheng, Jinsha, and Duyun counties are
in Guizhou Province and Huanjiang, Mashan, and Longzhou counties are in
Guangxi Autonomous Region. The mean annual temperatures (MAT) within the
Guizhou and Guangxi sampling sites were 14.6 °C and 20.9 °C,
respectively. The mean annual precipitation (MAP) levels in Guizhou and
Guangxi were 1,100 mm and 1,260 mm, respectively. Further details of the
climate conditions, geographic locations, and
dominant plant species are listed
in Table S1.
Soil sample collection was performed between August and October 2018.
There were 54 soil samples (three land use types × six sampling counties
× three repetitions). In each 30 m × 30 m sampling field, twenty
38-mm-diameter soil cores were collected at 15 cm soil depth and pooled.
Each soil sample was passed through a 2-mm sieve to remove rocks, roots,
and organic debris. Physicochemical properties were determined for
air-dried soil. Genomic DNA was extracted from soil stored at -80 °C.
2 .2 Soil
physicochemical property measurements
The
pH of 10 g soil suspended in 25 mL water was measured with a glass
electrode pH meter. Exchangeable calcium (Ca) and magnesium (Mg) were
displaced via compulsive exchange in 1 M ammonium acetate (pH 7). Soil
organic carbon (SOC) was measured by dichromate redox colorimetry.
Particulate organic C (POC) and mineral-associated organic C (MOC) were
isolated by size fractionation (> 53 μm and < 53
μm) and analyzed in the same way as the bulk SOC. Total N (TN) was
determined using an elemental (CN) analyzer (Vario MAX CN; Elementar,
Hanau, Germany). Ammonium N (NH4+) and
nitrate N (NO3−) were extracted from
10 g fresh soil by extraction with 0.5 M KCl and measured with an auto
flow analyzer. Total phosphorus (TP) and available phosphorus (AP) were
determined by digestion in a solution containing
H2SO4 and HClO4 and by
molybdenum blue colorimetry, respectively (Hu et al., 2021; Xiao et al.,
2019, 2020).
2.3 Soil DNAextraction and amplicon
sequencing
The
total DNA was extracted from 0.5 g soil samples using a FastDNA Spin Kit
(MP Biomedicals, Santa Ana, CA, USA) according to the manufacturer’s
instructions. DNA quality and quantity were determined with a NanoDrop
2000 UV-vis spectrophotometer (Thermo Fisher Scientific, Wilmington, DE,
USA). Bacterial 16S rDNA gene was amplified using the primer pairs 338F
(5’-ACTCCTACGGGAGGCAGCAG-3’) and 806R (5’-GGACTACHVGGGTWTCTAAT-3’).
Fungal ITS and protist 18S rRNA genes were amplified with the primer
sets ITS1F (5’-CTTGGTCATTTAGAGGAAGTAA-3’)/ITS2R
(5’-GCTGCGTTCTTCATCGATGC-3’) and TAReuk454FWD1
(5’-CCAGCA(G/C)C(C/T)GCGGTAATTCC-3’)/TAReukREV3 (5’-ACTTTCGTTCTTGAT
(C/T)(A/G)A-3’), respectively. Samples were distinguished by attaching
barcode sequences to each primer.
PCR was performed in a 20 µL mixture containing 4 µL of 5× FastPfu
Buffer, 2 µL of 2.5 mM dNTPs, 0.8 µL of 5 μM forward and reverse
primers, 0.4 µL TaKaRa rTaq DNA polymerase (TaKaRa Bio Inc., Shiga,
Japan), 0.2 µL BSA, 10 ng template DNA, and 20 µL ddH2O.
The thermal cycles for the 16S
rDNA and ITS genes were, 95 °C for 3 min, 27 cycles and 35 cycles for
16S rDNA and ITS genes, respectively, 95 °C for 30 s, 55 °C for 30 s, 72
°C for 45 s, and a final extension at 72 °C for 10 min. The PCR program
used for protist 18S rRNA gene was, initial denaturation at 95 °C for 5
min, 10 cycles of denaturation at 94 °C for 30 s, annealing at 57 °C for
45 s, extension at 72 °C for 60 s, 25 cycles of denaturation at 94 °C
for 30 s, annealing at 45 °C, 47 °C, 48 °C, and 49 °C for 45 s each,
extension at 72 °C for 60 s, and a final extension at 72 °C for 2 min.
PCR was performed at Shanghai Majorbio Bio-Pharm Technology Co. Ltd.
(Shanghai, China) on an Illumina MiSeq PE300 platform (Illumina, San
Diego, CA, USA).
2.4 Sequence analysis
Bioinformatic analyses of all
three amplicons were performed with QIIME 2 2021.2 (Bolyen et al.,
2019). The raw reads were imported into the QIIME 2 environment and
denoised with the DADA2 algorithm via the q2-data2 plugin (Callahan et
al., 2016). Operational taxonomic units (OTUs) were selected at 97%
sequence similarity in the q2-vsearch plugin (Rognes et al., 2016). OTU
taxonomies for bacteria, fungi, and protist were assigned using the
q2-feature-classifier (Bokulich et al., 2018) vs. the SILVA_v132 (Quast
et al., 2012), UNITE_v8.3 (Nilsson et al., 2019), and PR2 (V4.13.0)
(Guillou et al., 2012) reference sequences, respectively.
Alpha-diversity metrics (richness and Shannon index) were calculated in
the “RAM” package of R (R Core Team, Vienna, Austria).
2.5 Statistical analysis
Data were processed in R v. 3.6.3 (R Core Team). Significant differences
in microbial diversity, microbial groups, and environmental factors
among vegetation types were analyzed by one-way ANOVA and Duncan’s post
hoc test (p < 0.05). Three-way nested PERMANOVA was
used to determine the main and interactive effects of MAT, MAP, and
vegetation type on microbial diversity and community composition (Xiao
et al., 2019). Non-metric
multidimensional scaling (NMDS) was performed to reveal dissimilarities
in microbial, fungal, and protist community composition among
treatments. Pearson’s correlation coefficients and the Mantel test were
applied to determine the relationships among environmental factors
(geographical factors, climatic conditions, and soil properties),
microbial diversity, and community composition.
To determine the structure and identify interactions among microbial
groups, bacterial, fungal, and protist community networks were
constructed across all samples and vegetation types under low (GZ) and
high (GX) climate conditions. Spearman correlations were calculated to
clarify the associations among the microbial taxa detected in
> 1/3 of all samples per network. Significant Spearman
associations (ρ > 0.75; 1,000-fold permutations; p< 0.05) and the corresponding OTUs served as network nodes and
edges, respectively. Positive and negative associations suggested
cooperative and competitive relationships, respectively. Network
properties such as node, edge, density, and degree were calculated in
each network using the “igraph” package in R (R Core Team) (Csardi and
Nepusz, 2006).
Structural equation modeling (SEM)
based on SPSS Amos 21 (IBM Corp., Armonk, NY, USA) was conducted to
evaluate the causal relationships among MAT, vegetation type, soil
property, and microbial community composition (Hu et al., 2021; Yang et
al., 2017). To simplify the SEM, a principal component analysis (PCA) of
the first component (PC1) was applied to represent soil properties and
microbial community compositions based on bacterial, fungal, and protist
NMDS axes 1 and 2. Soil properties included pH, soil exchangeable
Ca2+ (Ca), soil exchangeable Mg2+(Mg), soil organic carbon (SOC), particulate organic carbon (POC),
mineral-associated organic carbon (MOC), total nitrogen (TN), ratio of
soil carbon to soil nitrogen (C/N), ammonium nitrogen
(NH4+), nitrate nitrogen
(NO3−), total phosphorus (TP), and
available phosphorus (AP).
3
RESULTS
3.1Bacterial, fungal, and
protist diversity and community composition
Only
bacterial diversity was affected by temperature (Table 1). Bacterial
richness and the Shannon indices for plantation forest and shrubland
were lower in Guangxi (high temperature) than Guizhou (low temperature)
(Fig. 2a). Vegetation type significantly affected bacterial, fungal, and
protist diversity (Table 1). Bacterial richness and Shannon index were
lower under shrubland than cropland (except for richness in Guizhou)
(Fig. 2a). Fungal richness was higher in shrubland than cropland under
both climate conditions (Fig. 2b). At Guangxi, protist Shannon index was
lower under plantation forest and shrubland than cropland. Cropland
protist richness was higher than that of shrubland at Guangxi (Fig. 2c).
There were no differences between plantation forest and shrubland in
terms of bacterial, fungal, and or protistan diversity (Fig. 2).
Across all sites, phylum-level bacterial community composition was
dominated by Actinobacteria (29.3%), Proteobacteria (25.6%), and
Acidobacteria (16.7%) (Fig. 3a). The most
abundant fungal phyla were
Ascomycota (45.0%) and Basidiomycota (27.9%) (Fig. 3b). The main
protist phyla were Cercozoa (33.9%) and Apicomplexa (22.1%) (Fig. 3c).
NMDS ordination revealed that bacterial, fungal, and protist community
composition differed among vegetation types between Guizhou and Guangxi.
Bacterial, fungal, and protist community composition differed between
cropland and vegetation recovery but were similar between plantation
forest and shrubland (Figs. 3d-f). Temperature and vegetation type
significantly affected bacterial, fungal, and protist community
composition (Table 1).
The relative abundance of Proteobacteria (bacteria) under all vegetation
types was lower at Guangxi than Guizhou (Table S2, Fig. S2). The
relative abundance of Mucoromycota (fungi) was higher under plantation
forest and shrubland at Guizhou than it was at
Guangxi
(Table
S3, Fig. S3). The relative abundance of Cercozoa (protist) was higher
under plantation forest and shrubland at Guizhou than it was at
Guangxi (Table S4, Fig. S4).
Relative bacterial phyla Cyanobacteria, Gemmatimonadetes, and
Nitrospirae abundances were higher under cropland than plantation forest
or shrubland. Relative Proteobacteria and Verrucomicrobia abundances
were higher under shrubland than cropland (Table S2, Fig. S2). Relative
fungal phyla Ascomycota and Mucoromycota abundances were higher under
cropland than plantation forest and shrubland. Nevertheless, the
opposite was true for Basidiomycota (Table S3, Fig. S3). Relative
protist phyla Ciliophora, Lobosa, and Ochrophyta abundances were higher
under cropland than plantation forest or shrubland. However, the
opposite was true for Apicomplexa (Table S4, Fig. S4).
3.2
Co-occurrence networks amongbacterial, fungal, and
protist
taxa
Microbiome network complexity was greater at Guizhou (GZ) than Guangxi
(GX). There were more nodes and edges under GZ than GX. The most densely
connected network occurred under cropland (Table 2, Fig. 4). There were
fewer positive and negative correlation links among bacterial, fungal,
and protist taxa under vegetation restoration than cropland. Shrubland
had the fewest positive and negative correlation links among microbial
taxa (Fig. 4). Node connectivity was greater under cropland than
plantation forest or shrubland. Shrubland had the fewest node
connections. Under all vegetation and both climate types, bacteria had
the most nodes, followed by fungi and protists (Table 2, Fig. S5).
Acidobacteria (bacteria), Ascomycota (fungi), and Cercozoa (protist)
formed remarkably centralized networks and showed the most links with
other microbial groups (Table S5, Fig. 4).
Bacterial phyla Acidobacteria and
Proteobacteria were the most closely connected with Ascomycota and
Cercozoa. The strongest and most numerous links between fungal and
protist taxa were determined for Ascomycota and Cercozoa, followed by
Ascomycota and Lobosa (Table S5). There were more interactions between
bacterial and fungal taxa (edge links) than there were between bacterial
and protist taxa or between fungal and protist taxa. The proportions of
edges linking bacteria to protists were lower under plantation forest
and shrubland than cropland. Conversely, the proportions of edges
linking fungi to protists were higher under plantation forest and
shrubland than cropland (Table 2).
3.3
Relationships among environmental factors and bacterial, fungal, and
protist communities
The
various bacterial, fungal, and protist taxa had different effects on
soil carbon (SOC, POC, and MOC), nitrogen (TN,
NH4+, and
NO3−), and phosphorus (TP and AP).
More bacterial, fungal, and protist OTUs correlated with C and N than
they did with P (Fig. 5). For bacterial taxa, soil SOC, POC, and MOC
were most strongly correlated with Proteobacteria. Acidobacteria and
Proteobacteria had the most OTUs correlated with soil TN,
NH4+,
NO3−, TP, and AP (Fig. 5a). Soil
carbon (SOC, POC, and MOC), nitrogen (TN,
NH4+, and
NO3−), and phosphorus (TP and AP) were
strongly correlated with fungal phyla Ascomycota (Fig. 5b). Protist
phyla Cercozoa were strongly correlated with soil C, N, and P (Fig. 5c).
Pearson’s correlation coefficients showed that bacterial richness and
Shannon index were negatively correlated with MAT, MAP, and MOC but
positively correlated with latitude and altitude (Fig. 6a). Fungal
diversity was not significantly correlated with any geographical factor,
climate condition, or soil property (Fig. 6b). Protist diversity was
positively correlated with AP but negatively correlated with SOC, POC,
TN, and NH4+ (Fig. 6c). Bacterial,
fungal, and protist community compositions were strongly affected by
environmental factors (geography, climate, and soil physicochemistry)
(Fig. 6).
Structural equation modeling (SEM) evaluated the linkages among
environmental factors and microbial community composition, and ranked
the relative importance of each driver on the microbial communities.
Temperature and vegetation type directly influenced soil property and
indirectly affected bacterial, fungal, and protist community
composition. MAT had greater standard total effect on bacterial and
fungal community composition than it did on protist community
composition. Vegetation type had the highest standard total effect on
protist community composition followed by fungal community composition
and bacterial community composition. Bacterial and fungal community
composition was affected by MAT and vegetation type. Moreover,
vegetation type and soil property had greater influences than MAT on
protist community composition (Fig. 7).