Negative Density Dependence Regulates Two Tree Species at Later Life Stage in a Temperate Forest

Numerous studies have demonstrated that tree survival is influenced by negative density dependence (NDD) and differences among species in shade tolerance could enhance coexistence via resource partitioning, but it is still unclear how NDD affects tree species with different shade-tolerance guilds at later life stages. In this study, we analyzed the spatial patterns for trees with dbh (diameter at breast height) ≥2 cm using the pair-correlation g(r) function to test for NDD in a temperate forest in South Korea after removing the effects of habitat heterogeneity. The analyses were implemented for the most abundant shade-tolerant (Chamaecyparis obtusa) and shade-intolerant (Quercus serrata) species. We found NDD existed for both species at later life stages. We also found Quercus serrata experienced greater NDD compared with Chamaecyparis obtusa. This study indicates that NDD regulates the two abundant tree species at later life stages and it is important to consider variation in species' shade tolerance in NDD study.

A challenge for studies of NDD is caused by lagged effects [24] which are non-fatal effects on individuals during one life stage that influence the growth and mortality rate during the following life stage [25]. Mortality patterns in seedlings can normally be analyzed through direct observation, due to their high susceptibility to natural enemies and environmental stress. However, for established trees that have lower mortality rates, long time-scales observation is required to analyze mortality pattern because NDD may not be strong enough to induce immediate mortality, but instead cause limited growth over short time-scales [24]. This may explain why few studies have been implemented to test NDD on trees at later life stages [4,[26][27][28][29]. Peters [4] found patterns consistent with NDD for saplings and trees of .75% of the species tested at sites in Pasoh, Malaysia, and BCI, Panama. Lan et al. [29] found saplings of 83.2% of species have a aggregated pattern, whereas adults of 96.2% species have a random distribution, implying that NDD can make the spatial pattern of tropical trees more regular with time. Evidence have shown that different mechanisms shape plant communities at different life stages [30][31][32], therefore, it is important to evaluate NDD at later life stages. One robust way to indirectly examine NDD on the long history of a forest is to analyze the changes in spatial patterns of trees at different life stages [20,33,34]. This is possible because the regulating mechanisms that are operating in a forest should have left a detectable spatial signature [35]. If NDD is prevalent and works constantly during the whole life-cycle of trees, the following signatures should be left on the spatial pattern of trees: trees survive better with fewer conspecific neighbors and the clustering of conspecific trees declines with time [36,37]. However, the above predictions can be obscured by habitat heterogeneity [33], because the performance of individuals can be greatly influenced by the availability of environmental resources [38][39][40], and the functional traits that are involved in plant-enemy interactions might be altered by abiotic and biotic factors [41]. This change in performance with environment can affect a population's susceptibility to herbivory to better link NDD with habitat heterogeneity [18]. For example, one population may show NDD pattern if resource availability is low and reduces the ability of individuals to survive herbivore damage, whereas one population with ample resources may not show the same effect. Therefore, habitat heterogeneity must be considered in NDD study. However, this is made difficult by the fact that many environmental covariates are difficult to quantify [20]. Getzin et al. [33] used a simple method to solve this problem. By utilizing the spatial pattern of adult trees (dbh$10 cm), they factored out large-scale habitat heterogeneity and were able to detect NDD in western hemlock populations.
Species varies in the strength of NDD, arise from differences in morphological, physiological and allocational traits [42]. Shade-tolerant, slow-growing species may experience less NDD than shade-intolerant, fast-growing species [43,44]. One of the possible mechanisms behind the differences is that shade-tolerant species tend to have higher tissue density and carbohydrate storage for increasing tolerance to herbivores [43] and recover from damage compared with shade-intolerant species [45,46]. Previous studies have shown that seedling shade tolerance is negatively correlated with disease susceptibility [40] and depression of growth and survival [47]. Also, in a tropical forest study, Comita and Hubbell [7] found seedling survival tended to be positively correlated with species abundance before controlling for shade tolerance. After controlling for shade tolerance, they found a significant negative relationship between seedling survival and species' shade tolerance. Therefore, shade tolerance should be accounted for in NDD study. However, few NDD studies implemented at later life stages considered the variation in shade tolerance among species.
The aim of this study is to examine NDD at later life stage for the most abundant shade-tolerant and shade-intolerant species at Keumsan long-term ecological research plot (Keumsan LTER) which is located in South Korea. The study was implemented by testing whether the intensity of aggregation was decreased from saplings to juveniles. The large-scale habitat heterogeneity was accounted for by utilizing the spatial pattern of adult trees. We first determined whether saplings exhibit additional aggregated patterns relative to adults, then tested whether the additional aggregation decreases from saplings to juveniles (the definitions of ''sapling'', ''juvenile'' and ''adult'' are provided in the Methods). Finally, we compared whether the shade-intolerant species experiences stronger NDD compared with the shade-tolerant species across size classes.

Study site
Our 1-ha (1006100 m) study site, Keumsan LTER (34u309N, 127u599E), is located in the southern area of the Hallyeohaesang National park in Gyeongsangnam-do, South Korea. The area is located in a warm temperate forest zone. The annual mean temperature and precipitation were 14.0uC and 2,180.0 mm, respectively, and the mean monthly temperature ranged from 21.5uC (January) and 25.6uC (August) in 2011 [48]. The study site is a secondary forest which had suffered great damage during the Korean War (1950)(1951)(1952)(1953). Thereafter, the forest has been well conserved and underwent secondary succession [49]. The Keumsan LTER is located at a 360-430 m hillside elevation with slopes ranging from 12u-28u, and is influenced by habitat heterogeneity such as edaphic gaps (e.g. gravels or wet drainage sites) [50]. The most abundant species are Quercus serrata, Chamaecyparis obtusa, Styrax japonica, Acer pseudo-sieboldianum, Carpinus tschonoskii and Stewertia pseudo-camellia. The Hallyeohaesang National Park is owned and managed by the state and its government and the location including our study area is not privately-owned. No specific permits were required for the described field studies. The field studies did not involve endangered or protected species.

Data collection
The Keumsan LTER was established in 2000. All woody stems $2 cm dbh were mapped, measured, identified to species, and tagged. The plot was recensused in 2006 and 2011. In the 2011 census, we documented 2,412 free-standing live individuals $2 cm dbh belonging to 20 families, 22 genera and 35 species. In this paper, we used data on live trees $2 cm dbh from the 2011 census for spatial point pattern analysis.
We selected Chamaecyparis obtusa (CHOB) and Quercus serrata (QUSE) as focal species. Chamaecyparis obtusa, a slow-growing, late successional species [51,52], was the most abundant shadetolerant tree species at the site, and accounted for 32.5% of the total individuals and 11.3% of total basal area (Table 1); Quercus serrata, an early successional species [53], was the most abundant shade-intolerant tree species at the site, and accounted for 10.7% of the total individuals and 53.9% of total basal area (Table 1). Young trees of the two species were usually found under different light conditions: QUSE were often found in forest gaps, while CHOB were often found in shaded areas.
All individuals of the two species were divided into three size classes as an indication for life history stages: sapling, juvenile and adult (Note the terminology differ from their original meanings, they were used only for classifying different size classes of trees and facilitating spatial point pattern analyses). For CHOB: saplings with dbh ranging from 2 to 5 cm, juveniles with dbh ranging from 5 to 10 cm and adults with dbh $10 cm; for QUSE: saplings with dbh ranging from 2 to 20 cm, juveniles with dbh ranging from 20 to 30 cm and adults with dbh $30 cm. Dbh cut-offs were selected to ensure adequate sample sizes for the spatial point pattern analysis. Both species had more than 60 individuals in each size class (Table S1 and S2). Because using different cut-offs may lead to different results, we also used all other cut-offs (with an increment of 1 cm) (on the condition that more than 60 individuals are included in each life stage) (Table S1 and S2) to implement the analyses.

Data analysis
In recent studies, spatial point pattern analysis which compares patterns of trees in different size classes has been shown to be an effective approach for testing NDD on established trees [18,20,33,34,54], with the assumption that populations of established trees are in an equilibrium stage. We utilized the bivariate pair correlation g(r) function [55,56] to implement the analyses. The g(r) function is the probability density function of the broadly used Ripley's K(r) function [57], which is calculated as: is the cumulative distribution function of the expected number of points of a pattern within the whole circle of a given radius r around a typical point of the pattern divided by the intensity l (points per unit area) of the pattern. K(r) is a cumulative distribution function where K(r) is the expected number of other points of a pattern within the whole circle of a given radius r around a typical point of the pattern divided by the intensity l of the pattern. In contrast, the g(r) is a non-cumulative distribution function in which g(r) is the expected density of other points in a ring of a given distance r around a focal point divided by the intensity l of the pattern [55]. The g(r) function can be used to estimate the strength of aggregation at specific scales. If g(r) . .1, there are more points at scale r than expected under a random distribution, which indicates an aggregated pattern at scale r. For the correction of edge effects, the translation correction method was used in the analysis [55].
We used the method of random-labeling null model within a case-control design [33] to study NDD on the two focal species from saplings to juveniles, where saplings and juveniles were used as cases (pattern i ) and adults as controls (pattern j). In the casecontrol design, the control pattern is used to account for the largescale habitat heterogeneity and intensity of seed rain [33,34]. This method is based on two assumptions: 1) small-scale patterns are usually attributed to plant-plant interactions, whereas large-scale patterns are usually attributed to habitat heterogeneity [58] and 2) adult trees have undergone excessive thinning over time due to habitat heterogeneity and most likely represent those which lived in sites most favorable for the species. Therefore, the large-scale pattern (e.g. at r .10 m) of adults can be used as an indicator of environmental habitat preferences and can be used to control for the large-scale habitat heterogeneity [59,60]. Under the null hypothesis of random-labeling, labels (case or control) are assigned to points randomly, conditioning the observed locations of the points in the joined patterns of cases and controls [61]. If the null hypothesis is rejected, the case-control approach may identify specific factors (i.e. NDD) other than habitat heterogeneity that may influence the spatial pattern of trees, by comparing the differences between case patterns and control patterns of trees in different life stages [18,33]. The g(r) functions are invariant under random thinning of trees, hence we would expect g ii (r) = g jj (r) = g ij (r) = g ji (r). In this study, NDD is examined by using g ii (r) -g ij (r) as test statistics [33,34,62].
According to the prediction of NDD, tree survival would be low if individuals are growing in high density patches of conspecifics. Therefore, if there is strong NDD in saplings, the surviving juveniles would be less aggregated than saplings. To test this hypothesis, we studied the difference in the degree of aggregation between saplings and juveniles. We used A i (r) = g ii (r) -g ij (r) as the test statistic to determine whether cases i show an additional pattern that is independent of the controls j. If A i (r) ..0, cases can be said to exhibit additional aggregated patterns relative to adults, irrespective of whether habitat heterogeneity is present or not [62,63]. The change in additional aggregation from saplings to juveniles at scale r can be expressed by the formula: T(r) = A saplings (r) -A juveniles (r). For a species, if A saplings (r) ..0 and T(r) ..0, we would infer that NDD existed in saplings. If this is the case, T(r) can be used to indicate the strength of NDD at different scales.
All spatial point pattern analyses were implemented with the ''spatstat'' package of R [64]. We focused on the scale of 0-20 m for the analyses above, assuming that tree-tree interactions could be efficiently indicated by this scale [65]. We performed 199 Monte Carlo simulations of the random labeling null model and used the 5th-lowest and 5th-highest values (i.e., extreme 0.25% simulated cases at either end) as simulation envelopes. However, this simulation inference yields an underestimated Type I error rate because the tests are performed at different concurrent scales [66]. We combined the simulation method with a goodness-of-fit test (GOF) [67]. Significant deviations were only determined for those data sets where the observed GOF's P value (P GOF ) were less than 0.025 [58,66].

