We analyzed the self-thinning relationship of Chinese fir-dominated mixed forests using National Forest Inventory (NFI) data for mainland China. Data for every 5 years of observations from 1991 to 2016 were obtained for Fujian, Jiangxi, Zhejiang, and Anhui provinces, which included 61.1% of areas with Chinese fir in China ( Table 1 ; Figure 1 ).

National Forest Inventory plots of Fujian (a 667 m² square), Jiangxi (a 667 m² square), Zhejiang (an 800 m² square), and Anhui (a 667 m² square) Provinces were established. The sample area was converted to units of hectares. The plots for this study were selected using the following criteria: (1) pure Chinese fir stands were unthinned, and seedlings were planted; (2) the volume of Chinese firs comprised more than 30% of the total stand; (3) the number of trees per hectare (TPH) continuously decreased with the quadratic mean diameter (QMD) every 5 years of observations according to the Andrews et al. (2018) method; (4) environmental disturbance (e.g., windthrow, pest, and disease) and artificial damage (e.g., illegal thinning) were absent; and (5) the regeneration of broad-leaved seedlings was natural. The selected plots were inevitably under self-thinning conditions. A total of 591 plots were selected, including 149, 146, 117, and 179 study sites in Fujian, Jiangxi, Zhejiang, and Anhui provinces, respectively ( Figure 2 ).

Measurements were taken from all trees above 5 cm in diameter at breast height (DBH) in every permanent plot; data collected included tree species, DBH, tree height, and tree status (living/dead). DBH was measured to the nearest 1 mm with diameter tape at a tree height of 1.3 m. Tree heights were measured (to the nearest 10 cm) with a Vertex (Haglöf Sweden). The recorded plot factors included soil type, latitude, longitude, altitude, slope aspect, slope position, slope gradient, soil thickness, soil texture, and humus thickness.

The number of living trees (DBH≥5 cm) per hectare is a basic index of stand density. QMD is a basic index that reflects the mean tree size in sample plots. It was calculated using (Equation 1).

where ( D B H ) i is the DBH of the i tree, and N is the total number of trees in the sample plot.

A site index (SI), representing site conditions, was estimated for each plot using the average height of the most dominant trees (95^th percentile and above) per hectare at the reference age of 20 years.

To determine tree species dominance in Chinese fir-dominated mixed forests, we calculated the ratio of the Chinese fir stand basal area to the total stand basal area (BA _{Chinese fir} /BA _total) in sample plots (Sun et al., 2011). When BA _{Chinese fir}/BA _total<0.5, the broad-leaved small trees were dominant in stands; however, when BA _{Chinese fir}/BA _total >0.5, Chinese fir were dominant in stands. In addition, BA _{Chinese fir} /BA _total /(1- BA _{Chinese fir} /BA _total) was used to describe the mixing proportion (MP) of each plot.

In this study, we used latitude to represent precipitation and temperature in different provinces to clarify the effect of climate changes on maximum size–density relationships in Chinese fir-dominated mixed forests.

The correlation was calculated using the “corrplot” package in R software (R Core Team, 2023), and the “ggplot” package was used to plot the results. Stochastic frontier analysis was performed using FRONTIER 4.1 (Coelli, 1996). The adjusted coefficient of determination (R_adj ²), root mean square error (RMSE), and Pearson correlation coefficient between estimated and observed values (r) were used as metrics to assess the efficacy of the model’s fit.

Bi et al. (2000) introduced stochastic frontier analysis to estimate maximum size–density relationships of even-aged Pinus radiata stands. The generalized expression (Equation 2) of the self-thinning relationship is as follows (Bi et al., 2000):

where Y is the observed value of the dependent variable related to stand growth; X₁,…X_k are the independent variables affecting stand growth; β ₁, β ₂…β_k are the exponents of the independent variables; A is the parameter of the model to be estimated; and e ϑ  and  e − μ are the two error exponents. Taking the logarithm, the expression (Equation 2) becomes (Equation 3) (Bi et al., 2000):

