Changes in Phylogenetic Community Structure of the Seedling Layer Following Hurricane Disturbance in a Human-Impacted Tropical Forest

: Disturbance plays a key role in shaping forest composition and diversity. We used a community phylogeny and long-term forest dynamics data to investigate biotic and abiotic factors shaping tropical forest regeneration following both human and natural disturbance. Speciﬁcally, we examined shifts in seedling phylogenetic and functional (i.e., seed mass) community structure over a decade following a major hurricane in a human-impacted forest in Puerto Rico. Phylogenetic relatedness of the seedling community decreased in the ﬁrst ﬁve years post-hurricane and then increased, largely driven by changes in the abundance of a common palm species. Functional structure (based on seed mass) became increasingly clustered through time, due to canopy closure causing small-seeded, light-demanding species to decline in abundance. Seedling neighbor density and phylogenetic relatedness negatively affected seedling survival, which likely acted to reduce phylogenetic relatedness within seedling plots. Across the study site, areas impacted in the past by high-intensity land use had lower or similar phylogenetic relatedness of seedling communities than low-intensity past land use areas, reﬂecting interactive effects of human and natural disturbance. Our study demonstrates how phylogenetic and functional information offer insights into the role of biotic and abiotic factors structuring forest recovery following disturbance. Software, L.S.C.; Investigation, J.T., J.K.Z., J.F.-M.; Data Writing—Original Preparation, Writing—Review M.N.U., M.U., J.F.-M., W.J.K., N.G.S., J.T. J.K.Z.; Visualization, L.S.C.; Supervision, M.U., J.T. J.K.Z.; J.K.Z. and J.T.; J.K.Z., J.T., M.U. and W.J.K.


Introduction
Disturbance plays a key role in shaping ecological communities by increasing spatial and temporal variation in resource availability and providing opportunities for recruitment [1,2]. Species often differ in their responses to disturbance, employing divergent strategies that allow them to take advantage of heterogeneous environmental conditions [3,4]. This variation in life history strategy is hypothesized Puerto Rico changed over the decade following a major hurricane (Hurricane Georges in 1998). To better understand the abiotic and biotic interactions driving changes in phylogenetic community structure, we also examined community-wide shifts in seed mass values, as well as effects of local neighborhood composition on individual seedling survival. In addition, we took advantage of spatial variation in the intensity of past human activity within the LFDP to examine the signature of land use history on phylogenetic community structure. We focused on the seedling stage because this life stage responds strongly and rapidly to disturbance [29,38] and because processes occurring during early life stages ultimately shape community composition and play a key role in maintaining diversity in species-rich plant communities [3]. During hurricanes, tree mortality is often low [39,40], but damage to large trees can be extensive (i.e., defoliation, broken branches), altering conditions beneath the canopy and resulting in significant effects on seedlings growing in the understory [29,41,42]. In addition, in forests recovering from past human disturbance, studies of the seedling layer can offer insights into the successional trajectory of forest recovery [43].
In this study, we examined changes in seedling phylogenetic community structure during forest recovery from disturbance and test several hypotheses. After the hurricane, as the canopy recovered and understory light levels dropped, we expected a decrease over time in the abundance of small-seeded, shade-intolerant species in the seedling layer resulting in increasing phylogenetic and functional clustering. However, at the same time, we predicted an initial decrease in phylogenetic clustering as seedling densities thinned out due to non-random mortality with respect to the phylogenetic relatedness of neighbors (i.e., phylogenetic negative density dependence). Finally, since human disturbance generally results in more homogenous areas dominated by closely related, disturbance-adapted species, we predicted that areas with a history of high-intensity land use would exhibit higher phylogenetic clustering compared to adjacent areas that had been less disturbed by past human activity.

