Habitation Sites Influence Tree Community Assemblages in the Great Bear Rainforest, British Columbia, Canada

Schang, Kyle; Cox, Kieran; Trant, Andrew J.

doi:10.3389/fevo.2021.791047

ORIGINAL RESEARCH article

Front. Ecol. Evol., 04 January 2022

Sec. Population, Community, and Ecosystem Dynamics

Volume 9 - 2021 | https://doi.org/10.3389/fevo.2021.791047

Habitation Sites Influence Tree Community Assemblages in the Great Bear Rainforest, British Columbia, Canada

$\r\nKyle Schang,*$ Kyle Schang^1,2*

Kieran Cox^2,3 $Andrew J. Trant,\r\n$ Andrew J. Trant^1,2

¹School of Environment, Resources and Sustainability, University of Waterloo, Waterloo, ON, Canada
²Hakai Institute, Calvert Island, BC, Canada
³Department of Biology, University of Victoria, Victoria, BC, Canada

Identifying how past human actions have influenced their environment is essential for understanding the ecological factors that structure contemporary ecosystems. Intertidal resource use by Indigenous Peoples for thousands of years has led to habitation sites containing vast shell midden deposits and facilitating long-term impacts on soil chemistry and drainage. Here we examine how these shell middens have impacted various forest metrics, such as species diversity, community composition, canopy height, and regeneration recruitment to determine if forests on habitation sites differ from the surrounding matrix. We surveyed known habitation sites with archeological evidence indicating past year-round human occupation, within the Hakai Lúxvbálís Conservancy on Calvert and Hecate Islands within the Great Bear Rainforest along British Columbia’s Central Coast. Our results demonstrate that habitation sites exhibit lower tree species richness, less relative species abundances, as such, displayed lower Shannon diversity and inverse Simpson values. The composition of tree communities on habitation sites was statistically different, with western hemlock and western redcedar densities increasing on non-habitation sites. Conversely, regeneration diversity at habitation sites was more even and exhibited elevated Shannon diversity and inverse Simpson values. The community composition of regeneration was more consistent among habitation and non-habitation sites; however, western redcedar, western hemlock and Sitka spruce were more abundant at habitation sites. For all tree species, maximum height was higher within the habitation sites; however, this trend was the most notable in western redcedar and Sitka spruce, which increased by an average of 4.8 m relative to non-habitation sites. Collectively, our findings suggest that long-term habitation alters forest community compositions. The landscape alterations within habitation sites promote conditions needed to support diverse, even, and abundant regeneration communities and consequently increase the height of the dominant coastal tree species. Thus, our results offer evidence that long-term influence by Indigenous communities have a persistent influence on coastal forests.

Introduction

Human influence on global biodiversity and ecosystem processes have been well-documented, with post-industrial environmental degradation increasing as human populations continue to grow (Steffen et al., 2015; Watson et al., 2016). Contemporary human land-use patterns are often associated with negative ecological impacts such as habitat destruction (McIntyre and Hobbs, 1999; Hoekstra et al., 2005), soil degradation (Lamb et al., 2005; Lal, 2015), and global biodiversity loss (Butchart et al., 2010). However, there is a growing body of evidence that suggests alternative spatial-temporal models for how humans have influenced their environment. These include phenomena such as the use of fire to enhance resource availability (Bird et al., 2020), constructed clam gardens that increase shellfish productivity (Groesbeck et al., 2014; Deur et al., 2015; Smith et al., 2019), terrestrial deposition of shell middens, which enhance resource productivity (Trant et al., 2016), and the increase of food producing tree species near archeological sites in Amazonia (Levis et al., 2017). Developing a conceptualization of these ecological legacies, or effects of past human modification, is imperative for understanding current and future ecosystem functioning by providing a more accurate baseline to assess patterns and processes (Swetnam et al., 1999).