where α = l n A   , β is a vector of parameters. The error term, ε = ϑ − μ , is a compound random variable with two components, and each is assumed to be independently and identically distributed across observations.

Bi (2001) incorporated site productivity in the stochastic frontier function for a generalized expression (Equation 4) of the self-thinning relationship:

where N is the number of trees, S is the relative site index, A is the intercept parameter, and β ₁ and β ₂ are parameters. The variables ϑ and μ are two random variables.

Taking the logarithm, the expression (Equation 4) becomes (Equation 5) (Bi, 2001):

where α = l n A , the error term ϵ is a compound random variable, ε = ϑ − μ . We introduced the mixing proportion and latitude into the model in expression (5) to analyze their effects on the number of trees, the slope and intercept of the self-thinning relationship, and the correlations between variables and determine the maximum stand density of the mixed-species forests. The model (Equation 6) can be expressed as follows (Kimsey et al., 2019):

where TPH is the number of living trees per hectare; Q M D is the quadratic mean diameter of each plot; ( B A C h i n e s e   f i r B A t o t a l ) / ( 1 − B A C h i n e s e   f i r B A t o t a l ) represents the mixing proportion; latitude represents climate changes (temperature and precipitation); SI is the site index; α is the intercept parameter; β 1 , β 2 , β 3 , and β 4 are parameters; and ϵ is the error term.

These models (Equation 7), with similar structures but different variable effects, resulted in a single, parsimonious model adaptable to each species that could be generally stated as (Kimsey et al., 2019):

where the TPH for any given species is a function of QMD (PP1), species basal area proportion ( B A C h i n e s e   f i r B A t o t a l ) / ( 1 − B A C h i n e s e   f i r B A t o t a l ) (PP2), latitude (PP3), and SI (PP4). Species models (PP1–PP4) are additive in nature to clarify relative sequential effects on TPH.

On the log-log scale, the number of living trees per hectare (TPH) was negatively correlated with the quadratic mean diameter (QMD) of Chinese fir-dominated mixed forests in the different provinces ( Table 2 ). That is, as QMD increased, TPH decreased by 1.605. Compared with the theoretical self-thinning slope of -1.605 in an even-aged pure stand, the estimated slope did not deviate from -1.605, albeit its trajectory was flatter for Chinese fir-dominated mixed forests (PP1 in Table 2 ).

Maximum size–density relationships remained consistent as the mixing proportion, latitude, and site conditions increased. Compared with the self-thinning trajectory for TPH~QMD (PP1 in Table 2 ), the self-thinning slope was flat, even when the mixing proportion factor was included (PP2 in Table 2 ). It became steeper when latitude was included (PP3 in Table 2 ), and even more so when the site condition factor was included (PP4 in Table 2 ).

In the PP4 model, the slopes of the maximum size–density relationships were shallower in Fujian, Jiangxi, and Zhejiang provinces and steeper in Anhui Province than in even-aged pure stands, yet the slopes of the maximum size-density relationships for all provinces (PP4 in Table 2 ) were consistent with Reineke’s assumed universal slope of -1.605 at the 95% confidence level, with no deviation from the slope of the self-thinning relationship in even-aged pure stands. The inclusion of interaction terms (e.g., ln(QMD)*SI) did not significantly affect the model intercept or slope coefficients (Kimsey et al., 2019). Interaction terms were not included in the PP4 model.

The correlation coefficients were lower between the mixing proportion (or latitude) and TPH than between the other variables (e.g., QMD and SI) ( Figure 3 ). However, the number of remaining Chinese firs and naturally renewed broad-leaved tree species in different plots varied with the mixing proportion. Moreover, temperature and precipitation conditions at different latitudes in southern China vary, and this affects the self-thinning of Chinese fir and the regeneration of broad-leaved tree species. The addition of the mixing proportion and latitude variables to the model led to an increase in the R_adj ² of the model and a decrease in the RMSE ( Table 2 ). Therefore, the addition of the mixing proportion and latitude variables is important for clarifying maximum size–density relationships in Chinese fir-dominated mixed forests.

