Huangshan (HS) is a member of World Biosphere Reserve, Cultural Double Heritage, and Priority Areas for Biodiversity Conservation. This study was conducted in the Forest Dynamics Plot (10.24 ha) at HS, a subtropical evergreen broadleaved forest, Anhui Province, Southeast China, with the geographical coordinates 30°8′26″N, 118°6′38″E and an elevation of 400–600 m (Figure 1; Ding et al., 2016). The annual average temperature is 7.8°C, the average temperature in the hottest month (July) is 17.7°C, and the average temperature in the coldest month (January) is 3.1°C. The average annual rainfall and duration of the rainy season are 2,394.5 mm and 182 days, respectively (Ding et al., 2016). The forest dynamics plot was established in the summer of 2014 in western Huangshan, and a re-census was carried out in 2019. The forest dynamics plot is one of the important biodiversity observation platforms for subtropical evergreen broad-leaved forests in China. The leaf samples for this study were collected in August 2020.

We firstly selected Castanopsis eyrei as the focal species. It is a large, evergreen, subtropical tree that grows in broadleaf forests in southeastern China and is also a dominant species at HS (Xu et al., 2019; Staab et al., 2021). In addition, C. eyrei is also the species with the highest importance value in HS forestry dynamics plot (Ding et al., 2016).

Here, we measured SGDCs under different community structures based on species diversity indices in 256 plots. The 10.24 ha forest dynamics plot was divided into 256 subplots that were 20 × 20 m in size. The community assembly differences among subplots were assessed using the following methods. We first calculated three species diversity indices (the Shannon-Wiener index, the Gini-Simpson index, and the Pielou index) based on the 2019 re-census data. Based on the subplot species diversity indices, we then used a principal components analysis (PCA) to divide the 256 subplots into four groups (Supplementary Figure 1). Group I had a negative Pielou index value, Group II had positive Simpson index and Shannon-Wiener index values, Group III had a positive Pielou index value, and Group IV had negative Simpson index and Shannon-Wiener index values (PC₁ + PC₂ = 97.78%). Secondly, we considered the numbers of focal species (>5) in each plot and the distance (>20 m) among the subplots. After considering these factors, we randomly selected 24 subplots from the four groups among the 256 subplots at HS.

Each selected subplot was considered as a population of the focal species based on a previous study (Xu et al., 2016). Here, we considered the distance among the focal species in the subplots to eliminate any possibility that they might have come from the same plant material. Therefore, we did not choose all focal species in the selected plots. Each 20 × 20 m subplot was also divided into 16 5 × 5 m grids and 0–3 individuals were selected from these grids. A total of 601 trees of the focal species in 24 subplots were selected (Supplementary Table 1).

Three topographic variables were quantified for each 20 × 20 m subplot. These were mean elevation, slope, and convexity, and these variables were calculated from the elevation data. In each plot, the mean of four vertices was used as the elevation of the plot; the average of the angles between the plane formed by any three angles of the target plot and the horizontal plane was taken as the slope of the plot; and the elevation of the target plot minus the average elevation of the surrounding eight squares was used as the convexity of the plot (Liu et al., 2013). Nitrogen (N) and phosphorus (P) are generally considered to be the two most restricted elements in terrestrial vegetation, and they play a vital role in the composition of plant communities. Carbon (C) and potassium (K) are important nutrients for plants. pH is believed to change the physical and chemical properties of soil and affect plant growth. The soil sampling methods followed the procedures described in Sun et al. (2020), and a total of 72 soil samples were obtained from HS. In short, one sampling point at the center of the subplot, and the other four points at the middle point between the center and the four corners, respectively. After removing the litter on the surface of the soil, we used a soil auger to sample at depth of 0–15 cm. Five soil traits were determined for each sample. Soil organic carbon (SOC) was measured using the K dichromate volumetric method; total nitrogen (TN) in the soil was determined using elemental analysis (MD SpectraMax 190, United States); total phosphate (TP) was assessed by molybdenum-antimony colorimetry; and total potassium (TK) was analyzed by a flame photometer. The soil pH was determined using a glass electrode at a 1:2.5 soil/water ratio (Pansu and Gautheyrou, 2006). Each soil property was measured three times. For the first time, we determined the impact and importance of various environmental factors on SD and GD in the HS metacommunity.

