From Young to Over-Mature: Long-Term Cultivation Effects on the Soil Nutrient Cycling Dynamics and Microbial Community Characteristics Across Age Chronosequence of Schima superba Plantations

Abstract

Optimizing forest management requires a comprehensive understanding of how soil properties and microbial communities evolve across different plantation ages. This study examines variations in soil nutrient dynamics, enzyme activities, and bacterial communities in Schima superba Gardn. & Champ plantations of 10, 15, 27, 55, and 64 years. By analyzing soil from depths of 0–20 cm, 20–40 cm, and 40–60 cm, we identified significant age-related trends in soil characteristics. Notably, nutrient contents, including total organic carbon (TOC), total phosphorus (TP), total carbon (TC), total nitrogen (TN), and nitrate nitrogen (${\mathrm{NO}}_3^{-}$-N), as well as soil water content (SWC), peaked in 55-year-old mature plantations and decreased in 64-year-old over-mature plantations. Enzyme activities, such as urease, sucrase, and acid phosphatase, decreased with soil depth and exhibited notable differences across stand ages. Microbial community analysis indicated the predominance of Acidobacteria, Chloroflexi, Proteobacteria, Actinobacteria, and Verrucomicrobiota in nutrient cycling, with their relative abundances varying significantly with age and depth. Mature and over-mature plantations exhibited higher absolute abundances of functional genes related to methane metabolism, nitrogen, phosphorus, and sulfur cycling. Reduced calcium ion levels were also linked to lower gene abundance in carbon degradation, carbon fixation, nitrogen, and phosphorus cycling, while increased TOC, ${\mathrm{NH}}_4^{+}$-N, ${\mathrm{NO}}_3^{-}$-N, and AP correlated with higher gene abundance in methane metabolism and phosphorus cycling. Our findings suggest that long-term cultivation of Schima superba enhances soil nutrient cycling. Calcium ion was identified as a significant factor in assessing soil properties and microbial dynamics across different stand ages, suggesting that extended plantation rotations can improve soil health and nutrient cycling.

1. Introduction

Stand age is regarded as an important indicator of plantation development status and structure characteristics following afforestation [1]. Stands of varying ages exhibit distinct differences in understory vegetation, soil quality, and microbial communities [2]. Soil quality, which directly impacts timber yield and quality, is typically assessed through a combination of physical, chemical, and biological properties [3]. The feedback between plants and soils across different forest stand ages is well-documented [4]. During various stages of forest development, notable changes occur in vegetation productivity and canopy density, which in turn affect internal water and thermal conditions, as well as the composition and distribution of surface litter [5,6]. These factors substantially influence soil nutrient cycling, enzyme activities, and microbial communities by regulating understory vegetation composition and the deposition of litter and root exudates [7,8]. Understanding the physical and chemical properties of soil and microbial structures at different stand ages is crucial for accurately assessing nutrient requirements and developing effective and sustainable management practices for plantations.

Plant growth significantly alters both biotic and abiotic dynamics in the surrounding soils, altering soil structure and properties, which in turn affects plant suitability [9]. Stand age indirectly influences soil nutrient cycling, and extensive research has examined its effects on soil biological and abiotic properties across various crops and forest species [10,11]. For instance, in Chinese fir plantations, the carbon-to-phosphorus and nitrogen-to-phosphorus ratios peak in 24-year-old mature forests, while the carbon-to-nitrogen ratio remains stable, with significant variations observed primarily in the 20–60 cm soil layer [12]. Similarly, in Pinus tabulaeformis plantations, increasing forest age generally enhances soil moisture, water holding capacity, and porosity, while bulk density shows minimal variation. Organic matter, nitrogen, and phosphorus levels display an overall increasing trend [13]. Studies on soil carbon fractions in Cryptomeria japonica plantations reveal that these fractions initially increase with forest age before stabilizing, with soil organic carbon stability improving over time [14].

Soil microbiome is a vital indicator of plant–soil feedback due to its dynamic role in the soil ecosystem and its sensitivity to changes in the soil microenvironment [15]. Stand age significantly impacts microbial diversity and community structure. Zhang et al. examined the soil microbiome at four developmental stages of Chinese fir plantations and observed that, as stand age increased, the abundance of oligotrophic bacteria decreased while copiotrophic bacteria increased [16]. Simultaneously, fungal functional groups associated with plants declined, whereas saprotrophic fungi showed a significant increase. In Larix gmelinii plantations aged 10, 30, and 50 years, Li et al. observed that older forests promoted more specialized bacterial and fungal communities with stronger bacterial–fungal associations [17]. Similarly, Wang et al. reported that phosphatase activity and the abundance of the phoC gene were higher in 18- and 25-year-old Chinese fir plantations, while the abundance of the phoD gene gradually increased with stand age from 6 to 49 years [18].

Across a chronosequence of Chinese fir stands aged 7, 16, 29, 36, and over 80 years, Zhan et al. found that the stand developmental stage created a spatial gradient in gross ammonification and nitrification rates by altering the abundance of microbial functional genes, which in turn affected plantation productivity by modulating the supply of bio-available nitrogen [19]. Recent advances in high-throughput sequencing and functional genomics have greatly improved our ability to profile microbial communities and their functional potentials [20,21]. However, the dynamics of keystone taxa and their gene functional profiles across varying forest age stages remain poorly understood.

Soil enzymes, primarily secreted by microorganisms and plant roots, play a key role in multiple soil functions [22]. Soil microbial communities and plants significantly influence enzyme activity, which serves as an indicator of microbial function [23]. García-Saucedo et al. investigated soil enzyme activity in Spanish black pine forests across various stand ages, specifically 50, 70, 90, 110, 170, and 210 years [24]. The study found that dehydrogenase enzyme activity varied significantly with stand age, highlighting its potential as an indicator of forest ecosystem development. Catalase, urease, sucrase, and acid phosphatase are essential enzymes involved in soil carbon, nitrogen, and phosphorus cycles. Their activities are important indicators of soil nutrient dynamics, reflecting nutrient transformation and availability within the soil ecosystem [25]. Wang et al. reported that catalase and urease activities gradually increased from 3 to 38 years in Chinese fir plantations, suggesting that older stands may enhance soil properties, particularly enzyme activities and microbial community diversity [26]. Therefore, preserving functional diversity in the soil microbiome is essential for long-term forest productivity, as it supports soil enzyme activity and nutrient cycling.