Mortality rates declined as the mixing proportion of Chinese fir-dominated mixed forests increased (PP2 in Table 2 ). In the PP4 model, the mortality rate of Chinese fir-dominated mixed forests decreased by 0.006, 0.028, 0.051, and 0.044 in Fujian, Jiangxi, Zhejiang, and Anhui provinces, respectively, for every unit increase in mixing proportion at the same latitude and site conditions ( Table 2 ).

The PP3 model showed that Chinese fir-dominated mixed forests had a higher carrying capacity in lower-latitude provinces than in higher-latitude provinces ( Table 2 ). The mortality rate of Chinese fir-dominated mixed forests increased by 1.158, 1.220, 1.152, and 1.000 in Fujian, Jiangxi, Zhejiang, and Anhui provinces for every unit increase in latitude at the same mixing proportion and site conditions (PP4 in Table 2 ). As the latitude increased (from Fujian Province to Anhui Province), the slope of the maximum size–density relationship gradually became steeper, and the number of living trees per hectare changed in parallel (PP3 in Table 2 ).

Forest site had a negative effect on maximum size–density relationships (PP4 in Table 2 ). Regardless of the mixing proportion and latitude, the mortality rate of Chinese fir-dominated mixed forests in Fujian, Jiangxi, Zhejiang, and Anhui provinces increased by 0.030, 0.013, 0.041, and 0.015 for every unit increase in site conditions (PP4 in Table 2 ).

The magnitudes of the effects of mixing proportion, latitude, and site variables on maximum size–density relationships varied ( Figure 4 ). Latitude had a stronger negative effect on maximum size–density relationships than mixing proportion and site conditions. The slope of the maximum size–density relationship of the Chinese fir-dominated mixed forests did not deviate from the slope of the self-thinning relationship of the pure forests.

The ecological carrying capacity for any given forest is driven by size–density relationships (Shaw, 2006). The difference in the slope of the limiting relationship is an indication of differences in the potential of a stand to withstand an understocked site and crowding conditions. The slope of the maximum size–density relationship of Chinese fir-dominated mixed forests was flatter than that of even-aged pure Chinese fir stands (PP1 in Table 2 ), and changes in the number of living trees per hectare were mainly driven by the appearance of natural broadleaved tree species and the death of a few Chinese fir trees (Sun et al., 2011). Flatter slopes indicate a greater capacity for a stand to withstand low-density conditions as diameter growth increases, and a poor ability of trees to survive crowding conditions (Gadow, 1986).

Our study showed that the slope of the maximum size–density relationships of the PP4 model in each province was -1.605 (PP4 in Table 2 ), which did not differ from the slope of the self-thinning trajectory in the even-aged pure stands. The characteristics of the Chinese fir-dominated mixed forests, including the uneven age of the stands and the presence of multiple layers, especially the vertical crown layers, were consistent with the self-thinning trajectory. The lack of a difference in the slope of the maximum size–density relationships and the slope of the self-thinning trajectory of even-aged pure stands can be explained by a dynamic negative feedback mechanism that regulates stand density; the negative feedback between trees reflects the increasing dominance of the upper Chinese fir crown layer and the continuous death of lower broadleaved tree species during the self-thinning process. The absence of a slope difference can also be explained by the observation that light was the main factor directly affecting the self-thinning trajectory of mixed forests with high stand densities (Morris, 1999).