Study Site
The study was conducted in the 16 ha Luquillo Forest Dynamics Plot (LFDP) in Northeastern Puerto Rico (18 • 20 N, 65 • 49 W). Annual rainfall at the site averages 3500 mm, with all months receiving >200 mm of rain, and daily average maximum and minimum air temperatures are 25.2 • C and 20.5 • C, respectively [44]. All freestanding woody stems ≥1 cm dbh (diameter at 1.3 m above ground) in the LFDP have been tagged, mapped, identified to species, and measured at~five-year intervals since 1990.
The Luquillo forest has a well-documented history of hurricane disturbance. In the two decades before our study, two major storms struck the forest, Hurricane Hugo (category 4) in September 1989 [39] and Hurricane Georges (category 3) in September 1998 [45]. In this study we focused on regeneration following Georges, which caused significant damage to trees in the Luquillo forest, altering understory conditions. Specifically, the mean and heterogeneity of light levels on the forest floor were elevated for 1-3 years following Georges, resulting in high seedling densities that persisted for 4-6 years [29].
In addition to the natural disturbance, the LFDP forest has experienced anthropogenic disturbance in the form of clear-cut logging and agriculture that occurred before 1934 when the USDA Forest Service purchased the land. Selective logging occurred in the LFDP until 1953 [46]. This past land-use history has significantly altered the species composition of the forest [47][48][49]. For this study, we divided the LFDP into two areas based on past land use, using information from aerial photographs taken in 1936. The northern portion of the plot (~2/3 rds of the LFDP) was categorized as high land use history and the southern area was categorized as low land use history (see [29,46] for more details).

Seedling Census
To monitor seedling dynamics following Hurricane Georges, 150 permanent 1 × 2 m plots were established at~20 m intervals along six north-south running transects spaced 60 m apart [29,50]. A census was conducted in March-April 2000, during which all tree and shrub seedlings ≥10 cm tall and <1 cm dbh within the seedling plots were tagged, measured for height and root collar diameter, and identified (99.99% of individuals were positively identified to species). Palm seedlings were included if ≥10 cm tall and their newest frond arose <1.3 m from the ground. Recensuses of the seedling plots took place in April-June 2002, August-October 2004, February-March 2007, and January-February 2008. Mean daily rainfall in mm (±sd) in the four census intervals was 9.03 (±17.17), 10.80 (±18.28), 11.03 (±8.60), 9.98 (±15.27), respectively. From 2000 to 2008, 13,957 individuals of 78 tree, palm, and shrub species were measured and identified in the 150 seedling plots.