Marine-derived nutrients influence on near-shore vegetation is a well-documented phenomenon, often seen through nitrogen loading via salmon carcasses or other marine subsidies deposited within terrestrial ecosystems (Ben-David et al., 1998; Hilderbrand et al., 1999; Cox et al., 2020). Marine-derived resources can have multiple effects on riparian zones by enhancing productivity and altering plant community structure (Bartz and Naiman, 2005; Hocking and Reynolds, 2012). Similarly, anthropogenically mediated nutrient inputs, such as the accumulation of charcoal and shell middens, can alter soil chemistry through increases in nutrients such as calcium (Ca), phosphorus (P), and pH (Cook-Patton et al., 2014; Vanderplank et al., 2014; Hoffman et al., 2016a; Trant et al., 2016; Fisher et al., 2019). Within the Pacific Northwest, shell midden deposits, often found along the shorelines, can be used alongside other unique landscape features, such as clam gardens or estuarine root gardens, to locate previously occupied habitation sites – places of repeated, long-term, or high occupancy human settlement (McLaren et al., 2014; Jackley et al., 2016; Cox et al., 2019). Within the context of this study, shell middens are considered to be the accumulated remains of intertidal and subtidal taxa that were purposely deposited to alter soil characters or discarded after consumption. The deposition of shell middens by local Indigenous groups for terracing or refuse have had various documented impacts within the Pacific Northwest and can be substantial in size, exceeding 5-m in depth (McLaren et al., 2014). Shell middens in these locations reflect the past activities of occupants, being composed of resources found in the intertidal zone such as bivalve shells, fish bones and other materials such as plant remains, human remains, and rich organic soils. Their long-term residency within these soils has resulted in the deposition of large amounts of calcium (Ca) and phosphorous (P) in systems where they might not be readily abundant and have been shown to enhance forest productivity (Cook-Patton et al., 2014; Trant et al., 2016; Toniello et al., 2019). Increased levels of nutrients have varying ecological impacts depending on the ecosystem; locations such as wetlands or tundra communities have shown decreases in plant diversity (Theodose and Bowman, 1997; Bedford et al., 1999), whereas semi-arid environments can experience increases in plant diversity (Perroni-Ventura et al., 2006). These alternations to forest community composition can persist at the regional and landscape level even after intense resource-use on habitation sites has ceased (Trant et al., 2016; Fisher et al., 2019).

The Central Coast of British Columbia (BC) is a place of cultural and ecological importance to the Indigenous peoples who, in some cases, have been occupying the landscape for over 14,000 years (McLaren et al., 2018). The long-term intergenerational occupation of these lands has led to a variety of landscape level modifications. For example, the development of shell middens signifies Indigenous peoples’ dependence on intertidal resources, and stewardship of these landscapes (Sawbridge and Bell, 1972; Deur et al., 2015). While Indigenous groups still utilize the land today, they have not seen full-time occupation since the nineteenth century where widespread deaths associated with European contact and a smallpox epidemic occurred (Boyd, 1994; Jones, 2014). With the overall intensity of resource-use diminishing habitation sites, which were once mostly cleared of vegetation, have since been recolonized by western hemlock (Tsuga heterophylla) and western redcedar (Thuja plicata).

The coastal temperate rainforests of BC are some of the largest and most intact globally. Many of the habitation sites located along the Central Coast of BC are located within the hypermaritime temperate rainforests systems seen on the outer coast (Banner et al., 2005). These hypermaritime temperate rainforests systems are dominated by conifer species which see little turnover due to the lack of large-scale disturbances present (Daniels and Gray, 2006). Given the close proximity of habitation sites to the intertidal zone, many of these sites are located within the zonal forest, one of four major dominant vegetative regimes within this system (Banner et al., 2005). These zonal forests often contain an array of slopes and valleys on loose colluvial soil of various drainage regimes (Meidinger and Pojar, 1991). Of the few tree species present, western hemlock and western redcedar are the dominant species of the dynamic climax community in the zonal forests in the Great Bear Rainforest (Banner et al., 2005). This region is dominated by low productivity old-growth forests where gap dynamics, fine-scale canopy gap disturbances and fine-scale mortality events drive forest turnover (Daniels, 2003; LePage and Banner, 2014; Hoffman et al., 2021). Post-disturbance interspecies competition for these fine-scale canopy gaps are often dominated by western redcedar and western hemlock seedlings due to their efficiency in low light and nutrient limited environments (Messier, 1993; Brown et al., 1996; Lacourse, 2005; Coates et al., 2009). Local climatology and microclimates often limit recruitment as site specific topography is often complex, leading to high levels of variation between soil and air temperature within the understory (Banner et al., 2005). Coarse woody debris, henceforth referred to as nurse logs, has been documented to act as a plant nursery within the Pacific Northwest, providing seedlings with the necessary conditions to germinate without the cover of dominant understory vegetation such as Gaultheria shallon (Salal) (Harmon and Hua, 1991). The role of nurse logs has been explored in some areas of the Pacific Northwest (Graham and Cromack, 1982), yet their role in facilitation regeneration within the Great Bear Rainforest is unknown.