A total of four species diversity indices [Tree richness (S), the Shannon-Wiener index (SD_SW), the Gini-Simpson index (SD_SI), and the Pielou index (SD_J)] were calculated by the “vegan” package (Oksanen et al., 2017). All indices were based on the 2019 re-census data. Four indices were calculated as follows (Magurran, 1988):

Tree richness:

Where S is the species richness in each 400 m² subplot, and n is the number of tree species with a diameter at breast height ≥ 1 cm.

Shannon-Wiener index:

Where p_i is the relative abundance of species i, p_i = N_i/N₀, N_i is the abundance of the i species, and N₀ is the sum of the abundances of all species.

Gini-Simpson index:

Where p_i is the relative abundance of species i, p_i = N_i/N₀.

Pielou index:

Where H is the Shannon-Wiener index, S is the tree richness.

To our knowledge, few studies have been conducted to detect SGDC based on a biodiversity dynamics observation plot, and most previous studies used amplified fragment length polymorphism (AFLP) to detect a population structure (Xu et al., 2017). In this study, we combined the multiplexed amplification of SSRs with high throughput sequencing to obtain more genotyping (Yang et al., 2019b; Lepais et al., 2020). Compared with traditional methods that rely solely on fragment size for genotyping, using next-generation sequencing technology for SSR genotyping can obtain higher levels of genetic polymorphism results (Sarhanova et al., 2018).

Leaf samples belonging to C. eyrei were taken from the selected 24 subplots and 601 samples were immediately dried in a plastic bag containing silica gel. Ten to 15 leaves were taken from each individual. We tested 69 nuclear microsatellite loci developed for Fagaceae, from which 30 (Supplementary Table 2) were selected as they had suitable levels of polymorphisms. DNA was extracted from silica-dried leaves using a DNAeasy Plant Mini Kit (Tiangen, DP350). Primers were designed and synthesized, and their DNA quality was evaluated. Then, PCR amplification was performed for each SSR locus. Optimal primers with clear bands were selected for further multiplex PCRs. Up to 30 loci primers were pooled together and adjusted using the results to achieve equal amounts of the products for each locus, to optimize primer compatibility, and to avoid undesired primer pairing. All the PCR products were pooled prior to library preparation and sequenced on the Illumina HiSeq2500 platform (PE150 bp) by Genesky Biotechnologies, Inc. (Shanghai, China). The raw reads were analyzed by pipeline, and quality control was conducted using FastQC, read merging in FLASH, and by the construction of reference alignments using Blastn. The number of SSR alleles was calculated using aligned reads and the sequence data. The SSR motif and repeat number were listed after two step adjustments, including slippage adjustment and amplification efficiency adjustment. Data analyses were based on SSR genotyping results, GenoDive version 2.0b27 (Meirmans and Van Tienderen, 2004) was used to calculate the observed number of alleles Na, effective number of alleles Ne, Shannon’s information index I, observed heterozygosity Ho, expected heterozygosity He.

All the statistical analyses were conducted using R 4.0.3 (R Development Core Team, 2020). The correlation between species diversity and genetic diversity was determined by calculating Pearson correlation coefficients. The Chi-squared Test was used to test whether the variables were significantly correlated. The effects of soil and topography on the SD and GD were investigated using hierarchical partitioning to obtain the explanatory rate of each explanatory variable (Lai and Pedro, 2020). In the hierarchical partitioning process, pH, SOC, TN, TP, and TK were set to one group and mean elevation, slope and convexity were set to another group variables. A redundancy analysis (RDA) was employed across subplots from different groups to evaluate the relationships among the soil properties, topography, and plant community structures using the “vegan” package (Oksanen et al., 2017). Five soil properties (pH, SOC, TN, TP, and TK) and three topographic variables (slope, mean elevation, and convexity) were included in the RDA. In the RDA, we first determined the Hellinger value for plant richness to eliminate the influence of extreme values and reduce the weight of highly abundant species. The Ezekiel (1930) method was used to obtain the adjusted R². A permutation test (n = 999) was used to test the significance of the RDA and the constrained ordination axes in the reduced model (Borcard et al., 2018).