As a whole, soil quality is a function of complex interactions between soil physicochemical properties, soil enzyme activities, and soil microbial activities that influence nutrient cycling dynamics. For instance, low soil pH reduces the availability of essential nutrients such as calcium, potassium, and ammonium and the sequestration of soil nutrients like phosphorus, rendering it insoluble by fixing with cations [27,28]. However, some microorganisms have been shown to solubilize nutrients—such as Bacillus, which solubilizes phosphate [29]. Soil microbial activities are associated with soil enzymes, which are important for the decomposition of labile organic matter, thereby activating biogeochemical cycling. The enzymatic activities reflect the functional diversity of the microorganisms involved in decomposition processes and adding a considerable amount of carbon into the soil system [30]. This carbon is decomposed by carbon degradation-related enzymes in the soil, which in turn provide carbon and energy for microbial growth and metabolism [7]. Furthermore, soil microbes produce extracellular enzymes that are responsible for the mineralization and cycling of terrestrial nitrogen, phosphate, and carbon [30].

Schima superba, an evergreen broad-leaved tree belonging to the Theaceae family, is a dominant species in China’s subtropical forests due to its adaptability, fire resistant, medium growth rate, and high yield [31]. This species naturally occurs south of 32°N and east of 96°E in China, encompassing 13 southern provinces, including Fujian, Guangdong, and Zhejiang [32]. It thrives in warm and humid climates, particularly in regions with an average annual temperature ranging from 16 °C to 22 °C and annual precipitation between 1100 mm and 2000 mm. Notably, Schima superba exhibits remarkable soil adaptability, growing on both yellow and red soils with pH values between 4.5 and 6.0. It also flourishes on slopes prone to severe soil erosion and arid conditions. In recent years, Schima superba has been increasingly adopted as a replacement for Pinus massoniana due to the widespread infestation of pine wilt disease in China [33]. This transition has led to a significant increase in the establishment of Schima superba plantations, particularly for cultivating large-diameter timber through extended rotation periods [34]. Despite the growing importance of Schima superba, a significant knowledge gap remains regarding soil properties and microbial communities associated with plantations of varying ages. This study aimed to investigate the response of soil’s physicochemical properties, enzyme activities, bacterial diversity, and functional gene composition in Schima superba plantations aged 10, 15, 27, 55, and 64 years in southern China. The goal was to provide scientific guidance for maintaining soil fertility and managing these plantations effectively. Specifically, the objectives of this study were to: (1) examine the dynamics of soil’s physicochemical properties; (2) detect shifts in soil’s enzyme activities; (3) assess variations in bacterial diversity and functional gene composition; and (4) analyze the relationships between stand age, soil properties, enzyme activities, and microbial activity. By elucidating these interactions, the study aimed to develop effective forest management practices that optimize soil health and enhance the productivity and sustainability of Schima superba plantations.

2. Materials and Methods

2.1. Study Site Description

The study was conducted in a small watershed within the Pusang State-Owned Forest Farm (26°30′ N –28°20′ N, 117°00′ E –119°25′ E), located in Fujian Province, Southeast China. This region is characterized by low mountain and hilly terrain, south of the Wuyi Mountains, and features a typical maritime subtropical monsoon climate with a mean annual temperature of 18.5 °C, annual precipitation of 1756 mm, 305 frost-free days, and a relative humidity of 82%. The soil, classified as red soil according to the genetic soil classification of China, is acidic and originates from deeply weathered deposits of Cretaceous-aged granites and conglomerates [35]. We selected five representative stands of Schima superba, spanning a natural age gradient of 10 years (young-aged), 15 years (young-aged), 27 years (middle-aged), 55 years (mature-aged), and 64 years (over-mature-aged). For each stand, three replicated sites with similar altitude, topography and soil texture were chosen. A standard plot measuring 400 m² (20 m × 20 m) was established at each stand (n = 15 plots) for soil sampling. Detailed information on altitude, aspects, density, tree height, DBH, and main understory vegetation is provided in

Table 1. Basic situation of studied Schima superba stands.

Stand Age	Altitude/m	Aspects/°	Density/Stems·hm−2	Tree Height/m	DBH/cm	Main Understory Vegetation
10	210–232	25–28	2307	5.93 ± 0.32 a	6.97 ± 0.37 a	Indocalamus tessellatus, Spatholobus suberectus
15	200–245	26–29	1966	7.98 ± 0.13 a	10.09 ± 0.23 b	I. tessellatus, Lophatherum gracile, Cibotium barometz
27	194–225	25–29	405	15.23 ± 0.91 b	21.73 ± 1.27 c	I. tessellatus, S. suberectus, C. barometz
55	180–220	25–26	345	17.58 ± 1 c	22.13 ± 1.24 c	L. gracile, Phoebe zhennan, Tetrastigma hemsleyanum
64	190–215	26–28	210	22.71 ± 0.97 d	32.06 ± 1.23 d	Aralia spinifolia Merr, T. hemsleyanum, Phoebe bournei

Different lowercase letters indicate significant differences among tree height and DBH within the same stand age (p < 0.05) (Mean ± SE).

.

2.2. Soil Sampling

We collected five 5 cm diameter soil cores using the “S-sampling” method from each plot, with cores taken to depths of 0–20 cm, 20–40 cm, and 40–60 cm using a stainless steel cylindrical driller. The cores from each depth were mixed to form composite samples, which were then passed through a 2 mm sieve to remove plant debris, roots, and rocks. Each composite sample was divided into three subsamples: One was air-dried and ground for chemical property analysis. Another was stored at 4 °C in a portable refrigerator for soil enzyme activity assays and measurements of ammonium nitrogen (${\mathrm{NH}}_4^{+}$-N) and nitrate nitrogen (${\mathrm{NO}}_3^{-}$-N). The final subsample was stored in a dry ice incubator, transported to the laboratory, and frozen at −80 °C for analysis of soil bacterial communities and functional gene abundance. A 100 cm³ sample from the middle of each soil layer was used to measure soil bulk density (BD) and maximum water holding capacity (SMC), and the SMC was determined by oven-drying the samples to a constant weight at 105 °C.

