Near-natural transformation of Pinus tabuliformis better improve soil nutrients and soil microbial community

Research area

The study area was carried out at the Shahai experimental area (41°29′N, 119°27′E), which is located in the Chaoyang, Liaoning Province. The climate belongs to semi-arid monsoon climate zone, with an annual mean temperature of 7−10 °C, with an absolute maximum temperature of 40.7 °C and minimum temperature of −26.4 °C, and a frost-free period of 145–150 days throughout the year. The annual average precipitation is 400–600 mm, with 60%–70% concentrated in June–August, the annual evaporation is about 1,600–1,800 mm, which is 3 to 4 times the precipitation, and the annual average relative humidity is 50%–60%. This area belongs to the edge of North China flora, the transition zone between North China and Mongolia flora. Woody plants include P. tabulaeformis, Robinia pseudoacacia, Quercus mongolica, Prunus sibirica, Prunus dabian, Zizyphus jujuba var. inermis, Hippophae rhamnoides, Vitex chinensis, Zizyphus jujube var. spinosus and Ostryopsis davidiana. The herbaceous plants are mainly Arundinella hirta, Cleistogenes polyphylla and Aneurolepidium chinensis. The typical soils are classified as Eluvial brown soil according to the FAO Taxonomy. Large areas of P. tabuliformis forests were planted in the 1950s, the government implemented the Natural Forest Conservation Program (NFCP) over the last 20 years (Zhang et al., 2000). Currently, P. tabuliformis (PT) has largely been transformation to Pinus tabuliformis -Armeniaca vulgaris (PTAU), Pinus tabuliformis - Robinia pseudoacacia (PTRP), and Pinus tabuliformis - Vitex negundo L. var. heterophylla (PTVN) (Table S1).

Samples collection

Four forest types were selected in this study in August 2020. In order to reduce the spatial autocorrelation, we established four spatially interspersed stands with 0.06 ha (20 m ×30 m) at the distances >300 m as the repetition in each forest type. After removing the litter layer, we collected the 10–15 soil samples of topsoil (0∼10 cm) using soil auger with an “S” shape and mixed together to form a single composite sample for each triplicate plot, giving a total of 16 samples. All collected soil samples were sealed in plastic bags, and then put into icebox for transport to the laboratory. Soil samples for chemical analysis were initially air-dried, finely ground, and sifted through a 100-mesh (0.15 mm) sieve. The microbial DNA samples were stored at −80 °C for DNA extraction.

The determination of soil characteristics

The soil pH was assayed via a 1:2.5 soil: water (w/v) suspensions using pH meter (Mettler Toledo pH (FE20)). The concentrations of soil total carbon (TC) and total nitrogen (TN) were quantified with an elemental analyzer (Euro Vector EA3000). Soil total phosphorus (TP) and available P (AP) contents were determined following soil digestion with H₂SO₄-HClO₄ and soil extraction with 0.5 mol l⁻¹ NaHCO₃, respectively, and assessed by the spectrophotometer (UV-9000S). The concentration of available N (AN) was determined by the alkali solution diffusion method. The soil C/N ratio was calculated as their mass ratios of respective elements. The determination of available potassium (AK) was performed atomic absorption spectrometry after extraction with 1.0 mol L⁻¹ NH₄OAc.

Soil DNA extraction and Illumina sequencing

Briefly, genomic DNA was extracted from 0.5 g of soil using the FastDNA SPIN Kit (MP Biomedicals, Santa Ana, CA, USA), in accordance with the instructions provided by the manufacturer. The target bacterial V3-V4 regions of 16S rRNA gene were amplified using the forward primer 338F (5′-ACTCCTACGGGAGGCAGCAG-3′) and reverse primer 806R (5′-GGACTACHVGGGTWTCTAAT-3′). And the target fungal DNA gene ITS1 regions were amplified using the forward primer ITS1F (5′-CTTGGTCATTTAGAGGAAGTAA-3′) and reverse primer ITS2 (5′-GCTGCGTTCTTCATCGATGC-3′). The total PCR reaction were carried out with 25 µl mixture, including 2 µl of dNTPs (2.5 mM), 2µl (40–50 ng) of DNA Template, 8.75 µl of ddH₂O, 10 uM (1 µl) of forward primer, 10 uM (1 µl) of reverse primer, 5 U/µl (0.25 µl) of Q5 High-Fidelity DNA Polymerase, 5 µl of Q5 High-Fidelity GC buffer (5 ×), and 5 µl of Q5 reaction buffer (5 ×) (Deng et al., 2018). The following PCR thermal cycling condition consisted of an initial denaturation step at 98 °C for 5 min, followed by 25 cycles of denaturation at 98 °C for 15 s, annealing at 55 °C for 30 s, and elongation at 72 °C for 30 s, with final elongation at 72 °C for 5 min. The PCR amplicons were further purified by Agencourt AMPure Beads (Beckman Coulter, Indianapolis, IN) and quantified using PicoGreen dsDNA Assay Kit (Invitrogen, Carlsbad, CA, USA). The purified amplicons of bacterial and fungal genes were paired-end sequenced on Illumina NovaSeq PE250 platform at Shanghai Personal Biotechnology Co., Ltd, Shanghai, China. The high-throughput sequencing raw data of bacteria and fungi were uploaded in the NCBI database with the SRA accession number of PRJNA720506.