While studies within this region have explored forest productivity, understory composition, and stand structure/composition (see Trant et al., 2016; Fisher et al., 2019; Hoffman et al., 2021), the interactions between shell middens and forest composition/regeneration remains unexplored. In this paper, we quantify the effect of shell middens on the composition of zonal coastal forests on the Central Coast of BC by addressing the following three questions: (1) Has the presence of shell middens impacted the diversity of adult trees and regeneration? With shell middens increasing the availability of nutrients, we aim to see if their presence acts as an ecological release to the normally nutrient limited systems. (2) Do adult trees on habitations sites have an increased capacity for maximum growth? While previous studies have showed the maximum height of western redcedars were higher on habitation sites, we seek to explore if this trend is present for other species found (Trant et al., 2016). (3) Finally, do nurse logs facilitate regeneration within this system and are there differences between those on habitation to control sites? We hypothesized that long-term habitation and the associated accumulation of marine-derived nutrients transform forest communities, and ecological responses to this habitat will be evident within adult tree species and regeneration assemblages. However, the multifaceted interaction between available resources and forest community dynamics suggests that forest’s community metrics response to long-term habitation will not be consistent. Tree species, specifically western redcedar, on habitation sites are expected to have an increased capacity for maximum growth compared to those found in surrounding areas (Trant et al., 2016). Thus, we aim to see if this trend is seen in other species present on these sites. However, the implications of long-term habitation for the diversity and abundance of adult tree species are less predictable and could play a significant role in shaping future forest dynamics. Previous research on herbaceous layers and shrub diversity (Vockenhuber et al., 2011; Fisher et al., 2019) suggest that regeneration communities will likely increase in diversity and abundance, with nurse logs acting as the primary substrate for seedling germination (Harmon and Hua, 1991). Given the effects of these shell middens have on forest productivity, nurse logs present on habitation sites may better facilitate regeneration. This is, however, speculative in nature.

Materials and Methods

Study Area

All data were collected in the Hakai Lúxvbálís Conservancy on Calvert and Hecate Islands within the Great Bear Rainforest. The Great Bear Rainforest is a 6.4 million-hectare low-productivity coastal temperate (perhumid) rainforest along the coast of BC (N 51° 38′–W 128° 05′). The study sites are located within the Coastal Western Hemlock zone (CWH) which is characterized by a high annual rainfall (1,000–4,400 mm), a cool mesothermal climate and poorly drained organic soils (Gessel et al., 1973; Meidinger and Pojar, 1991; Prescott et al., 2000; Richardson et al., 2004). Within this region lies the Very Wet Hypermaritime zone (CWHvh), a sub-variant of the Coastal Western Hemlock zone (CWH) classified by its continentality and precipitation (Meidinger and Pojar, 1991; Banner et al., 2005). This zone is dominated by western redcedar and western hemlock stands that are situated on high-water tables with low nutrients and poor drainage (Messier, 1993; Banner et al., 2005). However, other co-dominant species such as yellow cedar (Cupressus nootkatensis), Sitka spruce (Picea sitchensis), Pacific yew (Taxus brevifolia), and shore pine (Pinus contorta) can also be found in these areas (Spies and Franklin, 1989; Klinka et al., 1996).

High levels of precipitation and poor soil drainage often limit the occurrence of natural understory fires, such that this system experiences few large-scale disturbances of this type (Daniels and Gray, 2006; Hoffman et al., 2016a). While lightning seldomly causes fires within the Pacific Northwest, there is evidence of controlled burning on habitation sites by Indigenous peoples that have led to an elevated level of charcoal within the naturally nutrient-poor soil (Hoffman et al., 2016b). Additionally, the Pacific Northwest often experiences high intensity windstorms via the Pacific Decadal Oscillation that can influence coastal tree growth (Knapp and Hadley, 2012).

Calvert and Hecate Islands contain four different vegetative types that vary altitudinally and are classified based on dominant species and productivity (Banner et al., 2005; Thompson et al., 2016). All habitation and control sites were established in the zonal forest, the highly productive riparian interface between the coast and the bog forest further upland. Given their close proximity to the intertidal zone, the influence of altitudinal parameters was kept relatively constant across study sites.

Field Sampling

Five habitation and four control sites (Figure 1) were chosen based on several criteria. Habitation sites were identified using a combination of archeological and ecological data that contain strong evidence of past Indigenous management and settlement. These data included soil augers containing shell middens and artifacts, presence of culturally modified trees (mostly western redcedars with persistent scars from bark stripping), presence of elevated charcoal levels in the soil, and presence of nearby clam gardens (McLaren et al., 2015; Hoffman et al., 2016a; Jackley et al., 2016; Trant et al., 2016). Control sites were chosen based on proximity and ecologically relevant parameters, such as aspect and slope, to reflect the conditions of the habitation sites while having no evidence of the previously described Indigenous management and settlement. Due to the cultural sensitivity of these lands, specific habitat site locations will not be disclosed, and instead an approximate location will be provided (Figure 1). Within each site, three 11.28 m radius plots, roughly 400 m², were established (Supplementary Figure 1), with the exception of two habitation sites that contained only two plots due to the smaller extent of the shell midden. The boundary of one plot was positioned at least 10 m from the edge of the nearest plot, and at least 10 m from the forest-intertidal edge to avoid any confounds from edge effects such as increased light penetration. Altogether, we sampled 13 plots nested within five habitation sites and 12 plots nested within four control sites. The initial plot locations at each site were chosen using stratified random subsampling with the use of GE path (v1.4.6) to superimpose a square grid onto imagery in Google Earth Pro (v.7.3.1). Grid lines were set at 2.5 m and were positioned over locations of known shell middens within habitation sites. Any 2.5 × 2.5 m quadrat that landed outside of the shell midden boundaries, or only partially covered the shell midden boundaries, were excluded from the list of possible sites. Plots were then randomly selected from all remaining quadrats. The Universal Transverse Mercator (UTM) at the center of each quadrat was extracted from Google Earth Pro. Within each quadrat, all living adult species with diameter at breast height [(1.3 m; DBH) > 7.5 cm] were recorded for height (Figure 2).

