2.5 Statistical analysis
All data analyses were performed in R v4.0.4 (https://www.r-project.org/ ). One-way Analysis of Variation and the Least Significant Difference multiple comparison test were used to investigate differences in rhizosphere soil variables.
Microbial α-diversity (OTU richness, Shannon index, and Simpson index) was calculated using the “vegan ” and “picante ” packages in R v4.0.4. Variation analysis was performed using “EasyStat ” and “dplyr ” packages, followed by visualization using the “ggplot2 ” package, in R v4.0.4.
Microbial β-diversity was estimated based on the Bray-Curtis distances between samples. Bray-Curtis dissimilarity index was used to assess differences in β-diversity, and Wilcoxon rank sum tests were performed to investigate significant differences in microbial β-diversity. Unconstrained Principal Co-ordinate Analysis (PCoA) and Permutational Multivariate Analysis of Variance (PERMANOVA) were carried out using the “adonis ” function in the “vegan ” package in R v.4.0.4, and the results were visualized using the “ggplot2” in R v4.0.4. Subsequently, the “phyloseq ” package was used to extract normalized abundance of the top 10 phyla and phylogenetic trees were constructed using the “hclust” function in the “ggtree” package in R v4.0.4 based on the Bray-Curtis distance. Stacked bar plots and trees were combined using the “faceted ” function in the “ggplot2” package. In addition, β-Diversity of the core microbiota and other indices were quantified using the first two axes (PCoA1 and PCoA2) of a PCoA plot.
Correlation analyses between microbial community composition (based on Bray-Curtis distance) and soil environmental properties were conducted using the “mental test ” function in the “ggcor ” package in R v4.0.4. In the Mantel test, microbial community structure and soil environmental properties were considered distance matrices; and the correlation between the two matrices was tested.
Forest age-sensitive OTUs (fsOTU) were identified according to the method of Harman (17). First, we used the “Indicspecies ” package in R v4.0.4 to identify the indicator species. Subsequently, we tested for differences in OTU abundances of bacteria and fungi between two or more forest ages using the likelihood ratio test with the “edgeR ” package in R v4.0.4. The OTUs with significantly different relative abundances between two or more forest ages (corrected false discovery rate p < 0.05) were considered responsive to forest age. Afterward, OTUs that were confirmed by both indicator species analysis and likelihood ratio test were defined as forest age-sensitive OTUs (fsOTU).
Microbial co-occurrence networks were constructed for bacteria and fungi separately and visualized using the “igraph ” package in R v4.0.4 based on the Fruchterman-Reingold layout. The trimmed mean of M-values was used to normalize relative abundance to counts per million. The Spearman correlation between every two OTUs was estimated. Robust correlations with Spearman correlation coefficients >0.7 and p values <0.01 were identified for use in constructing networks, in which each node represents one OTU (16S_OTU for bacteria and ITS_OTU for fungi) and each edge represents a strong and significant correlation between two nodes. The descriptive and topological properties of the networks, including node number, edge number, and co-occurrence degree were calculated.
Core microbiota were identified in rhizosphere bacterial and fungal communities in R. pseudoacacia plantations with different ages. In all the samples, the ubiquitous (present in ≥ 80% of samples) and highly abundant (top 10 in terms of relative abundance) OTUs were selected (5, 10, 25). For the core microbiota, co-occurrence networks were constructed for bacteria and fungi among all samples after removing OTUs with low abundance (<0.001%) from the dataset. Spearman rank correlations were conducted between OTUs, and significant positive correlations (r > 0.07, p < 0.001) were visualized. The node-level topological feature (degree) was calculated for each node.
The associations between core microbiota and rhizosphere soil multinutrient cycling were evaluated using random forest analysis (5). Microbial α-diversity and β-diversity served as predictors of the rhizosphere soil multinutrient cycling index for estimating the contributions of the microbiota to ecosystem functions (25). Subsequently, multiple regression modelling with variance decomposition analysis was used to estimate the importance of the major taxonomic categories (class) of the core microbiota in explaining rhizosphere soil nutrient properties, based on the backward filtering approach. The “rfPermute ” package in R v4.0.4 was used to estimate the significance of the influence of the environmental variables on core microbiota abundance.
Results