2.3. Soil Chemical Properties Analysis

Soil pH was measured using a 1:2.5 ratio of air-dried soil to distilled water. Total carbon (TC) and total nitrogen (TN) content were analyzed by combustion with an elemental analyzer (Vario Macro Cube, Elementar, Frankfurt, Germany), and their contents were indicated by the C/N ratio. The total organic carbon (TOC), ${\mathrm{NH}}_4^{+}$-N and ${\mathrm{NO}}_3^{-}$-N were extracted from 3 g of soil using 30 mL solution of a mixture of 2 mol·L⁻¹ KCl, respectively. The TOC concentration was determined by a total organic carbon analyzer (TOC-VCPH/CPN Analyzer, Shimadzu, Kyoto, Japan), the ${\mathrm{NH}}_4^{+}$-N and ${\mathrm{NO}}_3^{-}$-N concentrations were measured using a continuous flow analyzer (San++, Skala, Amsterdam, Netherlands). The total phosphorus (TP) and total potassium (TK) were determined using hydrofluoric acid–perchloric acid digestion followed by analysis with inductively coupled plasma emission spectrometry (PE OPTIMA 8000, PerkinElmer, Waltham, MA, USA). Calcium (Ca²⁺) was extracted using a mixture of HF, HClO₄ and HNO₃, and was subsequently analyzed with an inductively coupled plasma emission spectrometer (PE OPTIMA 8000, USA). The available phosphorus (AP) was measured using the HCl-NH₄F leaching method and quantified with a multi-function microplate reader (SpectraMax M4, MolecularDevices, San Francisco, CA, USA).

2.4. Soil Enzyme Activity Analysis

Urease (UE) activity was measured colorimetrically using sodium phenol and sodium hypochlorite. Sucrase (SC) activity was also determined colorimetrically, employing 3,5-dinitrosalicylic acid [36]. Soil catalase (CAT) activity was assessed by mixing 0.1 g of soil with 1 mL of phosphate-buffer solution to create an ice-cold homogenate and then centrifuging at 4 °C for 10 min at 8000 rpm to obtain the supernatant. Reagents were added according to the instructions provided by the kit (Beijing Solai News Technology Co., Ltd., Beijing, China), and CAT activity was measured using an ultraviolet spectrophotometer. Acid phosphatase (ACP) activity was determined by adding 0.05 mL of toluene to 0.1 g of air-dried soil, incubating the mixture at 37 °C for 24 h, and then centrifuging at 25 °C for 10 min at 8000 rpm. Reagents were added following the kit instructions (Beijing Solai News Technology Co., Ltd., Beijing, China), and ACP activity was measured with a visible spectrophotometer.

2.5. Soil Bacterial Community Analysis

Bacterial communities were analyzed and sequenced using the Illumina MiSeq PE300 platform (Beijing Ovison Bioinformatics Co., Ltd., Guangzhou, China). The V4 hypervariable region of the bacterial 16S rRNA gene was amplified by polymerase chain reaction (PCR) using the primers GTGCCAGCMGCCGCGGTAA and GGACTACHVGGGTWTCTAAT. For each sample, three replicates of PCR products were combined for sequencing. Paired-end sequencing was performed on the Illumina MiSeq platform. The resulting data were processed using QIIME 1.8.0 to filter, splice, and remove chimeric sequences [37]. Operational taxonomic units (OTUs) were identified and compared with the Silva database to obtain taxonomic information for each OTU. Alpha diversity was assessed using Mothur 1.40.5 [38]. Beta diversity analysis was based on the Weighted UniFrac distance, with cluster analysis carried out using the pheatmap package in R 4.3.3. The UniFrac algorithm was used to assess phylogenetic differences between samples. Soil microbial community OTUs were converted to relative abundances, and the microorganisms at the phylum level were analyzed using one-way ANOVA in IBM SPSS Statistics 19 to determine statistical significance and errors, with final plots generated using Origin 8.5.

2.6. Analysis of Microbial Functional Genes Associated with C, N, P, S Cycling

The total microbial DNA was extracted and purified from soil samples collected at a depth of 0–20 cm. To quantify functional genes associated with biogeochemical cycling of carbon (C), nitrogen (N), phosphorus (P), and sulfur (S), we utilized a high-throughput quantitative PCR-based chip provided by Guangdong MAGI Gene Technology Co., Ltd. (Guangzhou, China). This technique, known as Quantitative Microbial Element Cycling (QMEC), allows for the detection and quantification of 71 functional genes associated with these elemental cycles, specifically targeting 36 carbon cycle genes, 22 nitrogen cycle genes, 8 phosphorus cycle genes, and 5 sulfur cycle genes [39] (Table S1).

Purified DNA samples and qPCR reagents were loaded onto a 384-well plate. High-throughput automatic micro-sampling equipment was used to introduce the samples and primers into the nanowells of the high-throughput qPCR chip. The SmartChip Real-Time PCR system was then employed to conduct the qPCR reactions and detect fluorescence signals, which resulted in amplification and dissolution curves. Cycle threshold (Ct) values for each gene were determined using Canco software. Data were normalized against the 16S rRNA gene, which served as an internal control, to obtain relative quantification for each functional gene. Absolute quantification of the 16S rRNA gene was performed using a Roche instrument, and absolute quantities for the other genes were derived accordingly. A summary table of the Ct values was compiled and subjected to quality control. To calculate the relative copy number for each gene, the following formula was used: Relative copy number = 10^{(31−Ct)/(10/3)}. Absolute quantification of the 16S rRNA gene was obtained with the Roche instrument, and absolute quantification for other genes was derived by conversion [40]. Details on primer synthesis and validation are reported by Zheng et al. [39].

2.7. Statistic Analyses