Bioinformation analysis

The obtained raw 16S rRNA and ITS rRNA sequences data were processed (i.e., filtering, clustering, trimming and chimera checking) using the Quantitative Insights Into Microbial Ecology (QIIME) pipeline with default settings. And the operational taxonomic unit (OTUs) were classified at 97% confidence level and clustered by the UPARSE pipeline with a cutoff of occupancy (100%, detected at every location) and relative abundance (>1%). Taxonomic characterization of the representative sequences of bacterial OTUs was performed using the SILVA 16S rRNA database and Unite database.

Statistical analyses

The significant differences in soil physical and chemical properties, the microbial community diversity index, and soil microbial community composition of different forest types were examined with one-way analysis of variance (ANOVA) followed by the least significant difference (LSD) at the 95% confidence level. The shared and unique OTUs among different samples were calculated and visualized using the package of “VennDiagram” in RStudio (Zaura et al., 2009). Non-metric multidimensional scaling (NMDS) analysis of the phylogenetic of OTUs based on Jaccard distance was used to better visualize and investigate the distribution characteristics and the dissimilarity of soil microbial communities among different samples using RStudio with the package of “Vegan” (Lozupone et al., 2007). Permutational multivariate analysis of variance (permanova) was used to analyze the differences between groups using RStudio with the package of “Vegan”. The correlations between soil properties and soil community diversity were studied using Spearman’s rank correlations using the package of “VennDiagram” in RStudio (R Development Core Team, 2009). LEFSE analysis, namely LDA Effect Size analysis, achieved by the rank sum test of non-parametric factor Kruskal-Wallis, can be used to find the groups with significant differences in abundance among different samples (Segata et al., 2011). In order to explore the contribution of edaphic variables on the soil dominant communities at phylum level, redundancy analysis was carried out using CANOCO for Windows, Version 4.5 (Braak & Smilauer, 2002).

Soil characteristics among different vegetation types

Except for AP, soil pH (F = 43.38, P = 0.00), TC (F = 8.15, P = 0.00), TN (F = 5.72, P = 0.01), C/N (F = 18.76, P = 0.00), TP (F = 16.34, P = 0.00), AN (F = 5.68, P = 0.01), and AK (F = 4.47, P = 0.03) existed significant differences among different vegetation types. Compared to PT, PTAU, PTVN and PTRP could increase the soil pH, TC, TN, AN, and AK in different degrees, the most obvious in PTAU with 7.71, 20.10 g/kg, 1.48 g/kg, 108.38 mg/kg, and 102.80 mg/kg. Soil C/N showed the highest in PT with 14.76, followed by PTAU, PTVN and PTRP (Table 1).

Soil microbial community among different vegetation types

A total of 1,097,939 high quality bacterial sequences were obtained across all samples, with average 68,621 sequences per sample, which were clustered into 57,544 OTUs. The number of shared OTUs among PT, PTAU, PTVN and PTRP was 1,160, and the unique OTUs was 12,512, 12,506, 13,055, 10,420 in PT, PTAU, PTVN, PTRP, accounting for 69.40%, 71.79%, 68.66%, and 62.31% of their total OTUs, respectively (Fig. 1A). A total of 2,302,306 high quality fungal sequences were obtained across all samples, with average 115,115 sequences per sample, which were clustered into 3,492 OTUs. The number of shared OTUs among PT, PTAU, PTVN, and PTRP was 200, and the unique OTUs was 394, 1,027, 556, 516 in PT, PTAU, PTVN, and PTRP, accounting for 41.78%, 59.64%, 42.97%, and 43.69% of their total OTUs, respectively (Fig. 1B).