Analysis of Phylogenetic Community Structure
An ultrametric community phylogeny at the species-level was generated for the LFDP using a multi-locus DNA barcode library of three plastid markers, shown in Figure S1: rbcL, matK, and trnH-psbA. The molecular phylogeny contained 141 species and was constructed using a super matrix approach and constraint tree at the level of order. A detailed description of the construction of the phylogeny can be found in [37,51]. We calculated phylogenetic community structure of the seedling layer (all tree, shrub, and palm seedlings ≥10 and <130 cm tall) for each of the five seedling censuses. We excluded seedlings of four species (Casearia decandra, Lasianthus lanceolatus, Phytolacca rivinoides, and Psychotria maleolens) because they either do not grow to be ≥1cm dbh or because individuals ≥1 cm dbh were not found when tissue was collected to construct the phylogeny. Together these excluded species contributed only 0.3% of individuals to the seedling census, and are therefore unlikely to influence the results.
Various metrics for assessing phylogenetic community structure have been proposed [52]. We chose two of the most commonly utilized: Net Relatedness Index (NRI) and Nearest Taxon Index (NTI) [53]. NRI and NTI are calculated based on the standardized effect size of mean pairwise distances (MPD) separating taxa and mean nearest taxon distances (MNTD) in communities, respectively. NRI is more sensitive to tree-wide patterns of phylogenetic clustering and evenness, whereas NTI is more strongly influenced by patterns closer to the tips of the phylogeny (i.e., within families or genera). We calculated NRI and NTI for each individual 2 × 1 m seedling plot containing at least two species (n = 136-146 plots, depending on census year), using abundance data to assess phylogenetic relatedness among individuals within each plot. We used random draws from the pool of species occurring in the community phylogeny as our null model (i.e., all species ≥1 cm dbh in the 16 ha LFDP), which maintains the observed species richness but does not constrain species occurrence frequencies [54]. This null model is appropriate for examining changes in phylogenetic structure in our system, over the spatial and temporal scales of interest, since seedlings were expected to be a subset of the adult trees present in the 16 ha LFDP. Given that the phylogenetic signal in species abundance can influence the patterns of phylogenetic structure [55], we tested for, and found no evidence of, a phylogenetic signal in species abundance in any census (all values of Blomberg's K < 0.61 and not significant based on randomization tests with α = 0.05). Phylogenetic analyses were carried out using the 'picante' package for R [56].
To examine temporal shifts in phylogenetic structure of the seedling community over the decade following Hurricane Georges, we used non-parametric Kruskal-Wallis rank sum tests to determine if values of NRI and NTI varied significantly among census years. When significant variation among years was found, we then tested for significant differences between individual pairs of years using paired Wilcoxon signed rank tests with p-values adjusted for multiple comparisons using the False Discovery Rate (FDR) method [57].
To examine the effect of land use history on phylogenetic community structure, we compared the NTI and NRI values for seedling plots in the high-intensity land use portion of LFDP versus seedling plots in the low-intensity land use portion of LFDP. We analyzed each census year separately, using Wilcoxon signed rank tests with p-values adjusted for multiple comparisons using the FDR method.
The most dominant species in the Luquillo forest is the palm Prestoea acuminata. The other palm species in the LFDP, Roystonia borinquena, is also relatively common. As palms, they are distantly related to most other species in the forest, and may therefore have a disproportionately large impact on patterns of phylogenetic community structure. Therefore, we recalculated all NRI and NTI values and repeated the above analyses after dropping palm species from the data set.

Phylogenetic Neighborhood Effects on Seedling Survival
To examine whether the phylogenetic structure of local neighborhoods influenced seedling dynamics, we used generalized linear mixed models (GLMMs) with binomial errors to model individual seedling survival as a function of the density and phylogenetic relatedness of seedlings within the same 2 × 1 m seedling plot for each census interval. The phylogenetic neighborhood of each individual was quantified as , where S is the number of species present in the same seedling plot as the focal individual, N is the number of individuals of species i in that plot, and P is the phylogenetic relatedness of species i to the focal species j. P was calculated by computing the cophenetic distances between species and then scaling values from 0 to 1, with 1 indicating a conspecific. Therefore, the effect of a neighbor was strongest if it was the same species as the focal individual, and declined with decreasing phylogenetic relatedness as a function of the exponent d. We initially ran models using several different values of d (0.5, 1, 2, 3, 4), but found that these values all gave qualitatively similar results. Therefore, we only present results from models with d = 1 (i.e., a linear decrease with relatedness). In addition to the phylogenetic neighborhood, seedling height was included as a fixed effect. Species and seedling plot were included as random effects to account for species-specific variation in mortality rates and spatial variation in survival unrelated to neighbor density (i.e., spatial autocorrelation), respectively. We also ran the model including canopy openness (light) as a random effect that varied among species; however, it did not alter the effect of phylogenetic neighborhood on survival (data not shown). In order to determine whether significant effects of phylogenetic neighborhood on survival were predominantly due to effects of conspecific neighbors or total neighbor density, we used Akaike information criterion (AIC) to compare the phylogenetic neighborhood models to models that substituted the density of either conspecific seedling neighbors or all seedlings neighbors in the 2 × 1 m plot in place of the phylogenetic neighborhood. GLMMs were fit using the 'lme4' package for R [58].

Temporal Variation in Community Mean Seed Mass
We used data on species' mean seed mass to examine whether effects of disturbance on phylogenetic community structure were likely driven by species' differential life history strategies. As previously reported [22], closely related species in the Luquillo forest have significantly more similar values of seed mass than expected by chance (i.e., significant phylogenetic signal in seed mass). For each species, between 20 and 200 seeds (depending on seed size and species abundance) were collected from a network of 120 seed traps located throughout the 16 ha LFDP and taken back to the lab where they were dried at 60 • C for several hours before weighing. Seed mass varied over several orders of magnitude among species and was highly skewed, therefore values were log-transformed before use in analysis. Each individual seedling was assigned the mean seed mass of its species. We then calculated the mean and variance of seed mass across all seedlings in the 150 seedlings plots in each census to assess community-wide changes over time in response to post-hurricane canopy dynamics. In addition, we constructed a trait dendrogram, analogous to a phylogenetic tree, and calculated abundance-weighted NRI and NTI values for each seedling plot in each census to determine whether patterns of seed mass structure mirrored patterns of phylogenetic structure. The trait dendrogram was constructed by first transforming the seed mass trait data, then using hierarchical clustering (UPGMA) on the trait distance matrix. We tested for significant variation among years in seed mass NRI/NTI values using Kruskal-Wallis tests. When significant variation among years was detected, differences between individual years were then tested using paired Wilcoxon-signed rank tests with p-values adjusted for multiple comparisons (see above). Seed mass analyses were performed using the 'picante' package in R [58]. All analyses were run using the R statistical package [59].

Phylogenetic Community Structure Following Hurricane Disturbance
We found significant variation in phylogenetic structure among years, as seen in Figure 1 Figure 1a,b. Thus, over the study period, individuals within seedling plots tended to be most closely related in the first census following hurricane disturbance and least related~five years after the hurricane.

Phylogenetic Community Structure Following Hurricane Disturbance
We found significant variation in phylogenetic structure among years, as seen in Figure 1, for both indices (For NRI, Kruskal-Wallis chi-squared = 26.7238, df = 4, p-value < 0.001; for NTI Kruskal-Wallis chi-squared = 15.8741, df = 4, p-value = 0.003). Values of NRI tended to be negative in all years, while values of NTI tended to be positive. However, both indices showed the same pattern over time: Mean abundance-weighted NRI and NTI values were highest in the 2000 census, decreased from 2000 to 2004, and then increased between 2004 and 2008, as seen in Figure 1a,b. Thus, over the study period, individuals within seedling plots tended to be most closely related in the first census following hurricane disturbance and least related ~five years after the hurricane.
However, when palm species were dropped from the data set, this pattern changed substantially. Specifically, values of both NRI and NTI tended to increase slightly over time when palms were excluded, up until the final census when values dropped significantly, shown in Figure  1c,d.

Phylogenetic Community Structure and Land-Use History
Mean values of NRI and NTI within seedling plots differed significantly between the high-intensity and low-intensity land-use categories for two of the five census years (2004 and 2007; Figure 2). In those However, when palm species were dropped from the data set, this pattern changed substantially. Specifically, values of both NRI and NTI tended to increase slightly over time when palms were excluded, up until the final census when values dropped significantly, shown in Figure 1c,d.

Phylogenetic Community Structure and Land-Use History
Mean values of NRI and NTI within seedling plots differed significantly between the high-intensity and low-intensity land-use categories for two of the five census years (2004 and 2007; Figure 2). In those years, both indices indicated that seedlings were more closely related to each other in plots in the low-intensity land-use history area than in the area that had experienced high-intensity land-use prior to 1934. When palm species were excluded, significant differences between land-use history categories in NRI and NTI values were detected not only in 2004 and 2007, but also in 2002, as seen in Figure S2, again with higher phylogenetic clustering in the low-intensity land-use history area compared to the high-intensity land-use area. years, both indices indicated that seedlings were more closely related to each other in plots in the low-intensity land-use history area than in the area that had experienced high-intensity land-use prior to 1934. When palm species were excluded, significant differences between land-use history categories in NRI and NTI values were detected not only in 2004 and 2007, but also in 2002, as seen in Figure S2, again with higher phylogenetic clustering in the low-intensity land-use history area compared to the high-intensity land-use area.

Phylogenetic Neighborhood Effects on Seedling Survival
We found a significant effect of the density and phylogenetic relatedness of neighbors on seedling survival over each of the first three census intervals (p < 0.05 for 2000-2002, 2002-2004, 2004-2007), but not in the final census interval (2007-2008, Table 1). In the first three censuses, the probability of survival declined with increasing density and relatedness of neighboring seedlings, as seen in Figure 3). For the 2000-2002 interval, the model that included density and phylogenetic relatedness of all neighbors performed better at predicting seedling survival than did the model that only included conspecific neighbor density, shown in Table 1, even when palm species were excluded from the analyses, shown in Table S1. This indicates that seedlings were negatively affected not only by seedlings of their own species but also by seedlings of closely related species that were present in the same seedling plot. AIC values were similar for both the phylogenetic and conspecific models in  Table 1.