Physicochemical properties and enzyme activity were analyzed using one-way analysis and two-way analysis of variance (ANOVA) and the least significant difference (LSD) with a significance level set at p < 0.05. Prior to data analysis, the normal distribution and homoscedasticity assumption for ANOVA were checked using Shapiro–Wilk and Levene’s tests. Both tests turned out non-significant (p > 0.05); hence, the normal distribution and homogeneity of variance assumptions were satisfied. The Levene’s tests were performed using SPSS Statistical software, Version 23. The Shapiro–Wilk test was performed using R Package. Microbial diversity metrics, including Chao1 and Shannon’s diversity index, were calculated using QIIME version 1.8.0 (accessed on 30 December 2013). The impact of soil physicochemical properties on the absolute abundance of carbon, nitrogen, and phosphorus functional microbial taxa was visualized using the Tutools platform, with further correlations analyzed using OriginPro 2018. Redundancy analysis (RDA) of dominant bacterial phyla and genera with environmental factors and enzyme activity was conducted using Canoco 5. Additionally, Partial Least Squares Path Model (PLS-PM) was established using R Package to explore the relationships between soil physicochemical properties, enzyme activity, bacterial diversity indices, and microbial functional genes.

3. Results

3.1. Soil Physicochemical Properties

Bulk density in Schima superba plantation soils increased with depth, with the BD value in the 40–60 cm layer being higher than that in the 0–20 cm layer (Figure 1a). In the 10a stands, SWC increased with depth, while in the 27a and 55a stands, it decreased with depth, potentially reflecting differences in soil moisture retention or root water uptake dynamics across stand ages (Figure 1b). Soil pH decreased with depth in the 10a stands but increased with depth in stands of other ages (Figure 1c). There is an interactive effect of pH between stand age and soil layer (p < 0.01). Nutrient content was generally higher in the 0–20 cm soil layer, with concentrations of TOC, TP, TN, TC, ${\mathrm{NH}}_4^{+}$-N, AP and ${\mathrm{NO}}_3^{-}$-N all decreasing with increasing soil depth (Figure 1d–k), indicating a clear accumulation of nutrients at the surface. And there is an interactive effect of ${\mathrm{NH}}_4^{+}$-N and ${\mathrm{NO}}_3^{-}$-N between stand age and soil layer (p < 0.05).

Bulk density decreased progressively in the 10, 15, and 27a stands, while it increased in the 55 and 64a mature stands, suggesting a potential stabilization or compaction of the soil structure in older stands (Figure 1a). In the 0–20 cm soil layer, SWC, TOC, TP, TN, TC, AP, and ${\mathrm{NO}}_3^{-}$-N increased with stand age, reaching their highest values in the 55a stands (Figure 1b,d–g,i,k). Soil pH was highest in the 10a stands and decreased sharply to its lowest level in the 15a stands. It then showed a gradual increase in the 27 and 55a stands (Figure 1c). The TOC, TP, and ${\mathrm{NO}}_3^{-}$-N generally increased with stand age. Notably, these parameters significantly peaked in the 55a stands but decreased in the 64a stands (Figure 1d,e,k). Both ${\mathrm{NH}}_4^{+}$-N and TK exhibited a “V” shaped trend with increasing stand age and the lowest point in the 27a stands (Figure 1h,j). Specifically, TK was lowest in the 15a stands within the 0–20 cm soil layer. The concentration of calcium ions in the 0–20 cm surface soil layer was highest in stands aged 15 years and decreased with increasing stand age. Notably, the calcium ion concentration in 64a stands was significantly lower than in all other age classes in the 0–20 cm soil layers. Furthermore, in the 20–40 cm and 40–60 cm soil layers, the calcium ion concentrations were significantly higher in the 10 and 15a stands compared to the 27, 55, and 64a middle-aged and mature stands (Figure 1l). These findings underscored that soil nutrient distribution and concentration in Schima superba stands are significantly influenced by stand age. Notably, the 55-year-old mature stands generally exhibited superior soil health and nutrient status compared to younger stands.

3.2. Soil Enzyme Activities

The activity of CAT was not significantly influenced by soil depth. However, the activities of soil SC, UE, and AP were significantly affected by soil depth, with all these activities decreasing as soil depth increases (Figure 2). Significant differences in enzyme activity were also observed across different stand ages. AP activity showed a “V”-shaped trend with increasing stand age, initially decreasing and then increasing, reaching its minimum in the 27a stands (Figure 2a). CAT activity gradually increased with stand age, with the 55a and 64a mature stands exhibiting significantly higher activity compared to other age groups, while the 10a stands showed notably lower activity (Figure 2b). UE content was significantly higher in the 55a stands compared to the other four age groups at the same soil depth, whereas the 15a stands exhibited the lowest UE activity (Figure 2c). And there is an interactive effect of UE between stand age and soil layer (p < 0.05). The activity of the SC enzyme was the lowest in the 10a stands, with significantly higher levels observed in the 27a, 55a, and 64a stands compared to the 10a and 15a stands.

3.3. Bacterial Communities

Through the annotation of 26,472 bacterial OTUs, we identified classification information for 55 bacterial phyla and 855 bacterial genera. The number of unique OTUs in 10a, 15a, 27a, 55a, and 64a was 2040, 2002, 1838, 1972, and 1957, respectively, and 6142 bacterial OTUs were present in all stand ages (Figure S1). The results reveal significant differences in the relative abundance of dominant bacterial communities across plantations of varying ages. The principal bacterial phyla observed in different stand ages included Acidobacteria (39.08%–53.18%), Proteobacteria (12.29%–24.61%), Actinobacteria (2.67%–7.24%), Chloroflexi (6.75%–22.08%), and Verrucomicrobiota (3.51%–5.82%) (Figure 3). HSB_OF53-F07, Candidatus_Solibacter, Acidothermus, uncultured_Acidobacteria_bacterium, and Acidibacter were the dominant bacterial genera (Figure S2). Additionally, phyla such as RCP2-54, Planctomycetota, Myxococcota, WPS-2, and Crenarchaeota exhibited dominance across different aged stands, collectively accounting for over 90% of the total bacterial community (Figure 3). Moreover, variations in soil depth significantly influenced the relative abundance of these dominant bacterial phyla (p < 0.05). Specifically, Acidobacteria, Chloroflexi, and Actinobacteria showed significantly higher relative abundance in surface soils compared to deeper layers. Conversely, the abundance of Proteobacteria increased with soil depth. There was no significant difference in bacterial diversity among different soil layers but among stand ages (Table S2). In the 0–20 cm surface soil, the Chao1, Observed species, and Shannon index values were higher in mature (55a) and over mature (64a) stands compared to young stands (10a, 15a) and middle-aged stands; in particular, the diversity of 15a stands was the lowest (Table S2).