There was considerable variation of soil bacterial Pielou_e index (F = 13.00, P = 0.00), and Simpson index (F = 9.52, P = 0.00) among PT, PTAU, PTRP, and PTVN. Soil bacterial Pielou_e index and Simpson index appeared the highest in PTVN with 0.905 and 0.999, respectively, followed by PTRP and PTAU, and lowest in PT with 0.881 and 0.998. However, soil bacterial Shannon index (F = 2.35, P = 0.12) and Chao1 index (F = 0.55, P = 0.66) appeared no significant differences among different samples (Table 2). Furthermore, soil fungal Pielou_e index (F = 2.31, P = 0.13), Shannon index (F = 3.03, P = 0.07), and Simpson index (F = 1.01, P = 0.42) also existed no obvious difference among PT, PTAU, PTRP and PTVN. While, soil fungal Chao1 index was the highest in PTAU, followed by PTVN and PTRP, and lowest in PT with 409.33 (Table 2).

Soil bacterial Pielou_e index existed obviously positive relation with soil pH (r = 0.52, P = 0.05) and TN (r = 0.67, P = 0.01). Soil C/N showed dramatically negative correlation with Pielou_e index (r = −0.76, P = 0.01), Shannon index (r = −0.68, P = 0.01), and Simpson index (r = −0.92, P = 0.01) (Table 3). Nevertheless, no significant correlations of soil fungal community diversity with soil characteristics were observed (Table 3).

Soil microbial community composition among different vegetation types

Soil bacterial communities ranking at the top ten were Proteobacteria, Actinobacteria, Acidobacteria, Bacteroidetes, Gemmatimonadetes, Chloroflexi, Verrucomicrobia, Patescibacteria, Planctomycetes, and Firmicutes, accounting for 99.05%, 99.01%, 98.54% and 98.83% in PT, PTAU, PTRP, and PTVN, respectively (Fig. 2A). The relative abundance of Proteobacteria in PT was the highest with 41.67%, followed by PTAU, PTTRP, and lowest in PTVN. While, Actinobacteria and Gemmatimonadetes in PTVN appeared the highest with 34.48%, and 6.07%. What’s more, lefse with LDA of 3.21 indicated that Proteobacteria were the predominant group in PT, while, Gemmatimonadetes enriched in PTVN at the phylum level, (Fig. 3A). At the genus level, the average relative abundance of soil bacterial community more than 1% were Subgroup_6, RB41, Sphingomonas, Blastococcus, Pseudonocardia, Solirubrobacter, Gemmatimonas, 67-14, Bradyrhizobium, Subgroup_7, TRA3-20, Actinoplanes and Microvirga (Fig. S1A).

For soil fungal community, Ascomycota, Basidiomycota, and Mortierellomycota were the predominant groups at the phylum level, accounting for 99.03%, 99.14%, 98.61% and 99.37% in PTAU, PTRP, PT, and PTVN, respectively (Fig. 2B). PTAU hold the highest Ascomycota with the relative abundance of 84.77%, and the lowest Basidiomycota with the relative abundance of 11.37%. In addition, lefse with LDA of 3.06 demonstrated that Ascomycota was the predominant group in PTAU, while, Basidiomycota was the predominant group in PTRP (Fig. 3B; Fig. S2B). At the genus level, the average relative abundance with more than 1% were Knufia, Penicillium, Trechispora, Solicoccozyma, Tomentella, Gibberella, Talaromyces, Didymella, Mortierella, Cladophialophora, Geminibasidium, Fusarium, Sebacina, and Exophiala (Fig. S1B).

More interestingly, the NMDS based on Jaccard distance and the permanova demonstrated that soil bacterial community (stress = 0.08) (Fig. 4A; Table S2) and fungal community (stress = 0.13) (Fig. 4B; Table S3) among PT, PTAU, PTRP and PTVN showed clear difference. Soil microbial community in PTAU clearly separated from PTAR, PTVN, and PT along NMDS1, indicating that PTAU distinctly changed the soil microbial community structure (Fig. 4).

The relationships between soil microbial community and soil characteristics