Diameter distributions
The diameter distribution of CHOB was strongly right-skewed with 89% trees less than 10 cm dbh (Fig. 1a); this is a pattern typical of late-successional species that have a large number of suppressed young trees. In contrast, the diameter distribution of QUSE was symmetrical with a lack of small trees (Fig. 1b); this is a pattern typical of early-successional species that establish as approximately even-aged cohorts in large forest gaps. For both CHOB and QUSE, dead trees were concentrated in the smaller diameter classes (Fig. 1). Young trees of QUSE had a much greater mortality rate compared with CHOB (Fig. 1). For trees that were

Negative density dependence
For both CHOB and QUSE, saplings showed additional aggregated patterns relative to adults, as their test statistic A saplings (r) ..0 (P GOF = 0.005 for CHOB and QUSE) (Fig. 2a,   b). For CHOB, significant deviations were found at the scale of 0-3, 4-4.5, 5.5-6.5, 8.5, 14-16.5, 17.5-19.5 m (Fig. 2a); for QUSE, significant deviations were found at all test scales (Fig. 2b). In juveniles, these aggregated patterns diminished for both species, as their test statistics demonstrate (A juveniles (r)<0) (Fig. 2c, d). These results indicate that both focal species were influenced by NDD at the sapling stage. For both species, the decrease in additional aggregation was the greatest at the smallest scale and had a trend The test statistic A j (r) (A saplings (r) or A juveniles (r)) evaluates if there is an additional aggregation within cases (saplings or juveniles) that is independent from adults. A j (r) ..0 means there is an additional aggregation within cases relative to adults. A j (r) <0 means patterns of cases and adults are created by the same stochastic process. The scale r is the radius around a focal tree of the pattern. Dotted lines were the 95% simulation envelopes constructed using the 5th-lowest and 5th-highest value of 199 Monte Carlo simulations of the null model. Comparison of the strength of NDD effect experienced by the two focal species (e) and Comparison of the strength of NDD experienced by the two focal species using all dbh cut-offs (f). Dbh cut-off points (a, b) are the cut-offs for defining life stages of trees: sapling ($2 and ,a), juvenile ($a and ,b) and adult ($b). Units of cut-off points are centimeters. The cut-off point lists can be found in Table S1 and S2. *: significant difference by the goodness-of-fit test (P GOF ,0.025). doi:10.1371/journal.pone.0103344.g002 of decreasing with increasing scales (Fig. 2e), indicating that NDD occurred predominantly at close distances among neighbors. QUSE experienced stronger NDD in the sapling stage as T(r) is greater for QUSE compared with CHOB at all test scales (Fig. 2e). For all dbh cut-offs, the mean T(r) values were greater for QUSE than CHOB (Fig. 2f, Table S1 and S2), indicating the above finding was not influenced by the selection of different dbh cut-offs.