3.4. Links Between Soil Physicochemical Properties, Enzyme Activities, and Bacterial Communities

The redundancy analysis (RDA) indicated that soil physicochemical properties and enzyme activities account for a substantial proportion of the variation in soil bacterial communities at various depths: 87.50% at 0–20 cm, 62.49% at 20–40 cm, and 50.42% at 40–60 cm (Figure 4a–c). The RDA analyses and hot map linking soil bacterial communities to soil physicochemical properties and to soil enzyme activities are presented in Figure 5. Examination of specific soil layers elucidates the key factors influencing bacterial community composition and diversity. In the 0–20 cm soil layer, the bacterial community composition was predominantly influenced by AP, TK, pH, ${\mathrm{NH}}_4^{+}$-N and Catalase (Figure 4a). AP showed a significant positive correlation with Acidobacteria (p < 0.05) and a significant negative correlation with Proteobacteria (p < 0.01) (Figure 5a). TK was significantly negatively correlated with Crenarchaeota (p < 0.05). Soil pH had a significant positive correlation with Myxococcota (p < 0.01) and Verrucomicrobiota (p < 0.05), and it had a significant negative correlation with Crenarchaeota (p < 0.01). ${\mathrm{NH}}_4^{+}$-N and ${\mathrm{NO}}_3^{-}$-N exhibited significant negative correlations with Chloroflexi (p < 0.01). Additionally, Catalase showed a significant positive correlation with Actinobacteria (p < 0.05) and a negative correlation with Chloroflexi (p < 0.05). TOC showed a significant positive correlation with Acidobacteria (p < 0.05) and a significant negative correlation with Chloroflexi (p < 0.05) (Figure 5a). Ca²⁺ demonstrated significant positive correlations with WPS-2 (p < 0.05).

In the 20–40 cm soil layer, the bacterial community composition was predominantly influenced by pH, TC, ${\mathrm{NH}}_4^{+}$-N and Catalase (Figure 4b). Soil pH exhibited a strong positive correlation with Myxococcota (p < 0.01) and significant positive correlations with Verrucomicrobiota and Acidobacteria (p < 0.05) (Figure 5b). Conversely, soil pH showed strong negative correlations with Actinobacteria and WPS-2 (p < 0.01) and a significant negative correlation with Crenarchaeota (p < 0.05). TC was significantly positively correlated with RCP2-54 (p < 0.05). ${\mathrm{NH}}_4^{+}$-N had a significant positive correlation with Verrucomicrobiota (p < 0.05). Ca²⁺ exhibited significant negative correlation with Catalase (p < 0.05), while Catalase showed a significant positive correlation with Crenarchaeota (p < 0.05). And AP showed a significant positive correlation with Planctomycetota (p < 0.05). In the 40–60 cm soil layer, the bacterial community composition was predominantly influenced by soil pH and Urease (Figure 4c). Soil pH had highly significant positive correlations with Planctomycetota and Myxococcota (p < 0.01) and a significant negative correlation with WPS-2 (p < 0.05) (Figure 5c). Urease activity demonstrated a significant negative correlation with Chloroflexi (p < 0.01) and a positive correlation with Myxococcota (p < 0.01). And TOC and TP showed a significant negative correlation with RCP2-54 (p < 0.05). Ca²⁺ showed a significant negative correlation with Verrucomicrobiota (p < 0.05) and a positive correlation with Chloroflexi (p < 0.05).

3.5. Absolute Abundance of Soil Bacterial Functional Genes

Analysis of the absolute abundance of soil bacterial functional genes across different stand ages revealed several significant differences. Functional genes with absolute abundances greater than 10³ copies/ng and their abundances distributions across various stand ages, are illustrated in Figure S5. Specifically, the absolute abundance of 16S rDNA was significantly higher in 15a and 27a stands compared to other stand ages (Figure 6a). Carbon degradation functional genes were also notably more abundant in 64a and 27a stands, with the abfA gene showing the highest abundance. Additionally, the cex gene was only detected in 54-year-old and 64a stands (Figure 6b). Regarding carbon fixation, no significant differences in gene abundance were observed across different stand ages, although higher levels were noted in 64a and 15a stands. The genes with the highest abundance in this category were acsA (acetyl-CoA synthetase) and mct (malate: citrate translocase) (Figure 6c). In methane metabolism, functional gene abundance was significantly higher in 64a stands compared to other ages, with pmoA (particulate methane monooxygenase) being the most abundant gene (Figure 6d).

The abundance of functional genes involved in nitrogen, phosphorus, and sulfur cycling was significantly greater in 64a stands compared to other ages. Among nitrogen cycling genes, ureC (urease) exhibited the highest absolute abundance (Figure 6e). In phosphorus cycling, phnK (phosphonate transport system protein) was the most abundant gene (Figure 6f). In the sulfur cycle, apsA (adenosine-5′-phosphosulfate reductase), soxY (sulfur oxidizing protein), and yedZ (sulfite oxidase) were the three most abundant genes, with their levels increasing with stand age (Figure 6g). Finally, while the total abundance of all functional genes did not show significant differences across stand ages (Figure 6h), the overall abundance was higher in 64a stands compared to other ages.

3.6. Correlations Analysis of Soil Physicochemical Properties and Bacterial Functional Genes Associated with C, N, P, S Cycling