Phylogenetic Neighborhood Effects on Seedling Survival
We found a significant effect of the density and phylogenetic relatedness of neighbors on seedling survival over each of the first three census intervals (p < 0.  Table 1). In the first three censuses, the probability of survival declined with increasing density and relatedness of neighboring seedlings, as seen in Figure 3). For the 2000-2002 interval, the model that included density and phylogenetic relatedness of all neighbors performed better at predicting seedling survival than did the model that only included conspecific neighbor density, shown in Table 1, even when palm species were excluded from the analyses, shown in Table S1. This indicates that seedlings were negatively affected not only by seedlings of their own species but also by seedlings of closely related species that were present in the same seedling plot. AIC values were similar for both the phylogenetic and conspecific models in  Table 1.

Changes in Seed Mass Values Following Hurricane Disturbance
In all years, the seedling community was dominated by individuals of relatively large-seeded species, as demonstrated by the left-skewed distributions of seed sizes shown in . This suggests that the heterogeneous conditions immediately following the hurricane permitted the recruitment of species with a wide variety of life history strategies, but the effect did not persist. This pattern was largely a result of an elevated number of seedlings of species with low seed mass in the 2000 seedling census, as seen in Figure 4.
Within seedling plots, patterns of structure in seed mass values varied significantly among census years for NRI (Kruskal-Wallis chi-squared = 56.7, df = 4, p < 0.001), but not NTI (Kruskal-Wallis chi-squared = 7.81, df = 4, p = 0.10). Values of NRI and NTI for seed mass steadily increased over each census interval, as seen in Figure 5, indicating that seedling plots were increasingly occupied by seedlings of species with more similar seed mass values. In all censuses, seedlings had more similar values of seed mass than expected (compared to the null model), as indicated by positive NRI and NTI values shown in Figure 5. Patterns of NRI and NTI for seed mass were largely consistent when palm species were dropped from the analysis, as seen in Figure S3.