In total, 8 chemical attributes, including pH, TC, TN, C/N, TP, AN, AP, and AK, were selected from all the physicochemical attributes measured in this study (Table 1). The RDA demonstrated that the first and second axis total explained 98.2% and 82.7% of the variance in the species environment relationship at the bacterial phylum and genus level, respectively, and 99.7% and 67.1% of the variance could be explained by the first and second axes at the fungal phylum and genus level, respectively, which indicated that these axes reflected the relationship between soil microbial communities and the environment to a large extent (Fig. 5).

In addition, Proteobacteria significantly declined with TN (r = −0.63, P < 0.01) and AN (r = −0.53, P < 0.05), while increased with soil C/N (r = 0.77, P < 0.01). Gemmatimonadetes existed negative relation with the TP (r = −0.62, P < 0.05). Chloroflexi appeared to be the significant positive correlation with soil pH (r = 0.62, P < 0.05), TN (r = 0.80, P < 0.01), TP (r = 0.50, P < 0.05), AN (r = 0.74, P < 0.01), and AK (r = 0.59, P < 0.05), while existed obviously negative relation with C/N (r = −0.60, P < 0.05). With regard to soil fungi, the relative abundance of Mortierellomycota increased with TN (r = 0.59, P < 0.05) and AN (r = 0.65, P < 0.01) (Table 4).

Soil characteristics response to different vegetation types

In some areas, certain successes have been achieved in the near-natural transformation of artificial forests through selective logging and creating gaps (Luo et al., 2013). As anticipated, we observed that forest type conversion altered several of the measured soil attributes, including soil pH, TC, TN, C/N, TP, AN, and AK. What’s more, PTAU, PTVN, and PTRP could increase the soil pH, TC, TN, AN, and AK in different degrees, especially in PTAU. Previous study also showed that the soil pH, TC, TN, and TP concentrations in the P. tabuliformis forest were obviously lower than in a coniferous-broadleaved P. tabuliformis mixed forest (Chen et al., 2017). What’s more, the findings from Li et al. (2021a) demonstrated that soil pH in the monoculture P. massoniana forest plot existed abundantly lower than that in the broad-leaved tree species plantation. Furthermore, different aboveground vegetation types could affect the soil C/N ratio due to a change in litter composition and the microbial community structure were well documented (Zhang et al., 2013; Djukic et al., 2010). Our present study revealed that soil C/N in PT was significantly higher than in PTAU, PTVN, and PTRP, which was agreement with previous study (Monkai et al., 2018). All these findings demonstrated that trees can also change soil chemical parameters mostly via changes in litter quantity and quality and root exudates (Umemura & Takenaka, 2015).

Soil microbial community diversity and composition response to different vegetation types

Our study revealed the impact of near-natural transformation on the diversity and composition of soil bacteria and fungi. Consistent with our hypothesis, near-natural forest transformation could improve soil bacterial Pielou_e index, and Simpson index, as well as soil fungal Chao1 index, which was consistent with previous study (Ding et al., 2021). In addition, we indicated that soil bacterial and fungal community composition among PT, PTAU, PTRP and PTVN showed clear variations, and PTAU visibly altered the soil microbial community structure, and our findings corroborated the previous reports (Sun et al., 2020; Jin et al., 2019; Liu et al., 2018). Collectively these findings nicely demonstrated that distinct difference of soil microbial communities is usually found to be associated with tree species (Meng et al., 2019).