FIGURE 1

Figure 1. Map of Calvert and Hecate Islands showing habitation (red polygons) and control sites (blue hashed polygons). Red polygons are presented as larger fields as not to disclose specific habitat site locations. Site types are denoted by H, Habitation and C, Control.

FIGURE 2

Figure 2. Diversity profile curves plotting habitation (green) and control (pink) sites species diversity as a function of order q-values for (A) adult tree communities, and (B) regeneration communities. The black dots denote Hill numbers: species richness (q = 0), Shannon diversity (q = 1), Simpson inverse diversity (q = 2). Curve slope reflects evenness of relative species abundances (e.g., the steeper the curve the more uneven the distribution of relative abundances).

Five 1 × 4 m sub-plots were established within each plot to document seedlings and saplings, henceforth regeneration (DBH < 7.5 cm) (Supplementary Figure 1). At the center of each sub-plot, a photo capturing a straight vertical field of view was used to determine canopy cover. Each photo was taken with a Sony Cyber-shot DSC-RX100 Digital Camera with a 35 mm equivalent lens mounted on a tripod. All photos were taken 1 m off the ground and balanced using a hand-held level with the top of the frame facing uphill from the coastal line. Both the exposure and shutter speed settings were automated by the camera. All photos were recorded in the JPEG format with a resolution of 5,472 × 3,648 pixels. Analysis of photos was done in GIMP v 2.10 using the “Select by Color” tool to identify all between canopy gaps. Limited access to these remote sites resulted in some photos being taken without an overcast sky, thus possibly underestimating canopy cover.

Aspect and slope were obtained using LiDAR data captured by Terra Remote Sensing Inc., on Calvert and Hecate Islands in 2012. Climatic variables were not included as the high-resolution climate data required to accommodate the close proximity of our sampling sites is unavailable for our study location. Due to the cultural sensitivity of these sites, we did not have permission to collect soil samples. Consequently, the habitation site’s soil composition reflects a unique combination of purposeful landscape alternations, including controlled burns for forest management, intertidal resource deposits for terracing, and the cultivation of culturally important species (Hoffman et al., 2016a; Trant et al., 2016; Fisher et al., 2019).

Species diversity was characterized using three indices to capture any variation between habitation and non-habitation sites. Parallel analyses were used to evaluate the adult and regeneration communities separately. The first index used was species richness, which denotes the number of unique species present. The second index used was Shannon diversity, the most commonly used metric that explores the relationship between species richness and evenness, placing less weight on dominant species (Jost, 2006; Tuomisto, 2012; Zhang et al., 2012). The third index used was inverse Simpson, which considers both species richness and evenness while weighting for dominant species. Collectively, these metrics are Hills number and can be represented by q-value, which designates the sensitivity of the measurement to the relative abundances where q0 is species richness, q1 is Shannon diversity, and q2 is inverse Simpson (Hill, 1973; Chao et al., 2014b; Cox et al., 2017). In the present study, a value for each of the three indices as a function of order q, in addition to all other q orders until a maximum of five, was generated for each sample plot (n = 25) within habitation and non-habitation sites (Chao et al., 2014a,b; Hsieh et al., 2016). The indices and order q-values were plotted as diversity profile curves, which illustrate Hills numbers as a function of order q (Chao et al., 2014b). The slope of the resulting curve reflects the unevenness of species relative abundances in that the steeper a declining curve, the more uneven the distribution of relative abundances. If the curve is constant at the level of species richness, the relative abundances are completely even (Chao et al., 2014b).

