Dynamics of Soil Nutrients, Microbial Community Structure, Enzymatic Activity, and Their Relationships along a Chronosequence of Pinus massoniana Plantations

Pinus massoniana is the major afforestation and vegetation restoration tree in southern China, and it plays an important role in the sustainable development of plantations. However, longterm single planting of P. massoniana has resulted in the decline of soil quality and forest productivity, and a soil fertility assessment is urgently needed. We selected P. massoniana plantations of four age stages for plot investigation and sampling to determine the soil physicochemical properties, microbial diversity and composition, and enzyme activities at different soil depths. The results showed that soil total phosphorus (TP) and available phosphorus (AP) decreased with the increase of age, especially low C/N ratio and high C/P and N/P ratio in the 30-year and 36-year stands, leading to P limitation. Meanwhile, the bacterial Shannon index also decreased with the increase of age and was positively correlated with AP, NO3-N, and pH. However, the fungal Shannon index decreased first and then increased with the increase of age; soil acid phosphatase (S-ACP) and urease activities showed a similar trend. Correlation analysis demonstrated that the increase of total organic carbon (TOC) and total nitrogen (TN) promoted the increase of fungal Shannon index, which was beneficial to the secretion of more enzymes. We found that soil physicochemical properties, microbial diversity, and enzyme activity decreased simultaneously when soil depths increased. Moreover, Acidobacteria and Basidiomycota were the most abundant bacterial and fungal communities, respectively, followed by Proteobacteria and Actinobacteria for bacteria and Ascomycota for fungi, and these microbial taxa were significantly affected by soil water content (SWC), TOC, AP, and C/P. In conclusion, this work reveals the potential correlation among soil physicochemical properties, microbial diversity and composition, and enzyme activities, and revealed potential correlations among them which will help to improve understanding of soil conditions and provide a reference for rational management of soil resources.


Introduction
Pinus massoniana is the main tree species for afforestation and vegetation rehabilitation in Southern China because of its special advantages such as fast growth, high adaptability to drought and barren soils, and its capacity to retain water and nutrients. According to a national survey, the forest area of the dominant species of P. massoniana in China is about 10 million ha, accounting for 6.08% of the total forest area in the country and 20% of the plantation forest areas in Southern China [1,2]. Consequently, P. massoniana occupies an extremely important position in promoting wood production and sustainable development

Experimental Design
The sample plot setting was carried out on 20-25 August 2017. The four stands had similar environmental characteristics, such as altitude, slope aspect, slope, slope position, and soil texture. Three 20 m × 20 m typical plots were randomly set in different age stands (12 sample plots in total). The diameter at breast height (DBH), tree height, height under the branch, canopy density, altitude, slope aspect, slope, slope position, and soil parameters of P. massoniana in the sample plot were measured (Table 1).

Experimental Design
The sample plot setting was carried out on 20-25 August 2017. The four stands had similar environmental characteristics, such as altitude, slope aspect, slope, slope position, and soil texture. Three 20 m × 20 m typical plots were randomly set in different age stands (12 sample plots in total). The diameter at breast height (DBH), tree height, height under the branch, canopy density, altitude, slope aspect, slope, slope position, and soil parameters of P. massoniana in the sample plot were measured (Table 1).

Soil Samples Collection
Soil sample collection was carried out on 20 June 2019. Five sampling points were set on the diagonal of each 20 m × 20 m plot to collect soil samples at depths of 0-20, 20-40, and 40-60 cm, respectively. Then, the debris such as stones, litter, and roots were removed, and the five different sampling points with the same soil depth were mixed. Each plot obtains representative samples of three soil depths (using the quartering method, taking 500 g soil for each soil depth). A total of 36 soil samples (4 ages × 3 plots × 3 soil depths) were obtained in this study. Each soil sample was divided into two parts. One-part samples were placed in sterile bags, and these samples were cooled immediately until arrival at the laboratory, and then they were subsampled and frozen at −80 • C as soon as possible until they were used to extract soil microbial DNA. The other part was used to determine the physicochemical properties and enzyme activities of soil.