In present study, Proteobacteria was the most abundant phylum in all soil samples under different vegetation types, and similar results have been repeatedly observed in past studies (Zhang & Lv, 2021; Deng et al., 2018; Zhao et al., 2012). While, Liu et al. (2013) reported that Acidobacteria was the predominant flora of soil bacterial communities, and the relative abundance of Acidobacteria ranked third in our study. Bacterioidetes, Gemmatimonadetes, Planctomycetes, Verrucomicrobia and Chloroflexi were also rich in all soil samples, which were basically consistent with previous studies on soil bacterial community in the Loess Plateau (Du et al., 2020; Yang, Wang & Shen, 2013). Among the dominant phyla, the relative abundance of Proteobacteria was the highest in PT, and PTVN hold the lowest relative abundance, which was in line with previous finding (Deng et al., 2018) indicated that the higher soil nutrients may be related to the higher relative abundance of Proteobacteria. Proteobacteria have been recognized as copiotrophic taxa (Davis, Sangwan & Janssen, 2010), which could be used as an indicator of soil nutrient status. In general, compared to oligotrophic taxa, copiotrophic taxa have higher N demands (Fontaine & Barot, 2005; Fierer, 2017). However, our reports documented that PT with lower nutrients N and a higher C/N ratio was beneficial to the survival of Proteobacteria. Actinobacteria, as important decomposers of cellulose and chitin (Pankratov, Dedysh & Zavarzin, 2006), is found to be abundant in diverse environments (Urich et al., 2008), which play a fundamental role in the turnover and circulation of organic matter. In our study, no significant difference of Actinobacteria were obtained among different vegetation types, while, findings from Li et al. (2021a) reported that the relative abundance of Actinobacteria in the mixed broadleaf-conifer forest presented lower than pure P. massoniana plantation. With regard to soil fungi, Ascomycota, Basidiomycota, and Mortierellomycota were the predominant groups at the phylum level (Fig. 2B), which appeared to be consistent with previous report documented that Ascomycota was the predominant fungal phylum (Deng et al., 2020). Basidiomycota, as one of the important components of functional group, participates in decomposing recalcitrant plant materials as well as plays an important role in biomass conversion and nutrient cycling (Morrison et al., 2016). Ascomycota, as a decomposer in soil, is a common fungal community in forest ecosystem and involves in many degradation processes (Riley et al., 2014), supporting the findings of earlier study (Li et al., 2021b).

Intriguingly, we found that Subgroup_6, RB41, and Sphingomonas were the dominant bacterial genera, as previously reported in the literature (Shi et al., 2017). Sphingomonas, as one of the abundant microbial groups in soil environment, played a crucial role in the determination of soil nutrient recycling. In addition, Gemmatimonas was a dominant genus in soil in our study, which had the ability of weathering minerals and solubilizing phosphate in various environments (Wei et al., 2019). Moreover, some bacterial genera, such as Pseudonocardia, Biastococcus, RB41 and Sphingomonas were also notably different among different vegetation types. Previous investigations indicated that soil bacterial community were strongly influenced by different vegetation types (Mendes et al., 2015), which was confirmed again in our research.

The relationship between soil characteristics and soil microbial communities

Soil characteristics, including soil TC, TN and TP, alter soil microbial community have been well documented (Zhao et al., 2019; Deng et al., 2018; Mendes et al., 2015). In present study, mixed coniferous and broadleaf forest increased soil bacterial Pielou_e index by increasing soil pH and TN (Table 3), as previous studies have shown that mixed broadleaf-conifer forest increased soil bacterial community diversity by increasing soil N, P and other nutrients in comparison to the monoculture P. massoniana (Li et al., 2021a). While, soil bacterial Pielou_e index, Shannon index, and Simpson index increased by decreasing soil C/N (Table 3). Nevertheless, no significant correlations of soil fungal community diversity with soil characteristics were observed (Table 3). All of these reports suggested that aboveground vegetation strongly influenced belowground microbial community diversity via soil chemical properties (Osborne et al., 2011).

In our study, the RDA results also showed that the soil microbial community composition were significantly correlated with soil characteristics (Fig. 5). Our results were in partial agreement with the findings of Yang et al. (2018), Rousk et al. (2010), and Li et al. (2021a) demonstrated that soil pH was a factor controlling the bacterial community structure. In the analysis of the relationships between soil chemical properties and the dominant bacterial communities, we found that Proteobacteria significantly declined with TN and AN, while increased with soil C/N. Gemmatimonadetes existed negative relation with the TP, however, previous study demonstrated that the relative abundance of Gemmatimonadetes was positively influenced by TC and SOM (Liu et al., 2020), suggesting that relatively high soil nutrients might stimulate the increase of the relative abundance of Gemmatimonadetes (Rosling et al., 2011). Chloroflexi appeared to be the significant positive correlation with soil pH, TN, TP, AN, and AK, similar with existed observation (Liu et al., 2020), while existed obviously negative relation with C/N, which was consistent with the investigations of the study from Wan et al. (2015). What’s more, other investigations also indicated that soil pH, TC and TN have impacts on the composition of bacterial communities (Cho, Kim & Lee, 2016; Barcenas-Moreno, Bååth & Rousk, 2016; Dawud et al., 2016). These findings collectively supported the idea that soil microorganisms were strongly linked with soil nutrients caused by different tree species (Siles & Margesin, 2016; Mueller et al., 2014; De Vries et al., 2012). Our research results revealed the influence of tree composition and soil nutrient availability on soil microbial diversity and composition, provided management guidance for incorporating soil microbial diversity and composition in plantation forest management.