An analysis of the correlations between soil physicochemical properties and bacterial functional genes associated with C, N, P, S cycling across varying stand ages of Schima superba revealed significant patterns. The absolute abundance of functional genes related to carbon degradation, carbon fixation, nitrogen cycling, and phosphorus cycling was negatively correlated with calcium ion content (p < 0.05) (Figure 7). Specifically, the abundance of genes associated with methane metabolism and sulfur cycling showed a highly significant negative correlation with calcium ion content (p < 0.01). Additionally, the absolute abundance of methane metabolism genes was positively correlated with TC and TOC (p < 0.05) and exhibited a highly significant positive correlation with ${\mathrm{NO}}_3^{-}$-N, ${\mathrm{NH}}_4^{+}$-N, and AP (p < 0.01). Furthermore, the abundance of functional genes involved in phosphorus cycling was positively correlated with TOC (p < 0.05) and showed highly significant positive correlations with ${\mathrm{NO}}_3^{-}$-N, ${\mathrm{NH}}_4^{+}$-N, and AP (p < 0.01). These results suggested that lower levels of calcium ion content are associated with reduced abundance of genes involved in carbon degradation, carbon fixation, nitrogen cycling, and phosphorus cycling. Conversely, higher levels of TC, TOC, ${\mathrm{NO}}_3^{-}$-N, ${\mathrm{NH}}_4^{+}$-N, and AP are positively associated with the abundance of genes involved in methane metabolism and phosphorus cycling.

Partial Least Squares Path Modeling (PLS-PM) was employed to analyze the correlations among stand ages, soil physicochemical properties, enzyme activity, bacterial communities, bacterial diversity, and the abundance of functional genes associated with C, N, P, S cycling, as illustrated in Figure 8. The model demonstrated a goodness-of-fit (GOF) value ranging from 0.622 to 0.657, with path coefficients below 0.3 being excluded. The analysis revealed that stand age significantly and negatively impacted calcium ion content (p < 0.05). Furthermore, calcium ion content exerted a significantly negative direct effect on the abundance of functional genes related to methane metabolism (path coefficient = −0.425, p < 0.05) and sulfur cycling (path coefficient = −0.516, p < 0.05). Conversely, stand age had a substantial and positive direct effect on soil physicochemical properties (p < 0.001). Soil physicochemical properties, in turn, demonstrated a significant positive direct effect on methane metabolism (path coefficient = 0.943, p < 0.05) and phosphorus cycling (path coefficient = 0.934, p < 0.05). Stand age showed no significant direct effect on enzyme activity, bacterial communities, bacterial diversity, or the abundance of functional genes associated with C, N, P, and S cycling. As forest age increased, the data indicated that stand age significantly and directly influenced soil physiological properties, including SWC and nutrient levels such as TOC, TP, ${\mathrm{NO}}_3^{-}$-N, ${\mathrm{NH}}_4^{+}$-N, and AP—all of which exhibited a gradual increase with stand age, peaking at stands aged 55 years. Calcium ion content was identified as a significant factor in the variation among stand ages. Improved soil quality enhanced the activities of acid phosphatase, peroxidase, urease, and sucrase, which, in turn, promoted increased soil microbial activity, methane metabolism, and phosphorus and sulfur cycling.

4. Discussion

4.1. Soil Physicochemical Properties Across Different Stand Ages of Schima superba

The development of forest stands impacts soil physicochemical properties through changes in nutrient uptake and litter return. Younger stands, with lower canopy cover, typically have higher soil temperatures and more fluctuating moisture, affecting bulk density and porosity [41]. Conversely, mature stands with dense canopies and deeper root systems generally favor soil organic matter accumulation and improve soil structure [42]. In this study, the TOC, TP, TN, TC, ${\mathrm{NH}}_4^{+}$-N, AP, and ${\mathrm{NO}}_3^{-}$-N showed a clear decline with soil depth, reflecting surface nutrient accumulation across all stand ages of Schima superba plantations, in line with findings for Chinese fir and Douglas fir. Bulk density increased with soil depth but exhibited a V-shaped trend with increasing stand age (Figure 1a) [43,44], primarily due to stand density and artificial tending [45]. In the 0–40 cm soil layer, SWC generally increased with stand age, peaked in the 55a stands, but decreased significantly in the 64-year stands. This trend contrasts with findings for Eucalyptus, where SWC decreases with stand age [46], suggesting that Schima superba plantations exhibit better water retention with stand development and that water limitations to stand development are minimal. Additionally, while Wang et al. reported an increase in soil pH with stand age in 3, 12, and 38-year-old Chinese fir plantations, our study found the highest soil pH in 10a stands and the pH value drops sharply in young plantations [47]. This discrepancy may be attributed to the fact that broad-leaved Schima superba plantations are less prone to soil acidification compared to fast-growing coniferous stands such as Chinese fir. And the rapid growth of younger plantations may limit the growth of microorganisms due to nutrient competition or secretion.

In this study, TOC, TP, and ${\mathrm{NO}}_3^{-}$-N generally increase with stand age, with TN, TC, and AP peaking in the 55-year-old stands in the 0–20 cm soil layer. This was likely attributed to the rapid growth phase of Schima superba between 10a and 27a, during which the trees absorb substantial amounts of soil nutrients. Schima superba experiences peak growth and flourishing between 5 and 45 years, with DBH growth accelerating from 15 years and peaking at 65 years [48]. However, by 55a, plantation density decreases, growth slows, and nutrient absorption diminishes. Decomposing litter and plant residues contribute to soil enrichment. This is consistent with the finding of Wang et al., which gradually increased with the age of the forest due to litter content and plant nutrient uptake, TP, TK, AN and AK content. [26,47]. In 64-year-old Schima superba plantations, the growth in height and diameter slows, litter return diminishes, pest resistance weakens, and topsoil fertility begins to decline. He et al. observed similar trends in Pinus massoniana plantations, where nitrogen, phosphorus, and potassium levels initially increased and then decreased with age [49]. The decrease in Ca²⁺ in the 64a stands within the 0–20 cm soil layer may signal a critical shift in soil nutrient cycling, where over-mature plantations may face nutrient depletion due to reduced litter input and altered leaching processes. ${\mathrm{NH}}_4^{+}$-N and TK were lowest at 27a, a period of rapid growth for Schima superba, highlighting the need for timely ammonium nitrogen and potassium fertilization. Additionally, phosphorus deficiency was observed throughout all growth stages of the plantation [50], emphasizing the importance of phosphate fertilizer application to address this limitation.