Changes in Seed Mass Values Following Hurricane Disturbance
In all years, the seedling community was dominated by individuals of relatively large-seeded species, as demonstrated by the left-skewed distributions of seed sizes shown in This suggests that the heterogeneous conditions immediately following the hurricane permitted the recruitment of species with a wide variety of life history strategies, but the effect did not persist. This pattern was largely a result of an elevated number of seedlings of species with low seed mass in the 2000 seedling census, as seen in Figure 4.     Within seedling plots, patterns of structure in seed mass values varied significantly among census years for NRI (Kruskal-Wallis chi-squared = 56.7, df = 4, p < 0.001), but not NTI (Kruskal-Wallis chi-squared = 7.81, df = 4, p = 0.10). Values of NRI and NTI for seed mass steadily increased over each census interval, as seen in Figure 5, indicating that seedling plots were increasingly occupied by seedlings of species with more similar seed mass values. In all censuses, seedlings had more similar values of seed mass than expected (compared to the null model), as indicated by positive NRI and NTI values shown in Figure 5. Patterns of NRI and NTI for seed mass were largely consistent when palm species were dropped from the analysis, as seen in Figure S3.

Discussion
We used a well-resolved community phylogeny coupled with long-term seedling dynamics data to improve our understanding of forest regeneration processes following disturbance. We found evidence that the phylogenetic structure of the seedling community shifted over time as the forest recovered from hurricane disturbance. Our analyses revealed that the direction and pattern of these changes in phylogenetic structure depended strongly on whether palm species at the site were included in the analysis. We also found evidence that changes in phylogenetic community structure were shaped by both the abiotic environment and biotic interactions in the seedling community. Specifically, seedlings of small-seeded, light-demanding species, which had been able to recruit under the high light conditions created by hurricane disturbance, decreased in abundance over time as the canopy recovered and light levels dropped. At the same time, the high densities of seedlings present immediately post-hurricane led to significant conspecific and phylogenetic density-dependent mortality, reducing phylogenetic relatedness of co-occurring seedlings. In addition, our study revealed that the signature of past anthropogenic impacts were still present in the phylogenetic structure of the seedling community, even decades after human land use ceased.