Discussion
NDD has been shown to play an important role in the population dynamics of tree communities in numerous studies [1,10,31,[68][69][70]. Our study examined the spatial patterns and tested NDD for two abundant tree species with different shade tolerance at later life stages.
Plants are most susceptible to natural enemies and abiotic stresses in early life stages, therefore, most demographic winnowing of propagules takes place at early life stages [71][72][73]. Thus, one might expect NDD to be subdued at later life stages. However, recent studies pointed out that NDD also exists at later life stages [4,15,18,22,54,74]. For example, Zhu et al. [18] studied 47 common tree species (dbh $ 1 cm) and found 39 of them exhibited conspecific density-dependent thinning in a subtropical forest by comparing the strength of aggregation from saplings to juveniles. Piao et al. [54] studied 15 common tree species (dbh $ 1 cm) in a temperate forest and found 11 of them showed conspecific densitydependent thinning, also by comparing the strength of aggregation from saplings to juveniles. Carson et al. [75] suggested that ecologists should focus on NDD occurring at later life stages in order to get a complete understanding of the effect in promoting species coexistence. In our later life stage study, we found NDD regulated both focal species. We also found NDD was the greatest at small scales and decreased with increasing scales for saplings. These results may be due to the decline in aggregation of saplings with increasing distance from parent trees, which had likely resulted from dispersal limitation [76]. Those findings are consistent with the distance-dependent model which predicts that survival is negatively correlated with the distance from the conspecific adult [77]. Distance dependence and NDD leads to a decreased probability that a species will replace itself locally and therefore promotes species coexistence [5,77].
Shade-tolerant tree species have been shown to have better defenses against natural enemies than shade-intolerant tree species [42,43]. Therefore, one could predict that shade-intolerant tree species experience greater NDD than shade-tolerant tree species. For example, Augspurger and Kelly [11] found that seedling disease susceptibility was negatively correlated with shade tolerance. Comita et al. [7] found that variations in species' shade tolerance could mask the negative relationship between seedling survival and abundance. Our result is consistent with the above prediction: we found that QUSE (a shade-intolerant and early-successional species) experienced greater NDD compared with CHOB (a shade-tolerant and late-successional species). Because this result can be caused by the great differences in the size class cut-offs between the two species, we used all possible size class cutoffs to implement the analyses and therefore our results are not likely caused by the specific size cut-offs selected. Our observation and previous studies suggest that attempting to quantify the contribution of NDD should take into account the variation in species' shade tolerance.
The stronger NDD experienced by QUSE may play an important role in the replacement of QUSE by other latesuccessional species such as CHOB, as self-thinning of earlysuccessional trees could progressively open up the forest canopy and create open spaces that promote the survival of latesuccessional trees [20]. However, further studies should be implemented to examine whether the subsequent survival of late-successional trees such as CHOB are higher in those open spaces created by QUSE.

Conclusions
Our study tested NDD for two abundant species found in the Keumsan LTER. We found NDD on the two abundant tree species at later life stages. We also found the shade-intolerant species suffered stronger NDD compared with the shade-tolerant species. Our findings highlighted the importance of testing NDD at later life stage for improved understanding of species coexistence. We also highlighted the importance of considering species shade tolerance in NDD study.

Supporting Information
Table S1 List of dbh cut-off points for defining life stages of Chamaecyparis obtusa: sapling ($2 and ,a), juvenile ($a and ,b) and adult ($b).