Community composition was characterized using non-parametric multivariate analyses to quantify the differences between ecological communities present at the habitation and non-habitation sites. Similar to the diversity metrics, parallel analyses were used to evaluate the adult tree species and regeneration communities separately. Dissimilarities between both adult and regeneration communities were determined using a nested Permutational Multivariate Analysis of Variance (PERMANOVA) with the “BiodiversityR” package in R (Messier, 1993; Brown et al., 1996; Prescott and Blevins, 2005; Kindt and Kindt, 2019). A nested PERMANOVA was utilized to account for the hierarchical nature of our data. Bray-Curtis coefficients were used to calculate community density matrices of all sites (s = 6). A total of 999 permutations were performed with two fixed factors, treatment type (habitation/control) and site, included in the model. Separate models were generated for the adult and regeneration data. Community compositions were then visualized using non-metric multidimensional scaling (nMDS). Finally, SIMPER analysis was performed using the “vegan” package in R (Dixon, 2003) in which the contribution of species contribution to the community dissimilarity was determined. The average dissimilarity (Diss) between each species to the standard deviation (SD) of the dissimilarities were measured as (Diss/SD) to determine which species are consistently contributing to community dissimilarity (Clarke, 1993). The average dissimilarity values of each species were visualized to illustrate the extent to which communities differed among treatments.

For forest height data, normality of raw data was determined using a Shapiro–Wilk test. Forest stand structure data contained the singular response variable of height (m) for each adult tree sampled (n = 835). To test for maximum height differences amongst the largest trees growing on habitation and control sites, data were separated into four quantiles. By separating all species tallied into their respective quantiles we are able to analysis our focal targets (the largest adult trees present) more efficiently. However, trees growing on control sites occurred at low densities within the fourth quantile prohibiting further analyses. Therefore, the analysis was conducted on the combined third and fourth quantiles; see Supplementary Tables 1, 2 for test results regarding the three quantiles analysis which only considered quantiles one and two. Tree height data for five habitation sites and four control sites were modeled against treatment type (habitation/control), slope and aspect using Generalized Linear Mixed Models (GLMM), where plot was fitted as a random intercept using the “lme4” package in R (Bates et al., 2007). Western hemlock and western redcedar analyses were fit with a Gamma distribution, as model residuals demonstrated a positive right skew, while Sitka spruce data were fit to a Gaussian distribution. The aspect was transformed using a cosine transformation where, x = –1 $\cos [θ (\frac{π}{180})]$ , where θ is the aspect angle in degrees to account for the issues associated with comparing aspect in degrees (Whittington et al., 2005). After transformation, values closer to 1 represent aspects better aligned to receive maximum solar radiation while those closer to −1 represent areas of lower solar radiation (Whittington et al., 2005). A single model was developed for each species tested for height. Thus, no model selection was performed. All models had their residuals plotted against fitted and visualized using the Q-Q plot function.

Finally, if forest canopy regrowth is homogeneous or if nurse logs facilitate regeneration individuals were tallied and counted as either establishing on a nurse log or the forest floor. The amount of space occupied by the nurse logs within the 4 × 1 m plots was also determined as the percent of plot occupied by nurse logs. However, due to logistical constraints, these data were only collected at six of nine sites (three habitation sites and three control sites). Generalized Linear Models (GLM) fitted to a Poisson distribution were performed, where regeneration counts were modeled against location (nurse log or forest floor) and treatment type (habitation/control). All analyses were done using R statistical software v.3.6.5 (R Core Team, 2013).

Results

Community Analyses

Six species of adult trees (western redcedar, western hemlock, Sitka spruce, yellow cedar, Pacific yew, and shore pine), and five species of regeneration (western redcedar, western hemlock, Sitka spruce, Pacific yew, and shore pine) were encountered across all sites. Diversity profile plots illustrate that the species richness of adult tree species was lower on habitation sites as only four of the six species were observed (western redcedar, western hemlock, Sitka spruce, and Pacific yew). Whereas species richness of regeneration was similar among habitation and control sites as four of the possible five species were observed at both site types (Figure 2A). Yellow cedar regeneration was not observed on habitation sites, whereas Pacific yew regeneration was not observed on control sites. The Shannon diversity and inverse Simpson of adult communities were elevated on control sites. In contrast, the opposite pattern was observed within the regeneration communities with diversity metrics elevated on habitation sites (Figure 2B).

Plot level community composition varied within adult communities more so than regeneration communities; however, in each assemblage, composition and species-specific differences in abundance were observed. Adult communities were significantly different between habitation and control sites (nested PERMANOVA: F = 9.13, p = 0.008, 999 permutations). The nMDS plot for adult communities highlights the differences in community composition observed at each plot and had a stress of 0.07, indicating a good model fit between habitation and control sites (Figure 3A; Clarke, 1993). The SIMPER analysis revealed adult communities on habitation sites were 34.64% dissimilar to those on control sites, with all species exhibiting higher densities on control sites (Figure 4A and Supplementary Table 4). The largest contributors of this dissimilarity were western hemlock (Diss/SD = 1.55) and western redcedar (Diss/SD = 1.35). The average dissimilarities of adult western hemlock and western redcedar on control sites relative to habitation sites was roughly 13.65 and 12.65, respectively (Figure 4A). On average, habitation sites exhibited a roughly 2-fold increase in western hemlocks and a 1.5-fold in western redcedars, totaling over 75% of cumulative dissimilarity (Supplementary Table 4).