Phylogenetic Community Dynamics Following Hurricane Disturbance
Phylogenetic relatedness in the seedling community was highest in the first census following Hurricane Georges. Following that census, seedling phylogenetic relatedness declined, reaching a minimum~five years post-hurricane, before increasing again. Comparing patterns of phylogenetic community structure with and without palm species (Figure 1a,b vs. Figure 1c,d) revealed that the observed changes in phylogenetic relatedness were driven by changes in the abundance of palm species. In particular, the palm Prestoea acuminata doubled in relative abundance in the seedling plots between 2000 and 2002, becoming the most abundant species in the seedling layer. Because palms are distantly related to most other species in the Luquillo forest, this increase in palm seedling abundance caused a sharp decrease in mean phylogenetic relatedness in the seedling community, particularly across the first census interval. Although Prestoea acuminata clearly recruits well under the high light conditions created by hurricane damage, it is a relatively large seeded species. This explains why patterns of phylogenetic community structure did not match patterns of mean seed mass over time when palms were included in the analyses.
When palms were excluded, phylogenetic relatedness tended to increase over time, before dropping significantly in the final census interval (Figure 1c,d). This increase in phylogenetic relatedness in the years following the hurricane was likely driven by decreased abundance of light-demanding species that had recruited in response to higher understory light resulting from hurricane damage to the canopy. Understory light levels were nearly four times higher in the Luquillo forest six months after Hurricane Georges compared to during seedling censuses 3-9 years later [29]. Small-seeded, light-demanding pioneer species responded strongly to the high light conditions found immediately after the hurricane, recruiting in higher densities than in other years. For example, seedlings of the pioneer species Cecropia schreberiana were only recorded in the seedling plots in the initial census in 2000. Other small-seeded, light-demanding species, such as Psychotria brachiata and Schefflera morototoni, were relatively common in the 2000 seedling census but had declined sharply in relative abundance by the subsequent census. Although we lack data on seedling abundances prior to the hurricane, such a high abundance of small-seeded, light-demanding seedlings were not present beneath the largely intact forest canopy that existed before the hurricane (JT, personal observation). As seedlings of these light-demanding species either grew quickly into larger size classes or died out as light levels dropped due to canopy closure, the relative abundance of more shade-tolerant species increased. For example, several large-seeded, late-successional species, such as Dacryodes excelsa and Tetragastris balsamifera, had substantially higher seedling densities in later censuses compared to the 2000 census. This shift towards larger-seeded, shade-tolerant species can be seen in our community-level analysis of mean seed mass values, which revealed higher abundance of seedlings of small-seeded species in the initial census compared to later censuses (Figure 4). At the scale of the seedling plots, this led to higher similarity of seed mass values over time following the hurricane, at least for NRI (Figure 5a). Similarly, other studies in Puerto Rico showed consistent increases in functional clustering [60,61]. This increasing functional clustering in the seedling plots appears to contribute to the post-hurricane dynamics of phylogenetic community structure when palms are excluded (i.e., similar temporal patterns of phylogenetic structure in Figure 1c,d and of seed mass structure in Figure S3a). While we attribute these changes to the effects of increased canopy openness due to hurricane damage, we also acknowledge that temporal changes in phylogenetic structure could be influenced by additional factors, such as inter-annual climatic variation [62].
Our survival model results suggest that species interactions are also shaping the phylogenetic structure of the seedling community. Specifically, we found evidence of phylogenetic density-dependent seedling mortality, particularly in the first census interval following Hurricane Georges when seedling densities were elevated due to increased canopy openness [29]. Specifically, mortality was higher for seedlings surrounding by seedling neighbors of their own or closely related species (versus distantly related neighbors). The negative effect of closely related neighbors most likely results from the fact that closely related species tend to have similar resource requirements or be susceptible to the same natural enemies [32], leading to strong resource competition or enemy-mediated apparent competition among closely related plant species [63] (reviewed by [64]). Similar negative phylogenetic neighborhood effects have been reported in other tropical forests. For example, similar results were reported by Webb et al. [34], Metz et al. [35] and Paine et al. [36] for seedlings in Borneo, Ecuador, and French Guiana, respectively, supporting the idea that the phylogenetic composition of local neighborhoods influences regeneration dynamics in tropical forests (but see [65]). However, our results contrast with those of a study of tree growth and survival, also conducted in the LFDP, which found little evidence that the phylogenetic relatedness of neighbors had a strong impact on performance of saplings and trees ≥1 cm dbh [22]. Thus, phylogenetic neighborhood effects appear to be strongest at earlier life stages in the Luquillo forest. Interestingly, this result contrasts with findings from the tropical forest of Barro Colorado Island, Panama (which does not experience hurricane disturbance), where negative effects of neighborhood phylogenetic relatedness on survival were found for adult trees, but not seedlings [66]. In the present study, we detected strong negative phylogenetic neighborhood effects on seedling survival only in the first few years following the hurricane, indicating temporal variation in the strength of such effects. Thus, our study contributes to a growing consensus that patterns of phylogenetic density dependence can vary both within sites (i.e., with life stage and temporally), as well as across sites (e.g., [67]).