The effects of topography and soil properties on SD and GD were examined using the following model:

where Y is the index for SD or GD; X₁, X₂, X₃ represent the variables of soil and topography; and ε is the random sampling error. Linear analysis with a restricted maximum likelihood was conducted to measure the relationship between soil, topography, species diversity, and genetic diversity. Before we ran the model, we conducted data conversion on response variables Y so that the model residuals met the conditions for normality. The variables SD_J, SD_SI, SD_SW, I, and Ne were SD_J/1-SD_J, SD_SI/1-SD_SI, square, quadrillion, and cubed transformed to ensure that the residuals remained normal. All model residuals were checked for normality by the “performance” package (Nakagawa et al., 2017).

The “vegan” package (Oksanen et al., 2017) was used to conduct a PCA to condense the topography and soil properties into a small number of uncorrelated variables for further path diagrams. The influence of topography and soil properties on species diversity and genetic diversity was determined using path analyses conducted based on local estimation by the “piecewise SEM” package (Lefcheck, 2016). We first ran a meta model (Supplementary Figure 1) and Shipley’s test of directed separation (d-separation test) was used to find the missing path. After adding the significant path, we then generated the final model to determine relationships among the community assembly, population of the C. eyrei and environmental heterogeneity.

A total of 87 species with a diameter at 1.3 m height (DBH) ≥ 1 cm were observed in 24 subplots (Supplementary Table 3). We found that the index of different species diversity has a large range of variation. Species richness, the Pielou index (SD_J), the Shannon-Wiener index (SD_SW), and the Simpson index (SD_SI) values for SD were from 10 to 36, 0.610 to 0.829, 1.496 to 2.890, and 0.639 to 0.921, respectively (Supplementary Table 4). Apart from C. eyrei, the most frequent species were Pinus massoniana, Vaccinium mandarinorum, Symplocos setchuensis, and Rhododendron simsii (Supplementary Table 3).

A total of 1,438 alleles at 30 loci were revealed across 601 individuals from 24 C. eyrei populations. The observed number of alleles (Na) among subplots varied from 3.395 to 8.026; the effective number of alleles (Ne) among subplots varied from 2.648 to 4.189; the Shannon’s information index (I) among subplots varied from 0.987 to 1.451; the observed heterozygosity (Ho) among subplots varied from 0.517 to 0.568; and expected heterozygosity (He) among subplots varied from 0.547 to 0.644 (Supplementary Table 4).

The elevation and slope of the 24 subplots ranged from 442.46 to 548.62 m and from 7.79 to 47.04°, respectively (Supplementary Table 1). The convexity of the subplots was also highly variable. The means for each soil property were 72.7 ± 27.6 g/kg (SOC), 3.78 ± 0.801 g/kg (TN), 0.335 ± 0.06 g/kg (TP), 21.8 ± 1.69 g/kg (TK), 4.90 ± 0.21 (pH), 19.2 ± 5.32 (C: N) (Supplementary Table 6).

There was no significant correlation between SD and GD for C. eyrei (Figure 2). However, the internal indices within SD and GD were significantly correlated. We found that both topography and soil had significant effects on genetic diversity (Figure 3B and Table 1). However, soil had no effect on community composition (Figure 3A and Table 1).

After adjustment, R squared was 0.2954, and axes 1 and 2 explained 11.56 and 7.78% of the overall variance in the plant communities, respectively (both whole model and each canonical axis, p < 0.05). Castanopsis eyrei and P. massoniana were almost negatively correlated with all soil properties and richness of other dominant species; however, both species were positively correlated with topography, indicating that dominant species, such as C. eyrei and P. massoniana, were significantly influenced by environmental heterogeneity (Figure 4).