FIGURE 3

Figure 3. nMDS plots comparing control (pink polygons) and habitation (green polygons) sites community composition for (A) adult tree communities, and (B) regeneration communities. Each point represents a plot, with plots of similar shapes being from the same site. The stress of each plot is denoted in its right corner.

FIGURE 4

Figure 4. Species-specific average dissimilarity between control and habitation sites for (A) adult tree communities, and (B) regeneration communities. Pink dots denote increases in abundances at control sites relative to habitation sites, whereas green dots denote increases in abundances at habitation sites relative to control sites. Average dissimilarity was derived from similarity percentage analyses evaluating adult and regeneration communities’ compositions within habitation and control sites. Error bars represent the standard deviation of the species contribution to the between-group dissimilarity.

Limited differences were observed between habitation and control sites for regeneration communities. A nested PERMANOVA analysis did not detect significant difference between the regeneration communities at control and habitation sites (F = 1.55, p = 0.24, 999 permutations). The nMDS plot illustrates the high amount of overlap between the communities and had a stress of 0.03, indicating a good model fit (Figure 3B; Clarke, 1993). The SIMPER analysis determined that regeneration communities at habitation are 60% similar control site’s communities (Figure 4B and Supplementary Table 5). The largest contributors to dissimilarity in regeneration composition between site types were western hemlock (Diss/SD = 1.34) and western redcedar (Diss/SD = 1.50). The density of western redcedar regeneration was 21.92 stems per plot on habitation sites relative to 13.83 on control sites (Figure 4B). This elevated stem densities corresponded to 1.58-fold increase in western redcedar regeneration on habitation sites relative to control sites. Conversely, the density of western hemlock regeneration was 29.92 on control sites relative to 23.69 on habitation sites (Figure 4B). These densities correspond to western hemlock regeneration stems being 21% more abundant on control sites, which was associated with an average dissimilarity of 18.75. Cumulatively, western hemlock and western redcedar contributed over 85% of the total dissimilarity observed between regeneration communities. The average density of Sitka spruce regeneration was 358% higher on habitation sites. While a seemingly large increase, the actual density of Sitka spruce on habitation and control sites were relatively low (1.33 and 4.77, respectively).

Forest Structure

The average height of western redcedar, Sitka spruce and western hemlock, when considered within three height quantiles, varied among species and according to the quantiles considered. For all species, maximum height was higher on habitation sites; however, this trend was the most pronounced in western redcedar and western hemlock third quantiles (Figures 5A,C). Specifically, the mean height of the third quantile western redcedar was 22.38 m within habitation sites and 17.50 m within control sites, resulting in a statistically significant 4.88 m increase in canopy height (Figures 5A,C and Supplementary Table 5). Similarly, third quantile Sitka spruce on habitation sites was 4.8 m taller in canopy height relative to control sites, as mean height was 24.77 m and 19.97 m, respectively (Figures 5A,C). However, the increase in Sitka spruce height was not significant (Supplementary Table 5). The third quantile for western hemlock was only marginally taller on habitation sites as a 0.5 m increase in mean height was observed relative to control sites (Figure 5B).

FIGURE 5

Figure 5. Mean tree height within the control (pink) and habitation (green) sites for (A) western redcedar, (B) western hemlock, and (C) Sitka spruce. All heights were divided into four even density quantiles; however, due to a minimal number of control trees occurring the 4th quantiles, the 3rd and 4th quantiles were combined. * denotes significant treatment effect, determined by Generalized Mixed Model (GLMM) analysis considering treatment, aspect, slope (Supplementary Table 5).

Forest Regeneration

Nurse log coverage differed only slightly between treatment type, with control sites having 19.55 ± 2.14% SE coverage while habitation sites had 23.43 ± 3.33% coverage within sub-plot plots. The total amount of nurse logs encountered and their average DBH on habitation sites and control sites were 77 (55.22 cm) and 91 (35.69 cm), respectively. There was significantly more recruitment on nurse logs than on forest floors (GLM, Est = 2.13, Z = 22.86, p < 0.01). There was no significant difference between recruitment of regeneration on nurse logs within habitation sites compared to control sites (Figure 6) (Est = 0.09, Z = 1.50, p = 0.13). There were no statistical differences between canopy cover values between treatment type.

FIGURE 6

Figure 6. Boxplot showing the total number of regenerations found on nurse logs versus the forest floor. Pink boxes represent control sites while green boxes represent habitation sites. Outliers are represented as black dots and represent samples that are >1.5 times the interquartile range. There was significantly more recruitment on nurse logs than on the forest floor (GLM, Est = 2.13, Z = 22.86, p < 0.01). Nurse logs made up on average 23.43% ± 3.33% of the forest floors on habitation sites vs. 19.55% ± 2.14% on control sites.