Signature of Past Land-Use on Phylogenetic Community Structure
In addition to impacts of natural disturbance on phylogenetic community structure, we also found that the phylogenetic structure of the seedling community differed significantly between areas that had experienced different intensities of past human disturbance. Specifically, seedling communities in areas of high-intensity past land use were similar to or less phylogenetically clustered than those in areas of lower-intensity past land use (both with and without palm species included). These findings are opposite to findings of previous studies examining the impacts of past human disturbance on plant communities (e.g., [13,15,18,68,69]), which generally found that human disturbance tended to homogenize not only the phylogenetic composition but also the functional and taxonomic composition of tree communities [15].
The discrepancy between results of previous studies and our present study may be due to the interaction between land-use history and hurricane disturbance at our study site, shaping patterns of both disturbance and recovery. Earlier studies of the Luquillo forest have shown that the two areas differing in land use history in the LFDP respond distinctly to hurricanes [70][71][72]. Specifically, anthropogenic disturbance increases the proportion of secondary species, which are often more vulnerable to storm damage [73]. Thus the high-intensity past land use area of the plot experienced more canopy damage in Hurricane Georges than the low-intensity area, leading to higher and more variable understory light levels in the high-intensity area [74]. These differences likely influenced recruitment patterns and phylogenetic structure after the hurricane. Indeed, while the two areas showed qualitatively similar patterns of seedling phylogenetic structure over time, the magnitude of temporal change was different ( Figure 2). The greater hurricane damage in the high-intensity land use area may have sustained a mix of early successional and late successional species in the seedling community in the years following hurricane disturbance, which resulted in the higher observed phylogenetic overdispersion than in the low-intensity land use portion of the LFDP, which suffered less damage and thus had a more homogenous understory.

Conclusions
Our results reveal how disturbance can impact both the abiotic and biotic environment, causing changes in phylogenetic and functional community structure over space and time. Our study focuses on a single life stage but suggests that both human and natural disturbance influence the phylogenetic structure of ecological communities and that their effects are interactive. Understanding the impacts of disturbance on forest communities, as well as the underlying processes driving forest recovery, is critical given that the frequency of some types of severe natural disturbances (e.g., hurricanes, droughts) is predicted to increase due to climate change [75,76]. In addition, an increasing proportion of forested area is made up of secondary forest recovering from past human disturbance [77]. Understanding the implications of past human land use for phylogenetic diversity (and thus evolutionary history and potential), as well as how human disturbance interacts with natural disturbances, should be a key research priority. Recent disturbance events, such as Hurricanes Irma and Maria in the Caribbean [78] and droughts and fires in the Neotropics and Southeast Asia due to the severe 2015-2016 El Niño event [79], offer the opportunity for studies that will increase our understanding of forest disturbance and recovery in a changing world.
Supplementary Materials: The following are available online at http://www.mdpi.com/1999-4907/9/9/556/s1, Figure S1: DNA barcode phylogeny for the Luquillo Forest Dynamics Plot; Figure S2: Mean phylogenetic relatedness of seedlings after excluding palm species; Figure S3: Mean similarity of seed mass values of seedlings after excluding palm species; Table S1: AIC values for neighborhood model comparisons with palm species excluded.