Functional Ecology of Forest, Heath, and Shrub Savannah Alternate States in Eastern Canada

: In eastern Canada, alternation of wildﬁre regime due to ﬁre suppression creates alternate vegetation states converting black spruce forest to heath and shrub savannah (SS). We compared the taxonomic diversity (TD) and functional diversity (FD) of post-ﬁre forest, heath, and SS alternate states to determine if community FD can explain their persistence. We hypothesized that (i) species diversity (TD and FD) would be the highest in forest followed by SS and heath due to decreased interspeciﬁc competition and niche differentiation, (ii) differences between TD and FD indices would be greater in communities with high TD in forest due to high trait differentiation and richness, and (iii) changes in community trait values would indicate niche limitations and resource availability. We conducted this study in Terra Nova National Park, Newfoundland, Canada. We calculated functional dispersion (alpha FD), functional pairwise dissimilarity (beta FD), Shannon’s diversity (alpha TD), and Bray–Curtis dissimilarity (beta TD) from species cover. We used ﬁve functional traits (speciﬁc root length, speciﬁc leaf area, leaf dry matter content, height, and seed mass) related to nutrient acquisition, productivity, and growth. We found lower beta diversity in forest than heath and SS; forest also had higher species diversity and greater breadth in niche space utilization. SS was functionally similar to heath but lower than forest in functional dispersion and functional divergence. It had the highest functional richness and evenness. There was no difference in functional evenness between forest and heath. Functional beta diversity was the highest in forest, and did not differ between heath and SS. Resource acquisition and availability was the greatest in forest and the lowest in heath. We suspect that this might be due to forest having the highest functional trait turnover and niche utilization. We conclude that alternate vegetation states originating from alterations to the natural ﬁre regime negatively impact ecosystem function.


Introduction
Alterations to the natural disturbance regime of a region can compromise ecosystem resilience and create alternate vegetation sates [1,2]. Post-fire alternate vegetation states may have distinct vegetation structures and functions that differ from the pre-fire community [3]. The severity and frequency of forest fires have powerful effects on the colonizing species and alter ecosystem dynamics [4]. Persistence of these effects may result in communities that support less species diversity, with negative implications for community stability and productivity [5]. Persistence of alternate vegetation states is often facilitated by environmental filters that interact with the species present and influence their establishment [6]. Strong abiotic filters in a system may select for specific "response traits" that flourish in a narrow range of conditions [7], or strong biotic filters may prevent the establishment of competing species in a community, creating what is often called "effect traits" [8]. Biotic and abiotic filters produce stabilizing effects and reinforce the alternate state [9].
Taxonomic diversity has been shown to predict resource availability to plants [10]. However, measuring taxonomic diversity alone ignores the differences in niche utilization.

Study Site
The study was conducted in Terra Nova National Park, Newfoundland, Canada (48 • 31 33 N, 53 • 57 58.4 W). Although this region is in the boreal ecozone, the topography and surrounding cold ocean causes irregular patterns of vegetation compared to the boreal region of the adjacent land to the west [20]. The region receives 1801-1200 mm precipitation (~30% as snow) and has mean summer and winter temperatures of 13 and −5 • C, respectively [21]. The region has three distinct terrestrial community states: forest, heath, and shrub savannah. Forest has high stem densities of trees forming a predominantly closed canopy consisting of Abies balsamea, Acer rubrum, Betula papyrifera, Larix laricina, Picea mariana, Sorbus americana, Prunus spp., and Populus tremuloides. This canopy is underlain by a variety of ericaceous shrubs such as Kalmia angustifolia, Rhododendron groenlandicum, Vaccinium angustifolium, V. myrtilloides, V. uliginosum, and Rubus idaeus and herbs such as Maianthemum canadense, Cornus canadensis, Clintonia borealis, and Carex spp., with acrocarpous and pleurocarpous mosses on the ground level. Heath is predominantly open canopy communities dominated by ericaceous species such as Kalmia angustifolia, Rhododendron groenlandicum, Vaccinium angustifolium, and Rhododendron canadense. Other herbs and shrubs, such as Ilex mucronata, Cornus canadensis, Amelanchier interior, Vaccinium vitis-idaea, and Viburnum cassinoides, are common but in low abundances. The ground is sparsely inhabited by lichens (Cladonia uncialis, C. mitis) and Dicranum spp. and Polytrichum spp. mosses. The shrub savannah has two distinct components. The interspersed tree islands consist primarily of Picea mariana with discrete Larix laricina in low abundances. Within these islands, there are pleurocarpous mosses and a few herbs such as Cornus canadensis, Trientalis borealis, and Maianthemum canadense. The second, and more widespread, open canopy component of shrub savannah is dominated by Kalmia angustifolia, Rhododendron canadense, Vaccinium angustifolium, and Cladonia stellaris. Three sites were selected in each post-fire community based on a similar amount of time since the last fire disturbance (heath: 24-30 years since fire; forest: 12-44 years since fire) [22], with the exception of shrub savannah, because the specific time since fire was unknown but presumed to be >100 years based on several tree cores ( Figure 1).