Discussion

The rich cultural history of the Great Bear Rainforest, including its occupation since time immemorial, and the dynamic processes that structure nearshore forest communities throughout the Pacific Northwest of North America offer ideal conditions to evaluate local and landscape-level human-ecosystem interactions (Hoffman et al., 2016b; Trant et al., 2016; McLaren et al., 2018). Specifically, the presence of shell midden deposits provides a unique opportunity to understand how long-term Indigenous habitation, resource-use, and stewardship shapes coastal forest communities. Our results provide evidence that the continuous occupation of habitation sites have had a quantifiable influence on present and future forest composition. The landscape modifications introduced during habitation have altered the coastal forest environment providing conditions that may favor a unique species assemblage and may also allow for the hyperdominance of culturally important forest species, such as western redcedar.

Our diversity analysis demonstrated that diversity was higher on control sites for adults but higher on habitation sites for the regenerative community. Adult diversity, quantified by Hill numbers of the order q, being higher on control sites is in contrast to previous work within the Pacific Northwest where species richness was often higher on sites with shell middens compared to the surrounding matrix (Fisher et al., 2019). While the increase in regenerative diversity is in agreeance with the work of Fisher et al., 2019 this discrepancy in diversity between the woody and herbaceous layer on habitation sites is of particular interest. While the relationship between each vegetative layer and shell middens has been explored independently, perhaps a holistic approach could allude to the drivers of this relationship.

Other work outside of the Pacific Northwest also found that species richness was higher in the herbaceous layer in the presence of shell middens (Vanderplank et al., 2014). In studies of both woody and herbaceous composition, a decrease in woody species abundance was detected as a result in a shift toward a more diverse herbaceous composition (Cook-Patton et al., 2014). Increased forest productivity on shell middens resulting in reduced, or a unimodal pattern of diversity is documented among both the herbaceous layer (Gough et al., 2000; Rajaniemi, 2002), and the woody layer (Trant et al., 2016). However, the maximum species richness was only six and four for control sites and habitation sites, respectively, with the two additional species found being poorly represented within the population. While both yellow cedar and shore pine were found only on control sites, it is unlikely they played a large role in influencing overall community diversity given their small representation in the overall stem count. Instead, the increased abundance and evenness of western hemlock and western redcedar on control sites, represented through a steeper curve (see Figure 2A), likely explain the increased diversity.

Our comparisons of community dissimilarity showed that habitation sites were significantly dissimilar to control sites despite sharing nearly all species present between sites. A SIMPER analysis determined that over 75% of the cumulative community dissimilarity was attributed to changes in western hemlock and western redcedar densities. This was attributed to an increased density of western hemlock and western redcedar’s on control sites compared to habitation sites. However, overall community dissimilarity was low (∼36%) as western hemlock and western redcedar, the dominant species, were seen on all habitation and control sites. This in combination with the low species richness of the zonal coastal forest could likely produce the lower overall dissimilarity value observed in this study. Differences in community composition presented here between habitation and control sites could also be as a result of differences in stand age or successional stage; however, tree core samples taken from habitation sites in our study area have documented to be more than 150 years old (Hoffman et al., 2021), representing an old-growth forest. Therefore, the successional stage is similar to the surrounding forest matrix and that our results are due to nutrient inputs facilitated by the presence of shell middens. However, the models used to derive these inferences were of the simplest form, and lacked features, such as climatic variables, that may influence composition beyond the effects of shell middens.

Similar to adults, both western hemlock and western redcedar regeneration made up over 85% of the community dissimilarity from changes in density. Western hemlock regeneration was higher on control sites while western redcedar regeneration was higher on habitation sites. This inverse relationship between a reduced western redcedar adult density on habitation site and increased regeneration density is of particular note given their limited dispersal range compared to western hemlocks and Sitka spruce (Minore, 1990). While forest regeneration in this system is primarily driven by gap-forest dynamics, we found little difference in canopy cover between control and habitation sites. However, given the differences in stem density between habitation and control sites, complex canopy differences, such as angular light penetration, that were undetectable within our study design could be present. The ability of western redcedar regeneration to establish in higher density on habitation sites may be a result of these differences and could be the focus of future studies within this system. Much like the assessment of our adult community, overall species richness was even lower in regeneration communities and compositional differences were expressed most in changes of species dominance and evenness.