4.2. Soil Enzyme Activities Across Different Stand Ages of Schima superba

Soil enzyme production and distribution are mainly influenced by microorganisms and root systems, as they play a key role in nutrient acquisition [51]. In this study, the activities of SC, UE, and AP decreased with soil depth, which is consistent with previous findings [26,30,52]. This trend can be attributed to the greater abundance of soil microorganisms in the upper soil layers, enriched by understory vegetation and organic matter, which supports greater enzyme activity [53]. Elevated CAT activity in mature and over-mature stands suggests enhanced microbial metabolism and better soil ecosystem function, aligning with soil physicochemical properties and microbial community structure [54]. Enzyme activities of UE, SC, and AP are associated with nitrogen, carbon, and phosphorus cycling, and generally increase with stand age due to higher organic matter inputs and increased microbial biomass [55,56]. However, in over-mature stands, these activities may decline due to reduced litter input, decreased microbial activity, or changes in organic matter quality. Specifically, AP activity was lowest at 27 years, while SC and UE activities were minimal at 10 and 15 years, respectively. Older stands showed higher enzyme activities, indicating more efficient nutrient cycling. These trends were linked to soil properties such as TOC, TP, TC, TN, AP, ${\mathrm{NH}}_4^{+}$-N, and ${\mathrm{NO}}_3^{-}$-N, which were significantly correlated with enzyme activities.

4.3. Soil Bacterial Communities Across Different Stand Ages of Schima superba

Soil bacterial diversity and community composition typically change as plantations mature [57]. In this study, mature and over-mature stands exhibited higher bacterial diversity compared to younger stands (10a and 15a), indicating that bacterial communities become more complex and stable with soil maturation. This observation is corroborated by the fact that mature plantations of Schima superba generally demonstrate better soil fertility and organic matter content, which supports a more diverse bacterial community. Dominant phyla identified include Acidobacteria, Proteobacteria, Actinobacteria, Chloroflexi, and Verrucomicrobiota, with similar compositions across different plantation ages, suggesting broad ecological niches and adaptability, essential for nutrient cycling and ecosystem functioning [16,17,58].

Acidobacteria, recognized as one of the most diverse and widespread phyla in soil ecosystems, exhibit high adaptability and often dominate acidic soils [26]. In Schima superba plantations, Acidobacteria demonstrated the highest relative abundance and had a relatively stable amount, likely due to the favorable acidic conditions that support their survival and that may be related to the wide ecological niche of their and its relatively low impact on environmental factors. When the competition between Acidobacteria and other microorganisms reaches an equilibrium state, their amount will also tend to stabilize. Furthermore, their abundance positively correlates with soil pH in the 20–40 cm soil layer. Acidobacteria transformed nutrients into plant-usable forms that are significantly involved in the breakdown of cellulose, lignin, organophosphate compounds, and other organic matter [59]. These functions may influence soil oxidative stress responses and nutrient cycling.

Acidobacteria exhibited the lowest abundance in the 15a stands and the highest abundance in the 55a mature stands, highlighting age-related variations in their functional contributions to nutrient cycling. Proteobacteria, essential for decomposing organic matter and facilitating carbon release, were strongly associated with soil carbon and nitrogen cycling [60]. Their abundance was notably correlated with ${\mathrm{NH}}_4^{+}$-N levels, indicating their involvement in nitrogen transformation [61]. Actinobacteria abundance was influenced by soil pH, with a decrease in diversity and abundance observed in more acidic soils [62]. The lowest levels were found in the 15a stands, which had the most acidic pH. These bacteria play crucial roles in plant residue decomposition and nitrogen fixation, linking their abundance to soil pH conditions. Chloroflexi, involved in methanogenesis and nitrogen cycling [63], displayed a significant positive correlation with calcium ions and a negative correlation with ${\mathrm{NH}}_4^{+}$-N and UE activity. PLS-PM analysis further indicated a negative impact of calcium ions on methane metabolism-related genes. These correlations between bacterial communities, soil nutrients, and enzyme activities illustrate complex interactions that influence soil health and fertility. Overall, these findings underscore the dynamic relationships between bacterial communities, soil nutrients, and enzyme activities across different plantation ages, enhancing our understanding of soil ecosystem functioning in Schima superba plantations.

4.4. Functional Genes Abundance and Correlations Across Different Stand Ages of Schima superba

Soil microorganisms are crucial drivers of the carbon, nitrogen, phosphorus, and sulfur cycles within ecosystems; however, the study of functional genes related to nutrient cycling in broad-leaved forests of different ages is still limited [39]. Our study revealed genes related to methane metabolism, nitrogen cycling, phosphorus cycling, and sulfur cycling showed an overall increase with plantation age. The increase in plantation age has fostered the development of understory vegetation, resulting in a rise in both litter quantity and decomposition rates, which in turn have contributed a significant amount of soil carbon [30]. The higher abundance of these functional genes in older plantations suggests that these ecosystems develop enhanced nutrient cycling capabilities over time.

Additionally, starch degradation genes (amyA) were most abundant in older Schima superba plantations (55a and 64a), highlighting a stronger starch degradation potential in these mature stages. Zhao et al. observed that hemicellulases had a strong capacity for degrading recalcitrant carbon in broadleaf forests [64]. Our study corroborates this, with the highest absolute abundances of genes involved in hemicellulose degradation (abfA, manB, xylA) found in Schima superba plantations. Interestingly, these genes were more abundant in younger plantations (10a, 15a, 27a) compared to older ones, suggesting that younger and mid-aged forests have a greater potential for degrading recalcitrant carbon. Methane metabolism genes were most prevalent in the 55a and 64a plantations, indicating higher methane production and oxidation rates in these soils.