Study Design and Sampling Protocol
Three 50 × 10 m belt transects were placed in each of the sites. The starting positions of the belt transects were determined by randomly selecting coordinates then placing the belt transects along a north-south gradient. Each 50 × 10 m transect consisted of five continuous 10 × 10 m quadrats. Additionally, five 1 × 1 m quadrats were nested within each 10 × 10 m quadrat. Four of the five 1 × 1 m quadrats were placed 0.5 m from the 10 × 10 m quadrat edges, and one 1 × 1 m quadrat was placed in the center of each 10 × 10 m quadrat. Tree abundance was determined within the 10 × 10 m quadrats and represented by basal area percent coverage. The percent coverage of tree basal area for each species was determined by taking a total area occupied by a species basal area divided by the total area of the plot. Herb and shrub percent cover was collected by ocular assessment within the 1 × 1 m quadrats and averaged among the five nested quadrats within the 10 × 10 m quadrats. The sampling method was repeated in three replicates in each site and three sites per plant community of each vegetation state (forest, heath, shrub savannah). One of the forest sites was restricted to two 50 × 10 m belt transects due to the unavailability of a third appropriate study site (Table 1).

Figure 1.
Location of study sites in the forest (F, in blue), heath (H, in red), and wood savannah (WS, in green) communities in Newfoundland, Canada [23]. Table 1 of Siegwart Collier and Mallik 2010 [15] shows the fire history (date of fire, fire weather index, duration of fire, area burned etc.) of the study sites.

Plant Trait Matrix
The plant trait matrix was based on the 25 most abundant species consisting of >90% of the total cover on the sites. Functional traits were determined for each of these 25 species. The traits considered were specific leaf area (SLA), leaf dry matter content (LDMC), specific root length (SRL), seed size, and height (Appendix A, Table A1). SLA strongly and positively correlates with foliar nutrient concentration, particularly with N and P [24]. LDMC represents aboveground productivity and specifically indicates the conversion rates of CO2 and H2O to dry matter [25]. LDMC also relates to resource availability and growth rate. Belowground, SRL is strongly correlated to soil fertility, meaning that species with higher SLA are adapted to habitats with fertile soils [26]. Height relates to light capture in plants, and it is also strongly positively correlated with seed dispersal [27]. Seed mass is a particularly important plant characteristic. Increasing seed mass relates to an increased germination rate, increased longevity, and less seed production [28].

Plant Trait Matrix
The plant trait matrix was based on the 25 most abundant species consisting of >90% of the total cover on the sites. Functional traits were determined for each of these 25 species. The traits considered were specific leaf area (SLA), leaf dry matter content (LDMC), specific root length (SRL), seed size, and height (Appendix A, Table A1). SLA strongly and positively correlates with foliar nutrient concentration, particularly with N and P [24]. LDMC represents aboveground productivity and specifically indicates the conversion rates of CO 2 and H 2 O to dry matter [25]. LDMC also relates to resource availability and growth rate. Belowground, SRL is strongly correlated to soil fertility, meaning that species with higher SLA are adapted to habitats with fertile soils [26]. Height relates to light capture in plants, and it is also strongly positively correlated with seed dispersal [27]. Seed mass is a particularly important plant characteristic. Increasing seed mass relates to an increased germination rate, increased longevity, and less seed production [28].

