3.1. Soil properties
Most of the soil physiochemical properties showed significant differences between the rhizosphere and bulk soil along the chronosequence (Table S2). The contents of SOC, SMC and TN increased gradually along the chronosequence and peaked in the 60-year-old plantation. However, soil pH showed an opposite trend. The SOC, AP, and C:N ratios were higher in the rhizosphere soil than in the bulk soil. Additionally, compared to the bulk soil, the pH in the rhizosphere soil was lower and more stable. Two-way ANOVAs indicated that the interaction between the chronosequence and soil compartment was significant for soil pH, TN, and NO3−-N.
3.2. Microbial community composition and diversity
In total, 62 771 bacterial ASVs and 6 906 fungal ASVs were obtained. Across all samples, eight bacterial phyla/classes, namely, Acidobacteria (21.66%), Actinobacteria (19.21%), Gammaproteobacteria (11.81%), Chloroflexi (10.60%), Alphaproteobacteria (9.32%), Gemmatimonadetes (7.96%), Deltaproteobacteria (5.37%) and Rokubacteria (4.62%), dominated the bacterial communities (Fig. S1). The most abundant fungal class was Agaricomycetes (65.11%), followed by Sordariomycetes (3.60%), Pezizomycetes (3.44%), Archaeorhizomycetes (1.84%), Mortierellomycetes (1.52%), and Leotiomycetes (1.29%) (Fig. S1). The relative abundance of Alphaproteobacteria and Gammaproteobacteria in the rhizosphere soil increased gradually along the chronosequence and was significantly higher than that in the bulk soil. The abundance of all of the detected fungi except Agaricomycetes and Sordariomycetes in the bulk soil decreased along the chronosequence. Pezizomycetes, Mortierellomycetes, and Eurotiomycetes were more abundant in the rhizosphere soil than in bulk soil from 30 to 60 years.
For bacteria, the alpha diversity (Shannon and Chao1) showed no significant differences along the chronosequence but significantly (P < 0.05) increased from the bulk soil to the rhizosphere. The chronosequence and soil compartment had no effects on the Shannon diversity and richness of the fungal communities (Table S3).
PCoA showed that the bacterial and fungal compositions were clearly separated among the three forest ages (Fig. 1). The first two axes explained 41.21% and 41.62% of the total variation in the soil bacterial and fungal communities, respectively. Two-way permutational multivariate analysis of variance (PERMANOVA) revealed that the chronosequence had a stronger influence than the soil compartment on the bacterial and fungal community structures (Fig. 1). The Mantel tests indicated that the bacterial community structure had significant correlations with SMC, SOC and AP. SMC was the key factor correlated with the fungal community structure (Table S4). The dominant bacterial taxa and diversity appeared to have a stronger interrelationship than fungi with the soil properties (Table S5).
Significant TDRs were found for the soil bacterial community in the different soil compartments (Fig. 2). The slopes of the bacterial communities in the rhizosphere soil were significantly steeper than those in the bulk soil, revealing rapid community turnover in the rhizosphere soil. However, a significant TDR was observed for fungal communities only in the rhizosphere soil, which suggested that the turnover rates of fungal communities changed dramatically along the chronosequence.
3.3 Microbial community networks
Multiple network topological properties revealed that the co-occurrence patterns of the bacterial and fungal networks changed constantly along the chronosequence (Fig. 3 and Table S6). All networks showed strong power-law distributions of degrees, indicating typical scale-free degree characteristics (Fig. S2). For bacteria, the network diameter and average path length decreased along the chronosequence, while the number of edges, graph density, average degree, and clustering coefficient tended to increase, suggesting a more connected network structure with stand age (Table S6). In addition, the positive correlations in the bacterial network in 60-year-old plantations increased by 9.98% and 12.55% compared with those in 30-year-old and 15-year-old plantations, respectively. This result indicated that the conflicting interactions (the competition among the bacterial species) decreased gradually along the chronosequence (Table S6). For fungi, the clustering coefficient, graph density, and average degree increased from 15 years to 30 years and later decreased at 60 years. The trends of the positive correlations were the opposite of those found in bacterial networks, indicating that the increase in stand age might lead to more conflicting interactions (Table S6). In addition, the bacterial networks were similar between the rhizosphere and bulk soil. However, fungal networks were different between the different soil compartments. The rhizosphere soil had higher numbers of nodes (86) and edges (188) and a higher average degree (4.37) than bulk soil, suggesting more complex relationships among the fungal species.
There were a total of 321 connectors, 13 module hubs (network hub ASVs), and 2 network hubs identified in the bacterial networks (Fig. S3 and Table S7). However, the numbers and their taxonomic affiliations differed among the three stand ages. The greatest number of connectors was observed in the 15-year-old plantation, which mainly belonged to Acidobacteria, Actinobacteria, Chloroflexi, and Gammaproteobacteria. There were three module hubs belonging to Alphaproteobacteria (1 ASV) and Chloroflexi (2 ASVs). In the 30-year-old plantation, 83 connectors and 5 module hubs were detected among various taxonomic groups. The connectors were mainly related to Acidobacteria and Gammaproteobacteria. Module hubs were related to Acidobacteria (1 ASV), Actinobacteria (2 ASVs) and Alphaproteobacteria (2 ASVs). In the 60-year-old plantation, the 103 connectors were mainly assigned to Acidobacteria, Gemmatimonadetes, and Rokubacteria. The 5 module hubs were assigned to Actinobacteria (1 ASV), Chloroflexi (1 ASV), Deltaproteobacteria (1 ASV) and Gammaproteobacteria (2 ASVs). The only 2 network hubs detected at this site belonged to Actinobacteria and Alphaproteobacteria. In addition, the number of connectors increased from the rhizosphere (63) to the bulk soil (77). No module hubs or network hubs were identified in the rhizosphere soil. In the fungal networks, the number of connectors decreased gradually with stand age (Fig. S3). A higher number of connectors was observed in the rhizosphere (39) than in the bulk soil (32). These connectors originated from a variety of taxonomic groups and were mainly related to Saccharomycetes, Sordariomycetes, Agaricomycetes, and Mortierellomycota (Table S7). A module hub (1 ASV) was identified only at the 30-year-old site, which was related to Tremellomycetes. No network hub was detected in the fungal networks.
The correlation between soil properties and keystone species suggested that bacterial keystone species in the rhizosphere had significantly positive correlations with TP and NO3−-N, and fungal keystone species in the rhizosphere were positively and significantly correlated with SOC, TN, and TP (Table S8). The MRM models revealed that TP and NO3−-N had strong effects on bacterial keystone species in the rhizosphere, and fungal keystone species in the rhizosphere were strongly affected by SMC, TN, TP, and NO3−-N (Table 1).
Table 1
Multiple regression with distance matrices (MRM) to estimate the contribution of soil properties on the rhizospheric keystone species composition.
Soil properties | Bacterial keystone taxa R2 = 0.506 P < 0.05 | Fungal keystone taxa R2 = 0.642 P < 0.01 |
pH | -0.240 | 0.484 |
SMC | ND | -0.064* |
SOC | 0.003 | 0.009 |
TN | -0.013 | 0.220* |
TP | 0.210* | 0.590* |
NO3−-N | 0.032* | 0.150* |
NH4+-N | -0.031 | ND |
AP | ND | ND |
C/N | ND | -0.003 |
Note. ND, not determined. * P < 0.05, ** P < 0.01, *** P < 0.001. |
3.4. Microbial assembly processes
Figure 4 shows that stochastic processes, especially dispersal limitation, dominated both bacterial and fungal community assembly. However, the relative contribution of stochastic and deterministic processes in explaining the variation in microbial community assembly changed with stand age. The relative importance of dispersal limitation decreased by 33.33% from the 15-year-old to the 30-year-old plantation and then increased by 26.67% for the bacterial community assembly, while the reverse trend was observed for the fungal community assembly. Regarding the soil compartment, dispersal limitation explained 55.56% of the variation in the bacterial community assembly in the rhizosphere soil, which was the same as that in the bulk soil. The contribution of dispersal limitation to the bacterial community assembly in the rhizosphere soil was the same as that in the bulk soil. However, dispersal limitation explained more of the variation in the fungal assembly patterns, increasing from 77.78% in the bulk soil to 88.89% in the rhizosphere soil. Homogeneous selection, homogenizing dispersal, drift, and heterogeneous selection slightly affected the assembly of the microbial community. In addition, βNTI and soil properties had no significant correlation, indicating less environmental filtering in controlling microbial community assembly (Table S9).
The migration rates of bacteria were lowest in the 30-year-old plantation, and fungi showed the opposite trend, indicating that the bacterial community was more limited and that the fungal community was less limited by dispersal in the 30-year-old plantation than in the 15- and 60-year-old plantations. In addition, the bacterial migration rates tended to be higher in the rhizosphere than in the bulk soil. The fungal migration rates were contrary to the bacterial migration rates, which suggested that the fungal community was more limited by dispersal in the rhizosphere than in the bulk soil (Fig. 5). Similarly, the habitat niche breadth of bacteria and fungi exhibited the same trends as the migration rates in the different plantations and soil compartments.