The ureC gene, a key marker for urease activity, showed lower abundance and reduced urease (UE) activity in the 15a and 27a stands, which is consistent with previous findings [39]. With the development of Schima superba plantations, there was a concurrent rise in TN, ${\mathrm{NH}}_4^{+}$-N, and ${\mathrm{NO}}_3^{-}$-N contents, alongside an increase in the abundance of genes involved in nitrogen cycling, indicating enhanced nitrogen metabolism. Gao et al. previously found that an increased abundance of soil phosphorus cycling genes, such as phnK, phoD and ppX, accelerates the mineralization of organic phosphorus and polyphosphates into inorganic phosphorus [65]. However, with increasing stand age, phosphorus content accumulates, and acid phosphatase activity shows a “V” shape in Schima superba plantations. This study observed that the abundance of phosphorus-cycling genes, specifically phnK (involved in organic phosphorus mineralization) and ppX (involved in polyphosphate hydrolysis), increases with stand age. This supports previous studies showing that aging stands promote microbial growth, enhance substrate availability, and facilitate nutrient cycling. In Schima superba plantations, the abundance of sulfur cycling genes such as apsA, soxY, and yedZ increases with stand age, suggesting a strengthening of sulfur cycling processes as the stand matures. Notably, the abundance of sulfur-oxidizing genes, particularly yedZ and soxY, greatly exceeds that of the sulfur-reducing gene dsrA, indicating a stronger sulfur oxidation process compared to reduction [66]. Sulfur-oxidizing bacteria, which utilize sulfides or elemental sulfur as energy sources and have the capacity to oxidize both inorganic and organic sulfur compounds, are likely more prevalent in the 64a and 55a stands, thus enhancing sulfur removal from the soil.

Calcium is vital for microbial growth, playing a crucial role in maintaining cell wall integrity, enzyme activity, and intracellular signaling in bacteria [67]. It also acts as an important exchangeable base cation in soil, influencing soil quality, nutrient retention, buffering capacity, and pH regulation [68,69,70,71]. Soil acidification can lead to the leaching of calcium and other cations, such as magnesium, while increasing exchangeable aluminum, which can be toxic to plants and alter soil microbial communities [72]. In this study, we found a negative correlation between calcium ion content and soil pH values, with significantly decreased calcium concentrations observed from 15a to 27a stands. Reduced calcium ion levels were associated with lower gene abundance related to carbon degradation, carbon fixation, nitrogen, and phosphorus cycling. Wang et al. investigated the effects of stand age on calcium characteristics in different-aged Chinese fir plantations (3, 9, 17, and 26 years) [73]. Their findings indicated that soil Ca²⁺ content at a depth of 0–30 cm was prone to leaching and uptake by trees, which can lead to soil acidification as the plantation matures. Comparatively, as the age of Schima superba plantations increases, calcium ion concentration decreases while pH values rise from young (15a) to mature (55a), which helps mitigate soil acidification and enhances both soil bacterial diversity and the abundance of functional genes. Furthermore, studies have also indicated that certain bacterial phyla, such as Acidobacteria, Proteobacteria, and Actinobacteria, may respond differently to variations in soil calcium levels. While current findings highlight the significance of calcium, there was a call for more detailed research to elucidate the specific mechanisms by which calcium affects bacterial communities and to explore interactions with other soil properties.

4.5. Intearctions Between Soil Properties and Bacterial Community

The partial least squares pathway model revealed that soil physiological properties, including SWC and nutrient levels, such as TOC, TP, ${\mathrm{NO}}_3^{-}$-N, ${\mathrm{NH}}_4^{+}$-N, and AP, increased with stand age. Soil physicochemical properties, in turn, demonstrated a significant positive direct effect on functional genes related with methane metabolism and phosphorus cycling. The quality and quantity of organic matter accumulated on the soil surface increases with stand age. Furthermore, the microenvironment (soil moisture content, temperature) in the understory of the forest was favorable for the proliferation of soil bacteria. Depending on the type of litter accumulated on the surface, various enzymes are secreted to decompose and mineralize carbon, nitrogen, phosphorus, and sulfur [55]. This increased enzymatic activity favors the increase in microbial biomass and their functional diversity [55]. The decomposition of organic matter by soil microbial activities through increased secretion of enzymes improves the nutrient cycling, thus the soil quality was determined by the complex interactions between soil physicochemical properties, soil enzyme activities, and soil microbial community composition.

5. Conclusions

This study provides compelling evidence of the intricate interplay between stand age, soil nutrient dynamics, and microbial community structure in Schima superba plantations. Our findings reveal that the enhancement of soil nutrient cycling and microbial activity was significantly influenced by stand age, with mature plantations (55 years) demonstrating optimal conditions for nutrient cycling efficiency and microbial diversity. This underscores the crucial role that long-term cultivation plays in sustaining soil health and ecological balance in forest ecosystems. The observed increase in SWC, TOC, TP, and ${\mathrm{NO}}_3^{-}$-N with stand age highlights the progressive enrichment of soil fertility over time, supporting the hypothesis that older forests foster more robust soil ecosystems. Conversely, younger plantations (10 and 15 years) exhibited lower nutrient cycling efficiencies, indicating a pressing need for targeted soil management practices to support early-stage forest development. Our detailed microbial community analysis further illuminates the pivotal roles of Acidobacteria, Chloroflexi, Proteobacteria, Actinobacteria, and Verrucomicrobiota in nutrient cycling processes. These insights offer a novel perspective on how specific soil properties influence microbial communities and their ecological roles. The observed decline in nutrient dynamics and bacterial diversity in over-mature plantations (64 years) prompts a re-evaluation of forest management practices. These results suggest that extending rotation periods beyond the maturity stage may lead to diminishing returns in nutrient cycling efficiency and microbial diversity, necessitating adaptive management strategies to optimize forest health.

In summary, our study not only advances our understanding of the temporal dynamics of soil nutrient cycling in Schima superba plantations but also highlights the need for nuanced forest management approaches that account for age-related changes. We advocate for the integration of these findings into forest management practices to enhance soil health, sustain biodiversity, and improve the overall productivity of forest ecosystems. This research lays the groundwork for future studies aimed at optimizing forest management strategies and advancing our knowledge of forest–soil–microbe interactions.