Data Analysis
Alpha and beta functional diversity were represented by functional dispersion (FDis) and functional pairwise dissimilarity (Dpw), respectively. Alpha and beta taxonomic diversity were represented by Shannon's diversity and Bray-Curtis dissimilarity, respectively. Additionally, community-weighted mean trait values were determined in order to quantify the relative dominance of species throughout communities and average trait values favored within communities. Correlations between alpha and beta diversity measures were computed with the 'cor.test' function of the R Stats Package [29]. The 'dbFD' function of the FD package [30] in R was used to calculate FDis and community-weighted mean. Functional pairwise dissimilarity was calculated using the 'mpd' function of the 'picante' package [31]. One discrepancy was identified between mathematical reality and ecological reality when calculating the functional space [15]. In situations where the number of species present in the functional trait matrix is less than the number of functional traits, the equations cannot yield a value for the five-dimensional space, despite having one or several species present. To calculate the functional composition of these anomalies, one standard deviation of the mean of the site was added to the non-applicable values in order to avoid representing low diversity values as zero diversity, as per [18]. Meanwhile, sites not containing any of the species in the functional trait matrix were represented as zero to maintain orthogonality in design and represent not only the lack of functional contribution but also the lack of functional diversity within certain communities. Shannon's diversity was calculated using the 'diversity' function of the vegan package [32], and Bray-Curtis dissimilarity was calculated using the 'vegdist' function of the vegan package [33]. Comparisons of alpha and beta diversity between sites were conducted using a Dunn test of multiple comparisons with a Bonferroni approach with the 'dunn.test' r-package [34]. All statistical analyses were completed in R version 1.1.463 [29].