The results of our structural data demonstrated species-specific responses to tree communities on habitation sites. Interestingly, this signal of increased maximum height was only detected in western redcedars and not in western hemlock or Sitka spruces. The increased sizes of western redcedar on habitation sites compared to the surrounding matrix are in agreeance with work in this region (Trant et al., 2016). Variability in species-level responses may be due to interspecific competition for Nitrogen (N) and Phosphorus (P) despite soils on habitation sites being less nutrient limited. Soils on habitation sites are less acidic than those in surrounding areas, which have been shown to decrease soil C:N ratios and elevate rates of N mineralization, leading to an increase in the availability of useable N for above ground vegetation (Fisher et al., 2019). Western redcedar have been shown to grow taller and be more productive when exposed to higher levels of P and N (Messier, 1993; Brown et al., 1996; Prescott and Blevins, 2005), yet this does not fully explain our findings. Competitive exclusion through elevated resource acquisition of soil nutrients may be a likely driver of the patterns seen within this study.

As expected, nurse logs in coastal temperate rainforests are the primary substrate for recruitment and forest regeneration and are consistent with past work in the Pacific Northwest where nurse logs contained over 94% of Sitka spruce and western hemlock regeneration (Graham and Cromack, 1982). Despite nurse logs containing fewer nutrients than surrounding soils (Takahashi et al., 2000), they allow for a competitive advantage for light by avoiding shading by the herbaceous layer. However, the advantages provided by a greater access to light come at a cost, as previous work has documented an increase in mortality rates associated with establishment on nurse logs (Harmon and Hua, 1991). were no observable differences in recruitment between both regenerative pathways on habitation and control sites despite the compositional differences observed in adult species. Given the contrast in stem density between control and habitation sites, the lack of signal detected in regeneration density between sites is of note. While nurse logs are the primary substrate for recruitment, they occupy less than 30% of the forest floor within our study plots, similar to other studies within the same region (Graham and Cromack, 1982; Christy and Mack, 1984). The increased stem density on control sites appeared to not be reflective of an increased density of nurse logs found within our plots as seen by their percent cover. While more nurse logs were found on control sites in observed plots, they were on average of a smaller DBH which may explain the lack of a significant difference between site types. This is however, speculative and the results presented here only represent preliminary findings, as other components such as fungi composition, nurse log decay stage, and height off ground were not explored.

Our results contribute to a growing body of literature in which Indigenous communities have left long standing impacts on systems previously thought of as free of human influence. The scale and degree of these effects are likely varied amongst the many habitation sites that exist along BC’s Central Coast. While this study only focused on habitation sites within a localized area along the Central Coast, these long-term impacts are likely expressed in other regions that have experienced similar Indigenous stewardship. Within Calver and Hecate Island, forest stands on habitation sites saw changes in diversity, expressed almost exclusively by compositional changes through a reduced density of adult western hemlock and western redcedar and the increased maximum height of western redcedars. While these changes at the site level may appear minor, the Central Coast of BC is home to thousands of habitation sites, which likely demonstrate the same shifts in composition seen here. Given the persistence of shell midden deposits in these environments, the results presented here represent only a singular “snapshot” in time, by which the results of future studies may confer too. How these shell middens will continually modify above ground forest stands remains to be seen, but by accounting for the role that Indigenous communities played in managing these landscapes, we have a greater insight to the underlying ecological processes that have altered these forest communities. The quantification of these persistent effects allows researchers to generate a better understanding of the factors that may be governing their desired ecosystem. Here, we encourage researchers to carefully consider the effects that Indigenous groups may have had on the lands they wish to study to account for any persistent effects present within their system (Schang et al., 2020).

Data Availability Statement

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

Author Contributions

KS and AT contributed to conceptual development, experimental design, and data collection. KS led the writing of the manuscript. All authors contributed to data analysis and critically to the editorial process, and gave final approval for publication.

Funding

This research was supported by funds from the Hakai Institute, the Natural Sciences and Engineering Research Council (NSERC) Discovery Grant (AT).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

We would like to thank Eric Peterson, Christina Munck, Will McInnes, and the staff of the Hakai Institute who helped facilitate this work, and Heiltsuk Integrated Resource Management Department (HIRMD) and Wuikinuxv Resource Stewardship Committee. We would also like to thank Brad Fedy, Stephen Murphy, Alana Closs, Sara Wickham, Emma Davis, and Alexandra Johnson for their insightful feedback throughout the preparation of this manuscript, and Sharon Wilson and Jennifer Walkus for their assistance with research development and data collection.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fevo.2021.791047/full#supplementary-material

References

Banner, A., Moran, J., and Groot, A. D. (2005). Pattern, Process, and Productivity in Hypermaritime Forests of Coastal British Columbia. Victoria, BC: B.C. Ministry of Forests and Range Research Branch.

Google Scholar

Bartz, K. K., and Naiman, R. J. (2005). Effects of salmon-borne nutrients on riparian soils and vegetation in Southwest Alaska. Ecosystems 8, 529–545. doi: 10.1007/s10021-005-0064-z