In this study, Chinese fir-dominated mixed forests comprised even-aged pure Chinese fir stands and natural tree species. Chinese fir occupied the upper crown layer, and naturally regenerated broadleaved tree species were located in the secondary crown layer. The main competitive interactions between Chinese fir and natural tree species are associated with differences in tree size. As mixed forests age, the growth space and habitat resources required by trees increase continuously. Asymmetric competition among trees intensifies (Xue and Hagihara, 1999; Stoll et al., 2002; Ogawa, 2005), and the interactions between adjacent individuals are stronger than interactions between non-adjacent individuals. The larger Chinese fir trees outcompete the smaller broadleaved tree species for light; this eventually leads to deviation from the light compensation point and the death of broadleaved trees (Weiner and Thomas, 1986), which further enhances the dominance of Chinese fir.

The death of trees in mixed forests is caused by a combination of the lack of light and limited growth space (Long et al., 2022). When even-aged pure stands of Chinese fir are converted to Chinese fir-dominated mixed forests, self-thinning of Chinese fir occurs, which involves the death of less dominant Chinese fir trees and the formation of canopy gaps. Canopy gaps allowed some broadleaved tree species, especially evergreen shade-tolerant tree species (e.g., Schima superba Gardner & Champ, Cinnamomum camphora (Linn) Presl, and Phoebe zhennan S. K. Lee & F. N. Wei), to grow in the secondary layer.

Chinese fir is a light-demanding tree that occupies the upper canopy when other shade-tolerant tree species are present, and shade-tolerant broadleaved tree species are located in the lower canopy. Niche differentiation of the vertical space can improve the overall light use efficiency of forests (Dănescu et al., 2016). The canopy interception of Chinese fir decreased during the self-thinning process, given that the precipitation in the atmosphere can penetrate through the forest layer and directly into the soil via canopy gaps, which enhances the supply of soil nutrients and water.

Temperature and rainfall play a significant role in modifying both the slope and intercept of the maximum size–density relationship in a given species (Kweon and Comeau, 2017; Andrews et al., 2018; de Prado et al., 2020). Latitudinal differences are closely related to climate conditions, which have implications for the growth and distribution of tree species. The growth of trees is mainly limited by temperature in high-latitude areas. By contrast, the effect of water availability on the growth of trees becomes increasingly important as the latitude decreases (Babst et al., 2013). Chinese fir is a shade-intolerant tree species with a low growth rate, and self-thinning is delayed in high-latitude areas. Decreases in temperature impede the regeneration of broadleaved tree species in the lower canopy and reduce stand density (Wu, 1984).

When water availability is not limited, the survival probability of tree species is high in warm areas. Some studies have shown that warm, dry forest sites have lower stand densities than cool, moist forest sites, which have higher numbers of shade-tolerant species (Weiskittel et al., 2009; Ducey et al., 2017). However, warmer conditions can enhance Chinese fir growth and decrease tree mortality (Wu, 1984). In addition, warm climates likely promote more rapid increases in photosynthetic activity than in respiration rates, which enhances the carbon assimilation rate (Way and Oren, 2010; Lines et al., 2012).

Site conditions do not affect the maximum size–density relationships of mixed forests, but they have positive effects on average tree size growth as forests develop; stand density decreases as forests experience competition-induced mortality (Harms et al., 2000; Pittman and Turnblom, 2003; VanderSchaaf and Burkhart, 2008). Canopy closure can occur rapidly in forests, and competition among trees can be particularly intense at fertile sites, which can result in canopy gaps (Yan et al., 2021). However, living trees, especially small broadleaved tree species, can rapidly colonize canopy gaps at fertile forest sites.

As site conditions increased, the number of living trees per hectare decreased (PP4 in Table 2 ). Trees with high fertility allocated more biomass to crown growth, which accelerated the competition among tree species for light, and this resulted in the death of competitively inferior tree species (Morris, 2002; Hautier et al., 2009; Campoe et al., 2013). The asymmetric competition for light is intensified at fertile sites, which results in the death of small Chinese fir trees. Moreover, the differences in tree size between Chinese fir and broadleaved tree species are particularly pronounced at highly fertile sites.