A cross structural equation modeling (SEM) analysis was used to explore the relationship between topography, soil properties, community composition, and the genetic structure of C. eyrei. The best fit model showed that topography had an effect on both community and population. However, soil properties only influenced the observed number of alleles. A slight relationship between SD and GD was found by the SEM, and the Shannon-Wiener index (SD_SW) for community showed there was a positive relationship with observed heterozygosity (Ho) (p < 0.01, Figures 5, 6).

The SGDCs improve our understanding of the relationships between the intraspecific variation and species diversity, which can have an influence on the community assembly (Lamy et al., 2017). Before we discuss the mechanisms behind the correlation between the species diversity and genetic diversity, we should understand how the C. eyrei population and HS community respond to environmental heterogeneity.

We discovered that both soil properties and topography had effects on population structure of C. eyrei (Figure 3). Based on the linear regression, we further found that phosphorus plays a key role in the growth of C. eyrei population (Table 1). We also found that elevation is negatively correlated with phosphorus (Figure 5), and the richness of other species decreases, which means that higher elevation are more suitable for C. eyrei growth and they may become the top succession species in the area. Environmental filter theories predict that plant diversity changes along the environmental gradient, and this can be interpreted as the environment filtering out of some plants that are not suited to the local environmental conditions (Laliberte et al., 2014; Blanchard et al., 2020). This discovery is meaningful, and suggests that C. eyrei is the main constructive species in the Huangshan area. Therefore, we chose C. eyrei as the focal species. However, we found that the community only responds to soil nutrients (Figure 3), which means topography is the main factor shaping community assembly in HS.

It is theoretically believed that the parallel impact of SD and GD on the community and the population would lead to a positive SGDC (Vellend, 2005). However, in this study, no significant correlation between the two was found in the forest dynamics plot (Figure 2). It is not the first time that a mismatch SD and GD has been recorded. In an integrated study of flora, Lamy et al. (2017) found that only 11% of the studies reported that SD and GD were significantly related. There have also been related previous studies on SGDC in one location. Most studies have reported non-significant or negative correlations (Avolio and Smith, 2013; Reisch and Schmid, 2018; Wei et al., 2018). When comparing our results to these older studies (Vellend, 2004; Lamy et al., 2013; Frey et al., 2016), it must be pointed out that connectivity and area will have effects on SGDC. Discrete and continuous sampling methods often determine whether we can find a positive SGDC. No significant correlation was found in a forest dynamics plot and grassland station (Avolio and Smith, 2013; Reisch and Schmid, 2018). This result has been used to successfully account for the island biogeography theory. A delineated vegetation community is generally regarded as an island-like research system (Reisch and Schmid, 2018). A small area often leads to increased drift, decreased internal connectivity, and loss of alleles, which leads to a decrease in genetic diversity and species diversity. Therefore, a weaker or insignificant SGDC correlation was often recorded.

In line with previous studies, the effects of environmental heterogeneity on species diversity and genetic diversity were different (Figure 3; Avolio and Smith, 2013). In this study, we found that the total phosphorus had a negative effect on the C. eyrei population, which showed that low phosphorus levels restrict the growth of C. eyrei. As a result, the effective population decreases, resulting in a decrease in genetic diversity. A further analysis by SEM showed that the community assembly response to soil nutrition was significantly smaller than that of a single population to soil nutrition (Figure 6). It is worth noting that the Shannon-Wiener index was positively related to observed heterozygosity (Figure 6). The results show that community assembly affected the C. eyrei population.

When comparing our results to those of older studies, it must be pointed out that the niche of the focal species, C. eyrei, was not congruent with other species and this might be the reason why there was no relationship between species diversity and genetic diversity. In this study, we found that C. eyrei and P. massoniana richness was negatively related to most dominant species. In addition, we also found that the population richness was negatively correlated with almost all soil nutrients (Figure 5). This may be due to the fact that C. eyrei competition for resources, such as phosphorus, is more specialized than other species in the community (Lamy et al., 2017). The increase in the availability of phosphorus might increase the average population size of the focal species and decrease all other species in the community. Therefore, there was a non-significant relationship between species diversity and genetic diversity.