Correlations between Functional and Taxonomic Diversity
Functional and taxonomic alpha diversity showed a strong correlation between functional dispersion and Shannon's diversity across all sites and life forms (Pearson's correlation, r = 0.872, p < 0.001) (Figure 2a). Correlations between functional pairwise dissimilarity and Bray-Curtis taxonomic dissimilarity showed a weak and non-significant relationship between functional and taxonomic beta diversity (Pearson's correlation, r = −0.117, p = 0.187) (Figure 2b).

Data Analysis
Alpha and beta functional diversity were represented by functional dispersion (FDis) and functional pairwise dissimilarity (Dpw), respectively. Alpha and beta taxonomic diversity were represented by Shannon's diversity and Bray-Curtis dissimilarity, respectively. Additionally, community-weighted mean trait values were determined in order to quantify the relative dominance of species throughout communities and average trait values favored within communities. Correlations between alpha and beta diversity measures were computed with the 'cor.test' function of the R Stats Package [29]. The 'dbFD' function of the FD package [30] in R was used to calculate FDis and community-weighted mean. Functional pairwise dissimilarity was calculated using the 'mpd' function of the 'picante' package [31]. One discrepancy was identified between mathematical reality and ecological reality when calculating the functional space [15]. In situations where the number of species present in the functional trait matrix is less than the number of functional traits, the equations cannot yield a value for the five-dimensional space, despite having one or several species present. To calculate the functional composition of these anomalies, one standard deviation of the mean of the site was added to the non-applicable values in order to avoid representing low diversity values as zero diversity, as per [18]. Meanwhile, sites not containing any of the species in the functional trait matrix were represented as zero to maintain orthogonality in design and represent not only the lack of functional contribution but also the lack of functional diversity within certain communities. Shannon's diversity was calculated using the 'diversity' function of the vegan package [32], and Bray-Curtis dissimilarity was calculated using the 'vegdist' function of the vegan package [33]. Comparisons of alpha and beta diversity between sites were conducted using a Dunn test of multiple comparisons with a Bonferroni approach with the 'dunn.test' r-package [34]. All statistical analyses were completed in R version 1.1.463 [29].

Overall, Tree, Shrub, and Herb Diversity
Forest displayed the highest diversity for all indices except beta TD when all life forms were analyzed. Overall alpha FD was significantly greater in forest than heath and shrub savannah (Kruskal-Wallis (KW), x = 8.096, p < 0.001; KW, x = 6.885, p < 0.001, respectively), while heath and shrub savannah did not differ (KW, x = −1.248, p = 0.318) (Figure 3a). A similar trend was observed within alpha TD, with forest surpassing both heath (KW, x = 8.014, p < 0.001) and shrub savannah (KW, x = 6.397, p < 0.001). No differences were observed between heath and shrub savannah (KW, x = −1.666, p = 0.144) (Figure 3b). Overall beta FD and beta TD showed different hierarchical trends among the three vegetation states and differed from each of the alpha diversity analogues. Similar to overall alpha FD, forest had higher overall beta FD than heath (KW, x = 7.931, p < 0.001) and shrub savannah (KW, x = 5.727, p < 0.001). However, unlike alpha FD, shrub savannah had significantly higher beta FD than heath (KW, x = −2.272, p = 0.035) ( Figure 3c). Heath had greater overall beta TD than forest (KW, x = −3.581, p < 0.001) and shrub savannah (KW, x = 2.204, p = 0.041), followed by forest, which had higher overall beta TD than shrub savannah (KW, x = 5.963, p < 0.001) (Figure 3d).

Discussion
Our results show that functional and taxonomic diversity were almost consistently higher in post-fire forest than heath and shrub savannah. shrub savannah showed an intermediate level of diversity and retained similar average functional trait values compared to heath, with the exception of height and SRL. Despite strong correlations between overall alpha TD and alpha FD, we found discrepancies in the hierarchical ranking of forest, heath, and shrub savannah between the two measures of alpha diversity within tree and herb life forms. Here, differences between shrub savannah and heath were found solely for alpha FD. In the case of trees, the amplified differences in alpha FD compared to alpha TD are likely a result of greater functional trait differences between the tree species present [33]. One possible reason for the decreased tree diversity in heath than shrub savannah is high biotic filtering, where Kalmia dominance interferes with black spruce regeneration [18,34]. This explanation also provides a rationale for alpha FD capturing more differences than the taxonomic analogue, as it accurately describes the reduction of successful trait ranges that occurs via habitat filtering [35]. The herb subgroup showed more differentiation between heath and shrub savannah due to increased microhabitat conditions captured within the spruce islands, which has been shown to favor biodiversity due to fine-scale spatial variation [36]. However, differences in shrub and overall alpha FD and alpha TD appear equal, suggesting that the spruce islands and their associated microhabitats contribute minimally to community-wide biodiversity and that shrubs control the functional processes in these communities.
The differences between the hierarchical ranking of the three vegetation states between beta TD and beta FD are not surprising, due to the weak correlation observed among them. Within heath, beta TD suggests high levels of heterogeneity. However, beta FD shows that the dissimilarity of plant species in the heath varies minimally. Examining beta diversity from a functional perspective indicates that a narrow range of functional traits are capable of establishing themselves in the environmental constraints of heath, which perpetuates the system in a state of homogenization [37]. Consequently, we see that functional homogenization is the highest in heath, followed by shrub savannah then forest. The tree and shrub beta FD do not differ because heath and shrub savannah select for only two conifer species (Picea mariana, Larix laricina). These two conifer species are able to persist in heath and shrub savannah by avoiding the inhospitable seedbed through vegetative reproduction. It appears that the dispersal of herbs, primarily in the spruce islands, is the primary driver of beta FD in shrub savannah (Figure 6c). In contrast, beta FD is the highest in forest due to the divergence of traits into a wider variety of niches. As a result, the biodiversity and consequent increase in competitive interactions between plant species in forest creates a positive feedback loop that results in increased niche differentiation and productivity within the community [38]. Therefore, forests have the highest degree of niche variability, followed by shrub savannah and heath.
It is difficult to include the traits of all plants (including rare species) to calculate FD. However, studies have shown that patterns of FD are often driven by the traits of the dominant species according to the mass ratio hypothesis [13,39,40]. Nonetheless, we recognize that the omission of rare vascular plants and cryptogams in the calculation of FD may have resulted in the loss of some information [18]. Trait analysis of cryptogams is challenging and is still in its infancy [41].
Although levels of niche differentiation vary across the three vegetation states (forest, heath, and shrub savannah) each remains favorable to specific functional traits adapted to the biotic and abiotic conditions in each community. Average SLA, LDMC, and seed mass are similar between heath and shrub savannah. The relatively high LDMC and low SLA observed in heath and shrub savannah represents nutrient conservation, whereas the relatively high SLA and low LDMC in forest is more representative of the rapid production of biomass [42]. This is consistent with pyrogenic organic matter resulting from incomplete combustion of organic matter in heath and wood savanna leading to lower available N, which subsequently lowers the C:N ratio and immobilizes N [43]. Charcoal (pyrogenic carbon) from forest fires has been shown to change soil pH [44] and produce a significant long-term carbon sink [45,46]. Moreover, past research [47] has shown that the abundance of ericaceous plants, such as Kalmia, is negatively correlated with soil ammonium and phosphorus contents, leading to limitations in soil nutrient availability. Additionally, the phenolic compounds produced from the highly abundant Kalmia [48] are known to inhibit nitrogen mineralization. This interferes with nutrient cycling and the establishment of non-ericaceous species. Additionally, thick organic matter occupied by Kalmia creates an inhospitable seedbed for conifers but a suitable condition for vegetative regeneration [49]. Hence, heath and shrub savannah have low seed mass values, as primarily seed regenerating species (black spruce) are outcompeted by vegetatively regenerating ericaceous species. Not surprisingly, SRL differs between each of the three communities and is the greatest in the heath, followed by shrub savannah, then forest. SRL decreases with soil fertility [25]. Therefore, it is likely that the spruce islands present in shrub savannah may contribute to soil fertility, either through nitrogen mineralization via mycorrhiza or simply through limiting the abundance of Kalmia. Plant height was significantly greater in forest and shrub savannah compared to heath. shrub savannah surpassed heath in plant height due to the prevalence of layered spruce islands that evade the hostile seedbed and promote the competitive pressures for light acquisition, resulting in species that are adapted to shaded environments, such as Cornus canadensis, and species that are adapted to acquire light through a rapid increase in height, such as Larix laricina. Accordingly, the three communities represent differentiation of traits to better suit their respective environmental conditions.
The changes in taxonomic and functional biodiversity across the three alternate vegetation states provide insight into ecosystem health and function. From a taxonomic perspective, the increase in alpha and beta diversity observed in forest shows that the biodiversity of the system is increased, which predicts that the system may also have a higher degree of functioning and stability than heath and shrub savannah [50]. Yet, beyond that prediction, taxonomic diversity provides little information regarding the functional role of species in the community [51]. Hence, alpha TD found minimal differences between the herb and tree life forms that propelled alpha FD differences between heath and shrub savannah. Alternatively, this suggests that both a lack of fire severity and/or the absence of fire results in a system with less resilience to environmental change and more susceptibility to regime shifts [52], meaning that long-term alteration to the natural fire regime of eastern Canada can result in fragmented assemblages of alternate states.
Functionally, the consequence of reduced black spruce forests has serious implications for productivity and ecosystem function. The loss of functional diversity reduces the level of niche complementarity and variability in resource acquisition strategies [53]. One of the driving factors for loss of function is the biological homogenization observed in heath and shrub savannah. Biological homogenization at the community level alters ecosystem productivity and function through the loss of specialists in a community [54]. This is supported by the increased alpha diversity in this study showing high degrees of coexistence among specialists in the forest.
The functional and taxonomic diversity derived from this study are consistently lower than those reported by other research [18] that compared the functional and taxonomic diversity differences in heath and forest. The differences occurred as a result of methodological differences. The aforementioned study used random quadrats and ocular assessment of tree crown cover, whereas in this study, we used belt transects and basal area of trees to assess tree species cover. The differences of tree diversity measures arise from the relative differences of abundance values generated by crown cover and basal area of trees, where crown cover is far greater than basal area. However, assessing crown cover provides a less precise measure of cover and thus provides the potential to misrepresent the influence of trees on the system's biodiversity. Alternatively, measuring the basal area of tree species may underrepresent their contribution to biodiversity when used in conjunction with the ocular assessment of ground cover by other species. While methods for conversion from basal area to crown diameter exist, they are inapplicable to immature trees. Thus, there is a need to refine these conversion methods to include immature trees (saplings) in order to properly assess tree species diversity in recently disturbed sites. Independent of methodological differences between tree assessments, herb and shrub diversity values remain lesser when using belt transects. The use of belt transects are superior to quadrats when assessing environmental gradients [55]. Therefore, it is likely that the values derived via the belt transect method more accurately capture the functional characteristics of the system than random quadrats.

Conclusions
Forest, heath, and shrub savannah of eastern Canada provide the opportunity to explore the consequences of alternate states that arise from alterations to the natural disturbance regime of an ecosystem. Where the black spruce forest is considered the historic community, deviations from the natural disturbance regime shifted the community to the alternate states of heath and shrub savannah, which resulted in the loss of alpha FD and alpha TD as well as beta FD. Furthermore, the dominant trait values shifted with the creation of alternate states in response to changes in abiotic environmental filters, such as increased post-fire residual organic matter, and biotic filters, such as Kalmia dominance. The systems with less biodiversity are associated with less nutrient availability and lower levels of niche differentiation. This results in a consequent reduction of the establishment potential of species, which dampens the overall ecosystem function.