Soil Physicochemical Properties and Enzyme Activity Measurements
The collected soil samples were air dried and crushed in the laboratory to determine the soil physicochemical properties and enzyme activities. The total organic carbon (TOC) concentration of soil was determined using the K 2 Cr 2 O 7 oxidation method. After digestion with H 2 SO 4 -H 2 O 2 , the concentrations of total nitrogen (TN) and total phosphorus (TP) were determined by the Kjeldahl and Mo-Sb colorimetric assay methods, respectively [32]. Soil ammonium (NH 4 + -N) was extracted by 2 mol L −1 KCl and determined by the Indigo phenol blue method [33]. Nitrate nitrogen (NO 3 --N) was determined by ultraviolet spectrophotometry [34]. Available phosphorus (AP) was extracted by NH 4 F, followed by the colorimetric method [35]. The pH value of soil was determined by the glass composite electrode method (water soil ratio of 2.5: 1) [20]. Soil water content (SWC), temperature (ST), and electrical conductivity (SEC) were measured by TDR350 soil moisture, temperature, and electrical conductivity meter. Soil acid phosphatase (S-ACP) was determined by diphenyl sodium phosphate colorimetry, sucrase by 3,5-Dinitrosalicylic acid colorimetry, and urease by indophenol colorimetry [36]. Soil DNA extraction, PCR amplification, and Illumina sequencing were completed by Hangzhou Kaitai Biotechnology Co., Ltd. (Hangzhou, China). The soil sample weighed 0.25 g. FastDNA ® Spin Kit (MP Biomedicals, Santa Ana, CA, USA) was used to extract the total DNA in the soil sample according to the manufacturer's instructions. The DNA concentration and purity were detected using the Thermo NanoDrop 2000 UV spectrophotometer and 1% agarose gel electrophoresis. The tested samples were stored at −80 • C for further analysis. The v3 + v4 region of the bacterial 16S rDNA gene and V4 region of the fungal 18S rDNA gene were amplified by StepOne Plus fluorescent quantitative PCR. The bacterial primers were 341F (5 -CCTACGGGNGGCWGCAG-3 ) and 805R (5 -GGACTACHVGGGTWTCTAA-3 ) [37]. The fungal primers were EF4F (5 -GGAAGGGRTGTATTTATTAG-3 ) and EF3R (5 -TCCTCTAAATGACCAAGTTTG-3 ) [38]. The PCR reaction was performed in a 20.0 µL mixture, which contained 1.0 µL of DNA template (10 ng/µL), 0.4 µL of forward primer, 0.4 µL of reverse primer, 0.4 µL of ROX Reference Dye, 10.0 µL of 2 × SYBR ® Premix EX Taq™, and 7.8 µL of water without RNA enzyme (Rnase-Free Water). The PCR amplification procedure was as follows: per-denaturation at 95 • C, 3 min, 1 cycle; denaturation at 95 • C, 30 s, 25 cycles; annealing at 55 • C, 30 s; extension at 72 • C, 30 s; final extension at 72 • C, 5 min, 1 cycle; and 4 • C storage. Library quality inspection was performed by Thermo NanoDrop 2000 UV spectrophotometer and 2% agarose gel electrophoresis. After the quality inspection was qualified, the library was quantified by KAPA Library Quantification Kit Illumina ® Platforms (KAPA Biosystems, Inc., Wilmington, NC, USA). Finally, an equal amount of PCR product from each sample was put into the single tube, and then high-throughput sequencing was performed using Illumina Miseq platform of Hangzhou Kaitai Biotechnology Co., Ltd. (Hangzhou, China).

Sequencing Data Processing
The two ends of the library were sequenced on the Illumina MiSeq. Each sample generated two files (reads1 and reads2), which were stored in a low-quality fastq format. Conducting quality control and filtering on the original offline data is necessary to obtain high-quality sequencing data and improve the accuracy of biological information analysis [39]. The main steps were as follows: (1) the fastq_mergepairs command of Vsearch software was used to assemble the paired reads from double terminal sequencing into a single sequence; (2) Cutadapt software was used to remove PCR primers from the sequence; (3) the fastq_filter command Vsearch software was used to remove the base sequence containing N (fuzzy base), the same base number of more than 8, and length of less than 100 bp.

Operational Taxonomic Unit (OTU) Clustering and Notes
To reduce the amount of redundant calculation in the process of analysis, we counted the repetitive sequence and extracted the non-repetitive sequence. We used Vsearch software to remove the once-detected sequences (Singleton sequences). The sequences with more than 97% similarity were gathered into one operational taxonomic unit (OTU), the chimeric sequence was removed, and the OTU representative sequence was used for subsequent species annotation. We used the Ribosomal Database Project (RDP) [40] and SLIVA databases [41] to classify and identify bacteria and fungi, respectively. The raw sequencing reads are available at the Sequence Read Archive (SRA) database of the National Center for Biotechnology Information (NCBI) under accession numbers SRP291512: PRJNA675184.

Statistical Analysis
Based on the OTU information of soil samples, the alpha diversity (Shannon index, H ) of bacteria and fungi in soil samples was calculated by using the software of mothur (version v.1.30.1) [10]. The Shannon index was used to estimate the microbial diversity in the soil samples, and the larger the Shannon index, the higher the diversity. Two-way ANOVA was used to examine the effects of stand age and soil depth on soil microbial community composition and diversity and elemental stoichiometry. Significant differences were measured by one-way ANOVA followed by an LSD test and significance was determined at the 0.05 level. Pearson correlation was used to determine the relationship between soil microbial diversity, enzyme activity, and soil physicochemical properties. These statistical analyses were carried out on SPSS 21.0 software (SPSS, Somers, NY, USA) and the figures were mapped on Origin 2020 (Origin Lab Corp., Northampton, MA, USA). Redundancy analysis (RDA) was used to study the effect of soil element stoichiometric variables on soil microbial community composition and implemented on Canoco 5.0 software (Microcomputer Power, Inc., Ithaca, NY, USA). The formula of Shannon index (H ) is as follows:  The two-way ANOVA showed that different stand ages and soil depths had significant effects on soil TOC, TN, and TP concentrations and their stoichiometry ( Table 2, p < 0.05 or p < 0.01). TOC and TN concentrations in different soil depths showed different trends with the increase of stand age. That is, the 0-20 and 20-40 cm soil depths showed that 30-year stand was significantly higher than other stands, while the 40-60 cm soil depth showed that 20-year and 30-year stands were significantly lower than the other stands ( Figure 2  decreased with the increase of stand age and peaked in 20-year stand (0.62 ± 0.04 g kg −1 , 0.45 ± 0.02 g kg −1 , and 0.40 ± 0.03 g kg −1 , respectively). C/N ratio of 0-20 and 40-60 cm soil depths first increased, then decreased with the increase of stand age and peaked in 30-year stand (6.32 ± 0.14 and 7.65 ± 0.49, respectively), while C/N ratio of 20-40 cm soil depth tended to decrease ( Figure 2). The changing trend of C/P and N/P ratios in the three soil depths were opposite to TP concentration with the increase of stand age. Except for the C/N ratio, the soil nutrient concentration and their stoichiometry of 0-20 cm soil depth in four stands were significantly higher than other soil depths ( Figure 2, p < 0.05).
The two-way ANOVA showed that different stand ages and soil depths had significant effects on soil TOC, TN, and TP concentrations and their stoichiometry (Table 2, p < 0.05 or p < 0.01). TOC and TN concentrations in different soil depths showed different trends with the increase of stand age. That is, the 0-20 and 20-40 cm soil depths showed that 30-year stand was significantly higher than other stands, while the 40-60 cm soil depth showed that 20-year and 30-year stands were significantly lower than the other stands ( Figure 2, p < 0.05). TP concentrations of 0-20, 20-40, 40-60 cm soil depths increased first, then decreased with the increase of stand age and peaked in 20-year stand (0.62 ± 0.04 g kg −1 , 0.45 ± 0.02 g kg −1 , and 0.40 ± 0.03 g kg −1 , respectively). C/N ratio of 0-20 and 40-60 cm soil depths first increased, then decreased with the increase of stand age and peaked in 30-year stand (6.32 ± 0.14 and 7.65 ± 0.49, respectively), while C/N ratio of 20-40 cm soil depth tended to decrease ( Figure 2). The changing trend of C/P and N/P ratios in the three soil depths were opposite to TP concentration with the increase of stand age. Except for the C/N ratio, the soil nutrient concentration and their stoichiometry of 0-20 cm soil depth in four stands were significantly higher than other soil depths ( Figure 2, p < 0.05).   Different stand ages had significant effects on soil physicochemical properties (p < 0.05 or p < 0.01). Different soil depths had no significant effects on NH 4 + -N, NO 3 --N, pH, and SEC (p > 0.05), but had significant effects on AP, SWC, and ST (Table 2, p < 0.01). NH 4 + -N and ST of three soil depths tended to rise with the increase of stand age and reached the highest in 36-year stand (the mean values were 12.59 ± 0.58 mg kg −1 and 20.01 ± 0.37 • C, respectively). On the contrary, except for 40-60 cm soil depth, NO 3 --N, AP, and pH tended to decrease with the increase of stand age, and the maximum appeared in 10-year stand (the mean values were 3.51 ± 1.56 and 7.15 ± 2.06 mg kg −1 , and 4.60 ± 0.05, respectively) ( Figure 3). SWC and SEC showed irregular changes, and they decreased first and then increased in 0-20 cm soil depth with the increase of stand age, which was the lowest in 20 (20.00 ± 2.46%) and 30 years (32.00 ± 2.00 µs cm −1 ), respectively. AP and ST of four stands decreased significantly with the increase of soil depth, while SWC increased significantly ( Figure 3, p < 0.05).   Different stand ages had significant effects on soil physicochemical properties (p < 0.05 or p < 0.01). Different soil depths had no significant effects on NH4 + -N, NO3 --N, pH, and SEC (p > 0.05), but had significant effects on AP, SWC, and ST ( Table 2, p < 0.01). NH4 + -N and ST of three soil depths tended to rise with the increase of stand age and reached the highest in 36-year stand (the mean values were 12.59 ± 0.58 mg kg −1 and 20.01 ± 0.37 °C, respectively). On the contrary, except for 40-60 cm soil depth, NO3 --N, AP, and pH tended to decrease with the increase of stand age, and the maximum appeared in 10-year stand (the mean values were 3.51 ± 1.56 and 7.15 ± 2.06 mg kg −1 , and 4.60 ± 0.05, respectively) ( Figure 3). SWC and SEC showed irregular changes, and they decreased first and then increased in 0-20 cm soil depth with the increase of stand age, which was the lowest in 20 (20.00 ± 2.46%) and 30 years (32.00 ± 2.00 µs cm −1 ), respectively. AP and ST of four stands decreased significantly with the increase of soil depth, while SWC increased significantly ( Figure 3, p < 0.05).

The Variation of Soil Enzyme Activities with Stand Age and Soil Depth
The two-way ANOVA showed that different stand ages had significant effects on S-ACP, sucrase, and urease (p < 0.01). Different soil depths had significant effects on soil acid phosphatase and urease (p < 0.01), while they had no significant effects on sucrase ( Table 2, p > 0.05). On the whole, except for urease of 30-year stand at 40-60 cm, soil depth was significantly higher than other stands (p < 0.05); S-ACP and urease decreased first and then increased with the increase of stand age. S-ACP and urease of 30-year stand at 0-20 cm soil depth was significantly lower than other stands, and S-ACP of 20-year stand at 20-40 and 40-60 cm soil depths and urease of 20-year stand at 20-40 cm soil depth was significantly lower than other stands (Figure 4, p < 0.05). The sucrase of 0-20 and 20-40 cm soil depth tended to rise with the increase of stand age, but there was no significant difference between different stands at the 20-40 cm soil depth (p > 0.05). The sucrase of 20-year stand at 40-60 cm soil depth was significantly lower than other stands ( Figure 4, p < 0.05). Moreover, S-ACP and urease of 40-60 cm soil depth in four stands were significantly higher than other depths (p < 0.05).

The Variation of Soil Enzyme Activities with Stand Age and Soil Depth
The two-way ANOVA showed that different stand ages had significant effects on S-ACP, sucrase, and urease (p < 0.01). Different soil depths had significant effects on soil acid phosphatase and urease (p < 0.01), while they had no significant effects on sucrase ( Table 2, p > 0.05). On the whole, except for urease of 30-year stand at 40-60 cm, soil depth was significantly higher than other stands (p < 0.05); S-ACP and urease decreased first and then increased with the increase of stand age. S-ACP and urease of 30-year stand at 0-20 cm soil depth was significantly lower than other stands, and S-ACP of 20-year stand at 20-40 and 40-60 cm soil depths and urease of 20-year stand at 20-40 cm soil depth was significantly lower than other stands (Figure 4, p < 0.05). The sucrase of 0-20 and 20-40 cm soil depth tended to rise with the increase of stand age, but there was no significant difference between different stands at the 20-40 cm soil depth (p > 0.05). The sucrase of 20-year stand at 40-60 cm soil depth was significantly lower than other stands (Figure 4, p < 0.05). Moreover, S-ACP and urease of 40-60 cm soil depth in four stands were significantly higher than other depths (p < 0.05).

Alpha Diversity of Soil Bacteria and Fungi
The total quantity of bacterial and fungal OTU obtained from 36 samples was 78,132 and 13,230, respectively, and the number of bacterial and fungal OTU in each sample was 1780 to 2461 and 151 to 528, respectively ( Figure 5). Based on 97% similarity, the rarefaction curves for bacteria and fungi tended to saturate for each sample ( Figure S1), which showed that the sequencing data were reasonable, and the majority of sequence types could be detected. The two-way ANOVA showed that different stand ages had significant effects on bacterial and fungal OTU and Shannon index (p < 0.05). Different soil depths had significant effects on bacterial Shannon index and fungal OTU and Shannon index (p < 0.05), while they had no significant effects on bacterial OTU ( Table 2, p > 0.05). The bacterial and fungal Shannon index had the same trend with OTU, that is, the bacterial alpha diversity decreased with the increase of stand age, while the fungal alpha diversity first decreased and then increased ( Figure 5). The bacterial Shannon index of 0-20, 20-40, and 40-60 cm soil depths were the largest in 10-year stand, which were 6.32 ± 0.08, 6.20 ± 0.10, and 6.14 ± 0.03, respectively. The fungal Shannon index of 20-40 and 40-60 cm depths were the smallest in the 30-year stand, which were 1.14 ± 0.18 and 0.95 ± 0.28, respectively ( Figure 5). Furthermore, the bacterial and fungal alpha diversity of 0-20 cm soil depth was significantly higher than that of the other two soil depths ( Figure 5). rial alpha diversity decreased with the increase of stand age, while the fungal alpha diversity first decreased and then increased ( Figure 5). The bacterial Shannon index of 0-20, 20-40, and 40-60 cm soil depths were the largest in 10-year stand, which were 6.32 ± 0.08, 6.20 ± 0.10, and 6.14 ± 0.03, respectively. The fungal Shannon index of 20-40 and 40-60 cm depths were the smallest in the 30-year stand, which were 1.14 ± 0.18 and 0.95 ± 0.28, respectively ( Figure 5). Furthermore, the bacterial and fungal alpha diversity of 0-20 cm soil depth was significantly higher than that of the other two soil depths ( Figure 5).

Community Composition of Bacteria and Fungi
A total of 9 and 3 phyla of bacteria and fungi with a relative abundance of more than 1% were identified (phyla with more than 1% average abundance in at least one sample were considered as predominant) (Figure 6a,b). The average relative abundances of Acidobacteria (31.12%), Proteobacteria (16.33%), Actinobacteria (17.06%), Bacte-

Community Composition of Bacteria and Fungi
A total of 9 and 3 phyla of bacteria and fungi with a relative abundance of more than 1% were identified (phyla with more than 1% average abundance in at least one sample were considered as predominant) (Figure 6a,b). The average relative abundances of Acidobacteria (31.12%), Proteobacteria (16.33%), Actinobacteria (17.06%), Bacteria-unclassified (16.02%), and Chloroflexi (11.11%) were all higher than 10%, which were dominant in the soil bacterial community. The relative abundances of Planctomycetes, Verrucomicrobia, candidate division WPS-2, and candidate division WPS-1 were all higher than 1% (Figure 6a), and they were also important members of soil bacteria. The trend of different bacterial relative abundance with the increase of stand age and soil depth was different. Acidobacteria and Proteobacteria first decreased and then increased with the increase of stand age, while Actinobacteria and Bacteria-unclassified were the opposite. Acidobacteria showed no obvious rule with the increase of soil depth, Proteobacteria and Actinobacteria decreased, while Bacteria unclassified and Chloroflexi increased (Figure 6a). The average relative abundance of Basidiomycota (83.34%) in the soil fungal community was the highest, followed by Ascomycota (8.96%) and Eukaryota-unclassified (6.71%). The relative abundance of Basidiomycota increased with the increase of stand age and soil depth, but Ascomycota and Eukaryota-unclassified were the opposite (Figure 6b).
The heat map showed that bacterial and fungal communities as 20-year, 30-year, and 36-year stands were clustered together, clearly separated from the 10-year stand. The 20-40 cm and 40-60 cm soil depths were clustered together, clearly separated from 0-20 cm soil depth (Figure 7). ta-unclassified (6.71%). The relative abundance of Basidiomycota increased with the increase of stand age and soil depth, but Ascomycota and Eukaryota-unclassified were the opposite. (Figure 6b).
The heat map showed that bacterial and fungal communities as 20-year, 30-year, and 36-year stands were clustered together, clearly separated from the 10-year stand. The 20-40 cm and 40-60 cm soil depths were clustered together, clearly separated from 0-20 cm soil depth (Figure 7). Relative abundance (%) Figure 6. Relative abundance distribution of bacterial and fungal communities at the phylum level at different soil depths and stand ages. Phyla with more than 1% average abundance in at least one sample were considered as predominant. Figure (a,b) show the variation of relative abundance of bacteria and fungi, respectively. Other: the taxa of bacteria or fungi with relative abundance less than 1%; 10 Y: 10

Relationship between Microbes, Enzyme Activity, and Soil Physicochemical Properties
Pearson correlation analysis showed that soil enzyme activities and the bacterial and fungal Shannon index were significantly correlated with soil physicochemical properties ( Figure 8, p < 0.05 or p < 0.01). S-ACP was positively correlated with TOC, TN and N/P ratios, pH, AP, and ST, as well as negatively correlated with C/N ratio and SWC; sucrase was positively correlated with C/P and N/P ratios, ST, and negatively correlated with SC; urease was positively correlated with TN, pH, AP, and ST, and negatively correlated with C/N ratio (Figure 8). Bacterial Shannon index was positively correlated with pH, NO 3 --N, SC, and negatively correlated with NH 4 + -N, and SWC. Fungal Shannon index was positively correlated with TOC, TN, AP, and ST. Besides, bacterial and fungal Shannon index were significantly positively correlated with S-ACP and urease (p < 0.01), but not significantly correlated with sucrase (p > 0.05) (Figure 8). In conclusion, S-ACP, urease, bacterial, and fungal Shannon index were significantly positively correlated with AP (p < 0.01), but negatively correlated with SWC.
Redundancy analysis was conducted with the abundance of soil bacteria and fungi as response variables and soil physicochemical properties as explanatory variables (environmental variables). The results showed that the first two axes of the RDA plot explained 40.56% and 9.36% (49.92% in total) of the total microbial variation (Figure 9). A Monte Carlo permutation test showed that the SWC, TOC, AP, and C/P ratio were significantly correlated (p < 0.05) with soil microbial community composition in all environmental variables considered. In particular, the SWC was the first factor that contributed to 28.2% of total variation (Table 3). Figure 9 show that SWC was positively correlated with Planctomycetes, Chloroflexi, Verrucomicrobia, and Basidiomycota, while negatively correlated with Acidobacteria, Actinobacteria, Proteobacteria, Bacteria-unclassified, candidate division WPS-1, candidate division WPS-2, Ascomycota, and Eukaryota-unclassified. Conversely, TOC, AP, and C/P ratio were negatively correlated with Planctomycetes, Chloroflexi, Verrucomicrobia, and Basidiomycota, while positively correlated with Acidobacteria, Bacteriaunclassified, Actinobacteria, Proteobacteria, candidate division WPS-1, candidate division WPS-2, Ascomycota, and Eukaryota-unclassified. Heat map depicting the hierarchy cluster results for the abundance of bacteria and fungi at the phylum level. Phyla with more than 1% average abundance in at least one sample were considered as predominant. Figure (a,b) show the clustering of bacteria and fungi at different ages, respectively. Figure (c,d) show the clustering of bacteria and fungi at different soil depths, respectively. Red indicates high relative abundance and green indicates low relative abundance. 0-20 cm: 0-20 cm soil depth; 20-40 cm: 20-40 cm soil depth; 40-60 cm: 40-60 cm soil depth.

Relationship between Microbes, Enzyme Activity, and Soil Physicochemical Properties
Pearson correlation analysis showed that soil enzyme activities and the bacterial and fungal Shannon index were significantly correlated with soil physicochemical properties (Figure 8, p < 0.05 or p < 0.01). S-ACP was positively correlated with TOC, TN and N/P ratios, pH, AP, and ST, as well as negatively correlated with C/N ratio and SWC; sucrase was positively correlated with C/P and N/P ratios, ST, and negatively correlated with SC; urease was positively correlated with TN, pH, AP, and ST, and negatively correlated atively correlated with Acidobacteria, Actinobacteria, Proteobacteria, Bacteria-unclassified, candidate division WPS-1, candidate division WPS-2, Ascomycota, and Eukaryota-unclassified. Conversely, TOC, AP, and C/P ratio were negatively correlated with Planctomycetes, Chloroflexi, Verrucomicrobia, and Basidiomycota, while positively correlated with Acidobacteria, Bacteria-unclassified, Actinobacteria, Proteobacteria, candidate division WPS-1, candidate division WPS-2, Ascomycota, and Eukaryota-unclassified.

Effects of Soil Physicochemical Properties on Microbial Alpha Diversity
The results were in line with our first expectation that soil bacterial and fungal alpha diversity (Shannon index) showed different trends with the increase of the chronosequence of P. massoniana plantations (Figure 5), and the bacterial Shannon index was significantly correlated with soil AP, NO3 --N, pH, and SWC (p < 0.05), while fungal Shannon index was significantly correlated with TOC, TN, AP, and ST (p < 0.05) (Figures 8). Pre-

Effects of Soil Physicochemical Properties on Microbial Alpha Diversity
The results were in line with our first expectation that soil bacterial and fungal alpha diversity (Shannon index) showed different trends with the increase of the chronosequence of P. massoniana plantations (Figure 5), and the bacterial Shannon index was significantly correlated with soil AP, NO 3 --N, pH, and SWC (p < 0.05), while fungal Shannon index was significantly correlated with TOC, TN, AP, and ST (p < 0.05) (Figure 8). Previous studies have shown that soil physicochemical properties, root biomass, and chemical composition change constantly with the continuous development of forests, which inevitably leads to changes in the diversity and composition of the soil microbial community [6,21,42]. Our research revealed that the bacterial Shannon index tended to decline with the increasing age of the P. massoniana plantation (Figure 5), and a similar result has been reported in other studies [6,42]. In our study, the decline of bacterial Shannon index could be attributed to the gradual decrease of soil nutrient availability, especially AP and NO 3 --N during the aging process of the P. massoniana plantation. As N and P were essential nutrients for microbial growth, the change of their effectiveness affected microbial composition and diversity [43]. As we knew, P was a necessary component in the production of the ribosome, ATP, DNA, and RNA to maintain rapid growth [44]. In 0-20 cm and 20-40 cm soil depths, the TP concentration in near-mature and mature forests was significantly lower than that in young and middle forests, while C/P and N/P ratios were just the opposite, which may lead to P limitation and thus inhibit bacterial growth ( Figure 2). Moreover, fungi are more likely to thrive at a lower pH, while bacteria occupy the dominant position at a higher pH [45]. Our results showed that soil pH was significantly positively correlated with bacterial Shannon index (Figure 9), indicating that the decrease of pH in the near-mature and the mature forest was also the reason for the decrease of bacterial diversity. Generally, excessive water reduces soil porosity and restricts oxygen input and CO 2 emission in the process of aerobic respiration, and thus affects aerobic microbial activities [46]. Hence, SWC increases with the age of the P. massoniana plantation, which may inhibit aerobic bacterial growth and affect its diversity. However, we found that the fungal Shannon index decreased first and then increased with the increase of the P. massoniana plantation age ( Figure 5). This result indicates that the soil condition may be more helpful to the growth and reproduction of fungi with the increase of P. massoniana plantation stand age. It is generally believed that the availability and quality of soil C and nutrients regulate microbial biomass and activity, and the C/nutrients ratio can reflect the quality of soil microbial growth environment to some extent [20]. Correlation analysis showed that fungal Shannon index was significantly positively correlated with soil TOC and TN (Figure 8, p < 0.05), which confirmed the above conclusion. In the present study, the increase of TOC and TN in the near-mature and mature forests provided carbon sources and nutrients for fungi and promoted the increase of fungal diversity. On the other hand, the needles of P. massoniana contain refractory components (cellulose, lignin, tannin, etc.) and fungi can decompose the lignocellulosic matrix which other microbes cannot use to obtain more nutrients [47,48]. Thus, the increase of fungal alpha diversity may be related to the increase of litters of P. massoniana. Kubartova et al. [49] studied the fungal diversity of litter decomposition in temperate forests of southern France and found that the decomposition rate and fungal richness of the mixed litter were significantly higher than that of beech and spruce litters in the late stage of decomposition (24 months). This indicates that the mixed litter can provide diverse niches for microbes. Consequently, P. massoniana coniferous broadleaved mixed forest should be more conducive to the microbial decomposition of litters and improve soil fertility.
Different soil depths had significant effects on bacterial and fungal Shannon index ( Table 2), and they decreased significantly with the increase of soil depth ( Figure 5, p < 0.05), which agreed with our second hypothesis. Du et al. [50] studied variations in bacterial and fungal communities through soil depth profiles in a Betula albosinensis forest and similar results were obtained. The main reason is that the soil organic C, N, and AP decreased significantly with the increase of soil depth (Figures 2 and 3), and the direct carbon sources and nutrients required by soil bacteria and fungi were limited, resulting in the decrease of their diversity. Secondly, SWC increases with the increase of soil depth, which may change soil structure, affect soil respiration, and inhibit aerobic microbial activities, causing their diversity to decrease [46]. In addition, the input of litter and roots on the soil surface can provide more abundant carbon and nitrogen for the surface soil to increase the microbial diversity [14,51]. However, the nutrients of the deep soil depth are derived from the combinations of allocation of root secretions, translocations of mineral nutrients from top surface by leaching, and weathering of minerals [52]. We also found that the soil C/P and N/P ratios significantly decreased with the increase of soil depth, while the C/N ratio significantly increased (Figure 2), indicating that the deep soil microbial diversity may be limited by N and the surface soil (0-20 cm) microbial diversity may be limited by P [29]. With no P fertilizer application, the P source in the 0-20 cm soil is primarily from weathering soil containing P minerals. Even though the P content in soil minerals of the P. massoniana forest in this study area is high, the rate of P released from the weathering process is less than the rate of P demand by both trees and the microbial. In addition, phosphorus can easily be fixed by Al and Fe when soil pH is as low as in our studied sites (<4.8, Figure 3). Thus, more N fertilizer should be applied to the deep soil and more P fertilizers should be applied to the topsoil to promote the absorption of nutrients by plants and maintain microbial diversity.

Effects of Soil Physicochemical Properties on Microbial Composition
Cluster analysis showed that the bacteria and fungi were clearly separated at 10-year stand from the others, and they were clearly separated at 0-20 cm soil depth from the other soil depths (Figure 7). It indicates that the composition of soil microbial communities varied greatly with the change of stand age and soil depth. Our two hypotheses have been confirmed once again. RDA analysis indicated that SWC, TOC, AP, and C/P ratio were the main factors affecting the microbial community ( Figure 9 and Table 3), and different bacterial and fungal taxa had different responses to these soil physicochemical properties, which supported our hypothesis. Specifically, the abundance of Acidobacteria was the highest among all soil depths and stands, followed by Proteobacteria and Actinobacteria (Figure 6a). The same results were found in B. albosinensis forest by Du et al. [50]. However, Acidobacteria were positively correlated with soil TOC, AP, and pH in our results (Figure 9), which were different from most previous studies. Zhang et al. [53] studied the forest soil in Shennongjia and found that Acidobacteria_GP1, Acidobacteria_GP2, and Acidobacteria_GP3 were significantly positively correlated with soil pH, and other subgroups (excluding Acidobacteria_GP7) and were significantly negatively correlated with pH. Fierer et al. [54] found that the relative abundance of Acidobacteria and soil TOC and C/N ratio were significantly positively correlated. The above showed that different subpopulations of Acidobacteria had different adaptability to soil physicochemical properties, and not all Acidobacteria were acidophilic and oligotrophic. Proteobacteria and Actinobacteria decreased with soil deepening and were significantly positively correlated with soil TOC and AP (Figures 7a and 9), which was consistent with the changes in soil bacteria in the B. platyphylla forest [50]. In general, the distribution of Proteobacteria is affected by the content of organic carbon, which is related to the high utilization rate of carbon [54]. More importantly, Proteobacteria preferred the high P, because its rapid growth requires a large number of RNA to meet the rapid cell division [55]. Previous studies showed that Actinobacteria is a key litter decomposer in soil, and part species of this phylum can secrete extracellular enzymes to promote the degradation of lignin, polysaccharides, and phenolics [17,56]. Consequently, Actinobacteria in topsoil can promote nutrient cycling by decomposing litter. SWC was negatively correlated with Acidobacteria, Proteobacteria, and Actinobacteria, and was the main limiting factor (the explanation rate was 28.2%) ( Figure 9 and Table 3). As mentioned earlier, SWC mediates microbial activity by affecting soil aeration and nutrient availability, so high SWC may not be beneficial for microbial growth [13,46].
For fungal communities, Basidiomycota was the most abundant in all soil layers and stands, followed by Ascomycota and Eukaryota-unclassified (Figure 6b). We observed that the abundance of Basidiomycota increased as the soil depths increased, and similar results were reported in previous studies [50,57]. In this study, soil nutrients (e.g., TOC, TN, AP) tended to decrease with the deepening of soil, and they were negatively correlated with the abundance of Basidiomycota (Figures 2, 3 and 9). Therefore, the change of Basidiomycota abundance may possibly be attributed to other causes. In most cases, the majority genera of Basidiomycota formed a symbiotic relationship with plant roots, enabling them to obtain carbon sources directly from the roots [57]. Cheng et al. [58] found that the P. massoniana roots of 2-10 mm were mainly distributed in the 20-60 cm depth range. Thus, the Basidiomycota in deep soil (40-60 cm) can form mycorrhiza with P. massoniana roots, providing carbon and nutrients, which makes the abundance of Basidiomycota higher. Usually, Ascomycota is mostly saprophytes, which decompose litter to obtain nutrients, while only a few of them can form mycorrhiza with plant roots, so Ascomycota is mainly distributed in surface soil [59]. Moreover, the soil nutrient content decreased and SWC increased with the deepening of the soil (Figures 2 and 3), which may not be conducive to the survival of Ascomycota, resulting in the decline of its abundance. We also noticed that the abundance of Ascomycota in the topsoil (0-20 cm) increased with the age of P. massoniana plantation (Figure 6b). This may be related to the increase of undergrowth shrubs and herbs in the near mature forest and mature forest stage after thinning, which makes litter easier to decompose. However, the abundance of Ascomycota in 20-40 cm and 40-60 cm soil depths tended to decrease with age (Figure 6b), which could be explained by the increase of nitrogen concentration (especially NH 4 + -N) in these two soil depths with age, because high nitrogen was harmful to the phylum [6].

Relationship among Soil Physicochemical Properties, Microbes, and Enzyme Activity
The activities of several enzymes have different trends with the increase of the chronosequence of P. massoniana plantations, but they decreased synchronously with the soil depths, and were significantly correlated with soil physicochemical properties and microbial alpha diversity. These results are in accordance with our hypothesis. Generally, a high S-ACP level is essential for improving soil P availability [6]. Our study observed that the topsoil has higher P and S-ACP, and there was a strong positive correlation between them, which was consistent with the function of S-ACP. Strangely, S-ACP first decreased and then increased with increasing stand age while AP continued to decrease, which indicated that S-ACP was not only related to phosphorus but may also be influenced by other factors. Previous studies demonstrated that TOC and TN could provide sufficient substrates for enzyme synthesis [28], so the lower soil TOC and TN in the middle-aged forest stage and deep soil led to the decrease of enzyme activity. S-ACP and urease were positively correlated with soil TOC and TN, which confirmed the above results. At the same time, the improvement of soil pH and the increase of ST can also create a favorable condition for the enzyme to exist stably in the soil [9]. Correlation analysis showed that S-ACP and urease were significantly negatively correlated with soil C/N ratio and SWC; sucrase was significantly negatively correlated with SEC, which indicated that higher C/N, SWC, and SEC might not be conducive to enzyme activity. In addition, soil enzyme activity was not only related to soil physicochemical properties, but also positively correlated with bacterial and fungal diversity, because most exoenzymes and endoenzymes are mainly from microbes [27,29]. In particular, S-ACP was significantly negatively correlated with C/N ratio, positively correlated with N/P ratio, and positively correlated with bacterial and fungal Shannon index ( Figure 8). Consequently, microbes may consume more C and N to secrete enzymes related to phosphorus metabolism under the condition of phosphorus deficiency [9]. In conclusion, fungi may have more competitive advantages than bacteria on TOC and TN in P. massoniana plantations under low phosphorus conditions. Therefore, the recovery of TOC and TN leads to an increase in fungal diversity, thus improving the enzyme activity. In turn, increased fungal diversity and enzyme activity promoted nutrient cycling.

Conclusions
Our results showed that the nutrient availability of woodland decreased significantly with the increase of the chronosequence of P. massoniana plantations, especially AP. Meanwhile, soil C/N ratio decreased while C/P and N/P ratios increased, indicating that there was a more serious P limitation in the near-mature forest and mature forest stage. On this basis, coupled with the decrease of pH and NO 3 --N, these resulted in the decrease of soil bacterial alpha diversity. However, fungi can increase their diversity by consuming more TOC and TN under the condition of P limitation, which leads to the increase of enzyme activity (especially S-ACP). From different soil depths, soil physicochemical properties, microbial diversity, and enzyme activities decreased synchronously, and there was a close correlation between them. We also analyzed the microbial composition at different ages and soil depths and found that bacteria and fungi in the 10-year stand were significantly separated from other stands, and bacteria and fungi in the 0-20 cm soil depth were significantly separated from the other two soil layers. In all stands and soil depths, Acidobacteria was the most abundant bacterial community, followed by Proteobacteria and Actinobacteria; Basidiomycota was the most abundant fungal community, followed by Ascomycota. Furthermore, SWC, TOC, AP, and C/P ratio were the main factors affecting the microbial community; in particular, SWC had the highest contribution rate, and excessive SWC might inhibit the survival of most microbes. In short, these results reflect the potential relationship between soil physicochemical properties, microbial diversity and composition, and enzyme activities, which enhances understanding of soil fertility in different soil depths during the development of P. massoniana plantations and provides a reference for manual intervention.
Author Contributions: Q.G. planned and supervised the project, conceived and designed the experiments, reviewed drafts of the paper. J.P. performed the experiments, analyzed the data and led manuscript production. H.L., S.L., Y.Z., S.Y., X.S., and Y.Q. collected and analyzed the data. X.F. analyzed the data and made charts. All authors have read and agreed to the published version of the manuscript.