The influence of climate change, site type, and disturbance on stand dynamics in northwest British Columbia, Canada

Stand and disturbance dynamics are key processes that need to be assessed along with climate-species interactions if we are to better understand the impacts of climate change on species. In this study we investigated the biotic interactions (competition) between species, the influence of disturbance type, and changes in resource availability (moisture and light) on the response of six tree species to climate change in the northwest region of central British Columbia, Canada. Two ecological models were parameterized, linked together and coupled to climate change scenarios to explore the interactions between: (1) the response of species in the regeneration phase and (2) the role of disturbance, resource availability and competition on determining stand composition and productivity. Climate change was found to reduce soil moisture availability which resulted in a decline in regeneration potential for all species on dry sites and negative to neutral responses on sites with higher water availability. Following fire, stand dynamics and composition were modeled to undergo significant changes under the 2080s climate compared to current climate conditions on dry and mesic sites. Changes in stand dynamics under climate change were marginal following bark beetle disturbances. While significant changes to stand dynamics were found on dry sites, the presented results suggest that the sites with the highest moisture availability maintain the same general stand dynamics and composition following disturbances under climate change. This study highlights the need to consider species response to climate change in interaction with existing stand conditions, disturbance type, competition, resource availability, not just temperature and


INTRODUCTION
Predicted changes in climate are expected to have significant impact on forest ecosystems through increases in disturbance events such as drought, fire, insect and disease outbreaks.This will in turn facilitate the migration of species to new ranges, where possible, and extirpation or extinction where not, which may lead to changes in ecosystem composition (IPCC 2007).Shifts in species distributions are expected to occur along environmental gradients if the environmental tolerances of species are exceeded (Miller and Urban 1999).These shifts will be complex and it may be difficult to determine whether shifts will be due to direct responses to climate change or indirectly through interactions with other species (Brzeziecki et al. 1995).
Under stable or slowly changing climate conditions, it is typically assumed that species distributions are in equilibrium with climate.However, this assumption cannot be applied globally and must be considered regionally and on a species by species basis (Pearson andDawson 2003, Arau ´jo andPearson 2005).At the macro-scale, climate is the principal driver of species distributions (Huntley et al. 1995); however, autogenic, allogenic, and biogenic processes interact with species ecophysiology and resource availability to affect species distributions at finer scales (Pearson et al. 2004, Guisan andThuiller 2005).Autogenic processes are driven by competition which can have a significant influence on species distributions (Arau ´jo and Luoto 2007).Arau ´jo and Luoto (2007) argue for the need to consider these interactions, particularly competition, at the site -level to improve our ability to assess the impact of climate change on species distributions at the landscape-level.
Allogenic factors such as fire, drought, and frost as well as biogenic agents such as insect and disease epidemics all interact with resource availability to affect species distributions (presence and abundance) and growth (Pulliam 2000).These interactions can be both spatially and temporally complex with biotic interactions occurring in one time period or at one scale that may not hold true in another due to changes in climate (Arau ´jo and Luoto 2007).The type, severity and timing of disturbance events also influences species distributions which can lead to declines in abundance and range of one species and the expansion of others at the site and landscape-scale (see Allen and Breshears 1998).The magnitude and/or frequency of a disturbance interact with each species' unique physiology, demographics and life-cycle characteristics and cause divergent responses to disturbance characteristics (Walker 1989).Species also vary in their tolerance to environmental stresses in different parts of their life cycle with their highest sensitivity occurring during their regeneration phase (Walck et al. 2011).This creates a potential interaction between regeneration, climate and disturbance type that is important for understanding the effect of changing environmental conditions on forest ecosystems (Guisan andZimmermann 2000, Astrup et al. 2008).
Modeling is typically used to predict individual species, ecosystem or biome responses to changes in climate (Zolbrod and Peterson 1999).Biome modeling generally focuses on broad lifeforms or plant functional types versus individual species (Gavin and Brubaker 1999).This generally results in species-specific responses being masked by the broad scale of the analysis (Bakkenes et al. 2002).Individual species modeling is typically done using one of two techniques: statistical modeling (e.g., Thuiller et al. 2005) or mechanistic modeling (e.g., Miller and Urban 1999).Statistical models generally take the form of bioclimatic envelope models and typically assume that species are at equilibrium with climate and do not consider biotic and abiotic factors directly (Thuiller et al. 2005, Arau ´jo andLuoto 2007).These types of models generally reflect the realized niche (i.e., ecological) of species but provide no explanation of the physiological mechanisms that may drive change and are poor predictors of species dynamics under changing environmental conditions (Guisan and Zimmermann 2000).These models are also sensitive to the spatial extent and resolution of data used for their development, which limits the ability to extrapolate over space and time (McKenzie and Halpern 1999).Despite these limitations statistical models can be useful for modeling species responses at multiple scales if both direct climatic and resource factors are used but where allogenic (disturbance) and autogenic (competition) factors are not considered (Guisan and Zimmermann 2000).The incorporation of disturbance and dynamics is important for understanding the impacts of climate change on species response and requires the use of mechanistic modeling to provide more accurate predictions of individual species response to climate change (see Guisan and Zimmermann 2000).
In this study we use a modeling approach to investigate how disturbance type, site conditions (moisture and nutrients), competition and changing climate interact to produce variation in tree species productivity and community responses to climate change.We carry out our study in the sub-boreal forests of British Columbia (BC), Canada where local soil resource gradients lead to large local variation forest composition and dynamics.For this region, as in many parts of the world, species response to climate change has largely been evaluated using bioclimatic envelope modeling (Hamann and Wang 2006) which fails to illustrate the combined effect of disturbance type, site conditions (moisture and nutrients), competition and changing climate.Here we linked two ecological models to explore the interactions between climate change, disturbance, resource availability and competition on stand composition and productivity and evaluate the ecosystem response in terms of (1) regeneration potential of species, (2) the sensitivity of species to climate change, (3) the effect of disturbance type on the sensitivity of species

STUDY AREA
The study area is located in the Bulkley Valley of northwest British Columbia (see Fig. 1).The study area is approximately 600,000 ha in area and has an elevation range from 356 m to 2736 m.In the study area, the Sub-boreal Spruce (SBS) ecosystem is found from 500 m to 1100 m.The SBS ecosystem has a continental climate with warm and moist summers and cold winters with extended periods of snow (Banner et al. 1993).The SBS zone is considered a transitional zone between the montane forests of in southern BC, the boreal forests of northern BC and the subalpine forests that occur at higher elevations within the central interior (Pojar et al. 1982).The forests in the SBS ecosystem are dominated by interior spruce (Picea engelmanni 3 glauca), subalpine fir (Abies lasiocarpa) and lodgepole pine (Pinus contorta var.latifolia).Pine is often the dominant species to regenerate after wildfires in the region.Trembling aspen (Populus tremuloides) can also be common in early seral communities.Local occurrences of black cottonwood (Populus balsamifera ssp.trichocarpa) and western hemlock (Tsuga heterophylla) occur within the study region.Fig. 1 shows the location of two weather stations within the study area which represent the SBS ecosystem in the Bulkley Valley.The average climate conditions of the study area are presented in Table 1 and Fig. 2.

Linking ecological models
The ecological model TACA (Tree and Climate Assessment) is a mechanistic species distribution model (MSDM) that analyses the response of trees to climate-driven phenological and biophysical variables (Nitschke and Innes 2008).It assesses the probability of a species to be able to regenerate, grow and survive under a range of climatic and edaphic conditions.The modeling approach reflects the regeneration niche of a species, since presence is directly related to establishment (McKenzie et al. 2003).
The original TACA model developed by Nitschke and Innes (2008) was modified for this study.First, the model was changed to incorporate a frost free period mechanism as Hamann and Wang (2006) found that the annual number of frost days had a significant interaction with observed species ranges in BC (see Fig. 3a).Next, the phenology component of TACA was improved to increase the interaction between chilling, heat sum accumulation, frost, and budburst based on Bailey and Harrington (2006) (see Fig. 3b).The improved phenology component integrates the species chilling requirement with the accumulation of heat sum which then interacts with frost events that delay bud burst and/or causes frost damage after bud burst occurs.The soil moisture function was modified to incorporate the Hargreaves model of evapotranspiration (Hargreaves and Samani 1985) and estimates of direct and diffuse daily solar radiation based on equations from Bristow and Campbell (1984) and Duarte et al. (2006).In addition, snowfall, snow accumulation, and snowmelt (timing and quantity) were considered in this version of the model.Snowfall and accumulation are tracked in snow water equivalent (mm) while snowmelt utilizes the model developed by Brubaker et al. (1996) which models snowmelt as a function of temperature and daily solar radiation.Finally, the soil  1943-2000 1966-2000 v www.esajournals.orgcomponent of TACA was expanded to allow for three different soil types (texture, coarse fragment percentage, and rooting zone depth) to be run simultaneously allowing for the representation of multiple edaphic conditions across the soil resource gradient present in the study area.
To study stand-level forest dynamics we use the individual tree, spatially explicit mixed species forest model SORTIE-ND to study stand-level forest dynamics.SORTIE-ND is a descendent (restructured and reprogrammed in Cþþ) of the small-scale disturbance model SORTIE developed in the early 1990s in the northeastern US (Pacala et al. 1996) and SORTIE/ BC, a version of SORTIE further refined for the forests of northern BC (Coates et al. 2003).The ND signifies the model's focus on local neighborhood dynamics.SORTIE-ND has an intermediate position between pure empirical and process-based models (Messier et al. 2003).The model predicts individual tree growth, survival, dispersal and recruitment as a function of neighborhood conditions throughout the grow-ing season.SORTIE-ND also allows for the deterministic simulation of disturbance events.
In SORTIE-ND the forest is represented by a large collection of interacting trees that are followed both in time and space.Those trees are currently divided among seedlings, saplings, adult trees and snags.Population-level dynamics are simulated by summing the collective activities of numerous individuals.Each tree is a discrete object that is described with various attributes (size, growth rate, crown dimensions, etc.).Each tree's (individual) behavior is modeled with functions that describe the interactions with other individuals (e.g., effect of species and distance of neighbors on growth of individual trees) or its environment (e.g., growth of seedlings in relation to available light levels).Many of the interactions have non-linear relationships and/or have random events associated with them.The non-linearity of many interactions, the stochastic behavior of some objects and processes, and the large number of objects, functions and stochastic events makes SORTIE-  Growing season terminates at the cessation of growing degree accumulation (species-specific) and/or by user defined date to account for declining photoperiod length (October 15 in this study).Frost events following termination of the growing season do not cause mortality/ damage except for lethal frost events as recorded in (a).Texture class defines the field capacity and available water storage capacity (AWSC) of a soil in mm/m of soil.Percolation parameter determines if soil losses groundwater when AWSC is exceeded or receives groundwater from upslope sites; for mesic sites, percolation is set to zero.
ND a good example of a modeling approach aimed at being able to represent interactions and complex behavior in forests.
To incorporate climate change into SORTIE-ND the results from the TACA model were incorporated in two ways.First, the change in regeneration coefficients between current and predicted climate change was used to modify the establishment parameters for each species.The species response in the current scenarios was used as a baseline with the establishment parameters modified by the proportion of change that species exhibit under the climate change scenarios.Secondly, the TACA model tracks changes in soil moisture availability on a site as a function of actual to potential evapotranspira-tion (AET/PET).The AET/PET ratio is used in British Columbia to classify sites into an actual soil moisture regime (ASMR).The change in the AET/PET ratio between climate scenarios was used to determine if the ASMR at each modeled site type would shift along an edaphic resource gradient under climate change.Shifts in the ASMR were used to adjust the growth rates of each species based on the current relationships between ASMR and productivity (as measured by Site Index) that exist in the province of BC.The approach used to link TACA to SORTIE-ND is detailed in Fig. 4 and the productivity modifiers summarized in Table 2. Table 2. Percent change in species productivity for the three site types due to changes in soil moisture regime modeled to occur under 2080s climate change.

Site type
Interior spruce

Subalpine fir
Trembling aspen

Western hemlock
Black cottonwood v www.esajournals.org

Model parameterization
Site types and stand conditions.-Seventyone fixed radius sample plots were established in the study area around Smithers, BC.At each sample site, soil depth to parent material, rooting zone depth and soil texture was measured, and the soil moisture and nutrient class of each site was classified based on the presence and absence of tree, shrub, herb and moss species that are used as indicators for delineating distinct site series (combinations of soil moisture and nutrient classes) (after Banner et al. 1993).This information was used to categorize sample sites into three site types based on actual soil moisture and nutrient regimes that are characteristic of the study region (Kranabetter and Simard 2008).The first site type represents crest and upper slope positions that had loamy sand texture with an average soil depth of 0.75 m, rooting depth of 0.27 m and coarse fragment content of 52%.This site type represents xeric to sub-xeric (moderately dry and poor) edaphic conditions found in sub-boreal forests of the study area.The second site type represents mid-slope positions with a sandy loam texture, average soil depth of 2 m; rooting depth of 0.34 m and coarse fragment content of 36%.This soil type represents submesic to mesic (slightly dry and medium) edaphic conditions found across the study area.The final site type represents lower slope and level positions with a loam to silt loam texture with an average soil depth of 2 m, rooting depth of 0.33 m, and coarse fragment content of 27%.This site type represents mesic to sub-hygric (fresh-moist and rich) edaphic conditions found across the study area.
Initial stand compositions.-Sixspecies were identified in the field sampling plots and were included in the initial stand structures: lodgepole pine, interior spruce, subalpine fir, trembling aspen, black cottonwood and western hemlock.The stand structure (density of seedlings, saplings, and mature trees) and composition information from the sampling plots were used to randomly generate three stem maps for each site type which represented the starting stand compositions and structures.Each replicated stand had a unique initial stem spatial distribution as v www.esajournals.orgthe X and Y coordinates of each stem were the elements randomized.Fig. 5 summarizes the stand composition information used to generate the stem maps.The first site type was dominated by lodgepole pine; the second site type was characterized by a mix of lodgepole pine, interior spruce, subalpine fir and trembling aspen and the third site type was dominated by subalpine fir and interior spruce.
Modeling scenarios.-Themodeling analysis involved the assessment of stand dynamics following disturbance with the empirically based SORIE-ND parameters for three site types: (1) dry and poor, (2) mesic and medium, (3) moist and rich.For each site type (n ¼ 3), three different initial stand conditions were replicated (n ¼ 3) under two climate scenarios (observed and 2080s climate change; n ¼ 2).For each site type, four disturbance simulations (n ¼ 4) were applied: no disturbance, fire; mountain pine beetle (Dendroctonus ponderosae) (MPB) attack, and spruce bark beetle attack (D. rufipennis) (SBB).The 72 model simulations (3 site types 3 3 site type replicates 3 4 disturbances treatments 3 2 climate treatments) were run over a time horizon of 200 years.The four deterministic disturbance treatments were characterized by the following parameters: (1) for no disturbance, both stand and substrate parameters were not modified; (2) for fire, 95% of the stems of all species were killed at time step five of the simulation; in addition, substrate parameters were modified to reflect the effects of fire on forest floor with the proportion of mineral soil increased to 35%; (3) for MPB disturbance, a severe outbreak was simulated with 100% of the lodgepole pine trees larger than 10 cm in diameter killed at time step five of the simulation; and (4) for SBB disturbance, a severe outbreak was simulated with 90% of the spruce trees with a diameter .10cm and ,30 cm, and 100% of all trees with a diameter .30cm were killed.Substrate parameters were not modified for the MPB and SBB disturbances.Fire removes 15-100% of the forest floor (Seedre et al. 2011); in our study we assumed that fire consumed the forest floor but left the availability of other substrates (mounds and logs) unchanged resulting in the 35% mineral soil proportion.The 95% mortality rate was used to reflect that individual stems do survive fires within stands.The 10 cm diameter threshold for lodgepole pine mortality is based on Safranyik (2004) who identified that mortality of host trees in stands affected by severe MPB outbreaks occurs down to the 8 to 10 cm diameter class.The 10 cm diameter threshold for spruce mortality is based on Safranyik (1985) and Berg et al. (2006) who identified that spruce trees with a diameter less than 10 cm diameter typically survive severe SBB outbreaks.Berg et al. (2006) and Werner et al. (2006) also reported 70-100% mortality in larger diameter (.24 cm) classes in severely affected stands in Alaska which was used to prescribe the 90% mortality rates in the 10-30 cm diameter classes while trees greater than 30 cm in diameter typically suffer higher morality in severe SBB attacks (Hard et al. 1983, Safranyik 1985) so 100% mortality was prescribed.The disturbance treatments are simplifications of highly variable processes but were developed to represent severe disturbance events that clearly changes stand structure in order to illustrate the potential effects that disturbance in interaction with climate and species ecology may have on future stand composition.
For trembling aspen, the ability to regenerate vegetatively was enabled while for lodgepole pine the vegetative regeneration parameter was modified to represent its aerial seed bank (Lotan and Critchfield 1990) that is released following disturbance.For each lodgepole pine individual it was assumed that 144 seeds could be released following a disturbance.This was calculated based on the average number of serotinous cones and seeds/cone found on mature lodgepole pine (Lotan and Critchfield 1990).Interior spruce and subalpine fir lack the ability to re-sprout and do not form aerial seed banks but they do form large seedling banks in the understory (Alexander et al. 1990, Nienstaedt andZasada 1990); particularly on the moister and richer site types of the study area (Astrup et al. 2008).
We assumed that the species in the study area have the adaptive capacity and plasticity to shift productivity in response to changes in temperature and soil moisture availability.Strong clinal variation does exist between populations of these species which will likely result in adaptation lags that may prevent species from capitalizing on changes in climate (Aitken et al. 2008).For the purpose of this study we assume that species productivity will shift due to changes in soil moisture and temperature and that using species productivity from a warmer region with soil moisture regimes that are similar to those predicted in the future is a sufficient proxy to parameterize the growth rates of species in the future SBS.
In modeling climate, we used multiple scenarios based on historical climate data from local weather stations and global climate change model (GCM) predictions for the region.A direct adjustment approach (Hamann andWang 2006, Nitschke andInnes 2008) was used to integrate climate change scenarios into the historical climate records for the 2 climate stations that represent the SBS zone in the study area (Table 1).For each climate station 10 weather scenarios representing the 10th, 25th, 50th, 75th and 90th percentiles for mean annual temperature and annual precipitation were used to represent the observed climate scenarios assessed in TACA.Two GCMs, the Canadian GCM2 (Flato et al. 2000) and Hadley CM3 models (Johns et al. 2003) found by Bonsal et al. (2003) to be the best and second best, respectively GCM for predicting historic temperature and precipitation in BC.The regional climate change predictions for the SBS zone were obtained from the Pacific Climate Impacts Consortium (2009).The Intergovernmental Panels 4AR emission scenarios: A1B, A2x, and B1x were used to represent an ensemble of future 2080s climate conditions.The mean values of these projections are a 3.18C (range: 2.28C to 3.88C) increase in mean annual temperature and a 16% (range: 12-24%) increase in mean annual precipitation.

Analysis of simulation results
In order to better understand the outcomes of the many simulations we carried out a simple statistical analysis of the modeling results to determine if differences exist in regeneration, stand composition and stand productivity (as a function of basal area).Regeneration coefficients generated by the TACA model were logit transformed (Warton and Hui 2011) to meet the assumptions of normality required for ANOVA (Zar 2010).A two-way ANOVA on a completely randomized design was used to assess the effect of climate, site and climate 3 site on species establishment and soil moisture (n ¼ 6).To test for changes in productivity a one-way ANOVA was used to assess the impact of climate change on each species' productivity for each disturbance treatment.Mean annual increments (MAI) (m 2 ha À1 yr À1 ) for each species and each stand were used for this assessment.To assess changes in stand composition the proportion of a stand occupied by each species was assessed.Species proportions were Logit transformed to provide normal distributions and the Kruskal-Wallis test used to test if a species had a different distribution (Zar 2010) as a result of climate change for each site and disturbance type (n ¼ 15).Subsamples were selected from year 40, 80, 120, 160 and 200 of each replicated simulation to include temporal changes in stand composition across each simulation.These analyses should not be used to infer statistical significance for use in prediction but to help highlight the difference in species response and stand dynamics over time.

Regeneration responses
The results of the TACA analysis identified that all species exhibited significant model changes in their regeneration coefficients (RC) as a result of climate change ( p , 0.001; p 0.0776 for western hemlock).Significant differences were detected for all species in regeneration coefficients between sites ( p , 0.001).Lodgepole pine and trembling aspen were the only species to exhibit a significant climate 3 site interaction ( p , 0.001).Both lodgepole pine and trembling aspen incurred a significant decline in their RC on site type 1 (dry-poor) due to climate change.Black cottonwood, interior spruce and subalpine fir were all modeled to incur significant declines in their RCs across all sites.For black cottonwood, its site type 1 RC was modeled to decline to zero.Species responses are illustrated in Fig. 6.

Moisture availability
Climate change resulted in a change in actual soil moisture regimes (ASMR).The ASMR for all sites shifted significantly ( p , 0.001 in all cases) towards a drier ASMR.Site type 1 shifted from an ASMR of moderately dry to very dry, site type 2 from fresh-slightly dry to moderately dry, and site type 3 shifted from moist to fresh.The 2080s ASMR were similar to the ASMR gradients v www.esajournals.orgidentified in the southern interior region of British Columbia by Nitschke (2006).Species productivity changes based on the shift in ASMR are provided in Table 2.

Productivity
Site type 1: dry-poor.-Onsite type 1, a decline in productivity was modeled to occur under 2080s climate conditions compared to the observed climate following all disturbances (see Table 3).Lodgepole pine productivity was expected to decline by 24% (Table 2); however, its productivity declined by 47% following fire but increased following other disturbance types as the species benefited from the decline in site suitability of the other species due to climate change (Table 3).Following fire, lodgepole pine had to compete with trembling aspen which was able to increase its productivity on the sites 1.68fold.Interior spruce and subalpine fir exhibited declines in productivity following all disturbances (Table 3).Black cottonwood was unable to regenerate and therefore had zero productivity in the 2080s.
Site type 2: mesic-medium.-Onsite type 2, stand MAI declined by 1% due to the fire and climate change treatment but by 11-18% for the other disturbances (Table 3).Lodgepole pine productivity was modeled to change by 21% to À40% though the assumed change was a 1%  2).Subalpine fir experienced declines in productivity following all disturbances while interior spruce experienced declines following fire and SBB disturbances.Both species were expected to suffer larger declines in productivity (Table 2).Trembling aspen was expected to exhibit a 12% loss in productivity, instead an increase was modeled following fire (16%) but À6% to À15% changes following the other treatments.The 21% increase for lodgepole pine and 16% increase in trembling aspen productivity following fire highlight the potential increase in the competitive ability of both these species following fire on this site type under a 2080s climate.
Site type 3: moist-rich.-Onsite type 3, only small changes in stand productivity were detected (Table 3).Stand productivity remained relatively unchanged though a significant decrease in black cottonwood's productivity was modeled to occur by the 2080s following fire, MPB and SBB attack.On site type 3 a 16% increase in productivity was expected for lodgepole pine (Table 2); interestingly, following fire, a 27% increase in MAI was modeled (Table 3).For trembling aspen a 14% decline was expected; however, a À4% to 10% change in MAI was modeled.For interior spruce and subalpine fir, 3% and 11% increases were predicted but À11% and 2% changes were modeled following fire respectively.Interior spruce had 0 to 10% changes following other disturbances while subalpine fir exhibited a 3% increase following SBB and no change following the other treatments.Positive increases in productivity were modeled for western hemlock across most treatments in the 2080s.

Stand composition
Site type 1: dry-poor.-Atthe species-level, lodgepole pine was modeled to increase in its proportion of the stand in the absence of disturbance (ND) ( p 0.024) and following both MPB ( p 0.026) and spruce bark beetle attack ( p 0.071).Following fire, site type 1 was modeled to suffer a significant change in stand composition with a decline in lodgepole pine ( p 0.026) and an increase in trembling aspen ( p , 0.0001) in the 2080s climate compared to the observed climate Interior spruce loses abundance due to climate change across all treatments with all changes being significant ( p 0.05).Subalpine fir was modeled to incur significant declines following fire and climate change ( p 0.078).Black cottonwood declined to zero ( p , 0.01) in all treatments due to climate change while western hemlock incurred significant changes under climate change ( p , 0.05) with decreases in all treatments except following fire.Species responses are illustrated in Fig. 7.
Site type 2: mesic-medium.-Asignificant change in stand composition was detected following fire on this site with lodgepole pine significantly increasing ( p 0.02) and black v www.esajournals.orgcottonwood declining significantly ( p 0.003).Trembling aspen increased following fire but change was not significant.All other species suffered non-significant declines in abundance.No significant changes in stand composition were detected following the ND, MPB and SBB treatments.Black cottonwood did exhibit a significant decline in abundance following the SBB disturbance ( p 0.028).Species responses are illustrated in Fig. 8.
Site type 3: moist-rich.-Nosignificant changes in stand composition were detected following the ND, MPB and SBB treatments.A significant change in stand composition was detected following fire on this site with lodgepole pine increasing by 25% ( p 0.02) and black cottonwood declining by 86% ( p 0.038) under the 2080s climate scenarios.Black cottonwood also exhibited a significant decline in abundance within the moist-rich sites following a SBB disturbance ( p , 0.001).Western hemlock exhibited its largest response on this site type increasing from 3% to 5% of the stand composition following fire though the change was not significant.Species responses are illustrated in Fig. 9.

DISCUSSION
In British Columbia (BC), as in many parts of the world, species response to climate change has been largely evaluated using bioclimatic envelope modeling.In the Bulkley Valley region of northwest BC it was predicted by Hamann and Wang (2006) that the region's climate would transition to a climate that is currently typified by ponderosa pine (Pinus ponderosa) and Douglas-fir (Pseudotsuga menziesii ) dominated ecosystems by the 2080s.Hamann and Wang (2006) used a bioclimatic envelope modeling approach to determine individual species responses and assumed that species would maintain similar realized niches and associations under a changing climate.These approaches generally provide v www.esajournals.orgno explanation of the mechanisms of change and are considered poor predictors of forest growth and succession under different environmental conditions (Kimmins 2004).In this study, we used a mechanistic species distribution model, TACA, driven by empirically derived climatic and soil parameters to determine species response to climate and then linked the outcomes to SORTIE-ND, a spatially explicit individual tree model.Our results illustrate a much more complex interaction between species, climate, edaphic conditions and disturbance and highlight the need to consider more than just bioclimatic relationships when assessing the impacts of climate change on species distributions and abundance.Guisan and Zimmermann (2000) identified that to improving predictions of species response to changing environmental conditions competitive, biotic and abiotic agents must be considered; the findings of this study supports this assertion.

Impact of changes in climatic and resource availability on species regeneration potential
It was found that climate change had an impact on the regeneration potential of the studied species and these responses were not uniform among species and site types.The species-specific changes in regeneration potential and potential productivity suggest that species will respond individualistically to climate change.Individualistic responses to climate change are expected to result in changes in ecosystem composition and structure (Bartlein et al. 1997).Walck et al. (2011) stated that climate change may affect species regeneration ability which can lead to shifts at the population and community level.In this study changes during the regeneration phase did result in changes in stand composition on site types 1 and 2; particularly after a fire.

Interaction between regeneration, competition, and disturbance: composition
To understand if species will be sensitive to v www.esajournals.orgclimate change both competitive and biotic and abiotic agents need to be considered (Guisan and Zimmermann 2000).Linking TACA model results to SORTIE-ND allowed for the investigation of the multiple roles of disturbance, resource availability and competition on species response to climate change.The findings of this study support this mechanism and highlight the importance that edaphic conditions may have in mediating species response, even under competition.Stand conditions that provide wetter edaphic conditions are regarded as being important for maintaining species that cannot tolerate climatic change (Aide andRivera 1998, Burke 2002).For subalpine fir and interior spruce, the combination of a disturbance that maintains a forest canopy and soil properties that provide higher soil moisture availability are likely to be important for preventing a shift towards pine-dominated ecosystems in the future.
The findings also highlight that climate climax plant communities (site type 2) and edaphic climax communities (site types 1 and 3) are likely to behave differently under climate change.The interaction of physiological and edaphic factors that were modeled to occur on site type 3 sites are likely to play an important role in allowing the current sub-boreal spruce ecosystem species to persist under future climates as edaphic communities within a landscape of newly developing climatic communities (Theurillat and Guisan 2001) likely mediated by an increase in fire activity which is expected under climate change (Stocks et al. 1998).The responses of species on different site types is an important finding of this study and suggest that some species may be able to cope with climate change at the local-scale where suitable resources are maintained even following disturbance (Theurillat and Guisan 2001).
On the drier site type 1, a shift to lodgepole pine dominated forests occurred following fire and SBB attack.The dominance of pine following fire in the 2080s climate did differ significantly from the current climate which was reflected in v www.esajournals.orgthe change in stand composition which shifted towards more trembling aspen, very few spruce and no subalpine fir.Following all disturbances, productivity decreased; however, the largest declines occurred after fire and MPB attack due to the inability of lodgepole pine to regenerate under the low light conditions following the latter and outcompete trembling aspen following fire.The inability of interior spruce to regenerate on organic substrates, which are maintained following SBB and MPB disturbances, along with declines in climatic suitability results in the recruitment of interior spruce being severely restricted which in turn limits productivity under the 2080s climate.Astrup et al. (2008) identified that low light and substrate availability are important factors that are currently limiting recruitment of pine and spruce, respectively, in stands following MPB attack within the study region.Axelson et al. (2009) also found that regeneration of pine and spruce following MPB attack is a slow process that is affected by substrate suitability and the presence of grasses; as well as light availability for lodgepole pine.For subalpine fir, the decline in regeneration potential resulted in the absence of recruitment due to increased soil moisture stress, particularly after fire, and the reliance on advanced regeneration following MPB and SBB attack.Cui and Smith (1991) found that subalpine fir is sensitive to soil moisture stress which can prevent establishment on exposed sites.Ettl and Peterson (1995) also found that subalpine fir is sensitive to variation in soil moisture with species response to climate also being mediated on sites which maintain suitable soil moisture and exacerbated on sites where moisture is limiting.The response of subalpine fir in this modeling study reflects this behavior.
On site type 2, an increase in lodgepole pine and trembling aspen were modeled following fire with decreases in all other species suggesting that fire may be a driver for change on these sites.On site type 3, species exhibited the highest resilience to climate change and disturbance.On these sites, following MPB and SBB disturbance events, stands were modeled to maintain their current composition and structure of species and productivity.Slight changes in composition were modeled following fire but these changes were not significant.Subalpine fir and white spruce (Picea glauca) (interior spruce is a hybrid of white spruce) have both been documented establishing under sheltered conditions while failing to establish under exposed conditions due to moisture stress (Cui andSmith 1991, Man andLieffers 1999).Interestingly, subalpine fir, which is currently benefitting from a MPB epidemic within the study region (Astrup et al. 2008), could exhibit a similar response across all sites under the 2080s climate.Leonelli et al. (2008) predicted that trembling aspen in northeast BC will likely be negatively affected by climate change on sites where moisture resources are limiting but should also be able to maintain growth or increase productivity on sites where moisture is not consistently limiting.The response of trembling aspen in this study supports the hypothesis of Leonelli et al. (2008).
Black cottonwood was modeled to decline across all sites and disturbances due to climate change.Mild soil moisture deficits can prevent regeneration in this species while its productivity is sensitive to soil moisture availability (Krasny et al. 1988, Yarie 2008).Changes in soil moisture availability were modeled across all three site types which impacted this species regeneration and productivity significantly.This suggests that black cottonwood will become increasingly rare across the region and may be restricted exclusively to riparian areas with wet moisture regimes.Western hemlock exhibited a variable response across all sites and disturbances but did not increase its abundance to any significant degree as its regeneration coefficients remained lower than all other species except black cottonwood.The modeled increases in regeneration and productivity were not sufficient enough to give this species a competitive advantage over the other species under climate change.

Interaction between resource availability and competition: productivity
The potential change in productivity under climate change is a large unknown for all of the study species; however, for lodgepole pine, productivity across BC under climate change has been predicted by Wang et al. (2006) and O'Neill et al. (2008). O'Neill et al. (2008) predicted that lodgepole pine productivity (m 3 / ha) is expected to decline by 45-74% with lodgepole pine disappearing from portions of v www.esajournals.orgthe study region.Wang et al. (2006) predicted that productivity could increase in the region by 20-30% if seed source is optimized but would only increase by 6.9% under a 28C increase in temperature and may decrease by 9.1% with a 48C increase in temperature if seed source is not optimized.Interestingly in this study, lodgepole pine productivity was modeled to decline by 47% on dry sites but increase by 21-27% on the mesic and moist sites.This study illustrates that variations in edaphic conditions and disturbance types may affect productivity highlighting the need to incorporate edaphic variation in predictions of climate change impacts on species productivity.

The role of disturbance in increasing the sensitivity of species
Climate induced changes in establishment, growth, competitive ability and disturbance frequency and severity will all impact on species distributions and abundance at multiple scales which in turn will affect the ecosystems they define (Johnstone and Chapin 2006b).Structural and compositional adjustments of ecosystems to environmental change can either be delayed by community inertia or accelerated by large-or small-scale disturbances (Hofgaard 1997).Ecological resilience refers to the magnitude of disturbance that can be absorbed before the system is restructured with different controlling variables and processes (Gunderson et al. 2002).In this study, both species and ecosystem resilience were found to vary depending on climate, disturbance and site type.Changes in ecological resilience, determined by changes in stand composition, was found to be a timedependent process influenced strongly by natural disturbances which in turn was impacted by climate change affects on species regeneration and productivity that cascaded from the tree to the stand level.These findings highlight that a loss of resilience precipitated by changes in processes at the tree and stand level can have ramifications on species composition and community structure respectively (Wirth et al. 2008).
By far the biggest impact on stand composition across the site types was the occurrence of fire.Across all edaphic sites, fire resulted in an increase in lodgepole pine and trembling aspen which on site types 2 and 3 resulted in the maintenance of stand-level productivity.Fireadapted species such as lodgepole pine and trembling aspen are expected to expand in areas where future climate causes increased fire frequencies at the expense of other species such as subalpine fir (Bartlein et al. 1997).This process was apparent following the stand-level simulations where lodgepole pine increased in abundance on sites that were traditionally dominated by interior spruce and subalpine fir while the latter declined in abundance.Frequent fire occurrence is also expected to give fire-adapted conifers, such as lodgepole pine, a competitive advantage over many deciduous species (Bradshaw et al. 2000).This was not the case in this analysis where trembling aspen was able to gain the competitive advantage over every other species following fire on site type 1. Trembling aspen has been found to respond more positively to increasing fire severity and occurrence than its co-occurring boreal species (Johnstone and Chapin 2006a).Though aspen has the lowest seed survival rate compared to white spruce and lodgepole pine, it is able to release a larger number of seeds permitting high rates of seedling establishment following a fire (Johnstone and Chapin 2006a) while its re-sprouting ability allows it to dominate areas at the expensive of seed producing conifer species (Johnstone 2006).Dense suckering following fires can lead to early crown closure and the suppression of other species (Johnstone 2006).These responses were apparent in this analysis as fire increased abundance of aspen at the expense of the conifers.Lodgepole pine was able to compete to some degree due to its serotinous behavior (Nyland 1998); particularly on the mesic sites.Interestingly, trembling aspen was unable to outcompete interior spruce and lodgepole pine on the moist-rich sites which highlights the ability of these conifers to maintain their competitive advantage over aspen following fire under these edaphic conditions.

Conclusion
In this study, we considered the impacts that climate, disturbance, resource availability and competition will have on the response on the dominant tree species in the sub-boreal forests of northwest British Columbia.All species were modeled to incur changes with differential responses identified between species and between sites.The results highlight the need to consider local responses of species, not just regional or landscape responses that are commonly investigated by bioclimatic envelope models.Species were modeled to undergo the largest changes on dry sites with a shift to an increase in lodgepole pine and trembling aspen abundance and the loss of subalpine fir and decline in interior spruce.
The individualistic responses of species to climate, disturbance and resource availability in this study highlights the complex interactions that may occur under future climate change.The findings highlight the potential resilience that species may exhibit to climate change when biotic and abiotic factors are considered.As a result, model outcomes of species response to climate at larger spatial scales may not be representative of species response at smaller spatial scales and vice versa (Bugmann et al. 2000, Barrio et al. 2006).In this study, disturbance drove change which affected the resilience of the existing communities following fire; particularly on dry sites and to a lesser extent on mesic sites.On fresh to moist sites, resilience was maintained following all disturbances which suggest that changes in these forest communities may be a slow process that will allow for species to persist in this composition and structure under climate change even if fire occurs.These findings also highlight the need to consider management actions that reduce disturbance risk and promote healthy ecosystem functions in order to exploit the ecological resilience that seemingly exists in these systems (Noss 2001, Gunderson et al. 2002).Based on the results of this study, fire will be a large driver of change and change will be greatest on dry sites which may create the best opportunities for new ecosystems to form while moist sites may provide the best opportunities for conserving the ecosystem composition and structure that currently exists within the SBS zone under predicted climate change.

Fig. 1 .
Fig. 1.Sub-boreal spruce zone research area in the Bulkley Valley of Northwest British Columbia.Smithers Airport and Topley Landing weather stations are demarcated with triangles.

Fig. 2 .
Fig. 2. Average annual climate in the Bulkley Valley based on Smithers Airport and Topley Landing weather stations.

Fig. 3 .
Fig. 3. Model framework used to determine species establishment in TACA.(a) Habitat niche elements including new frost day parameter and expanded soil component.(b) phenology module depicting the relationship between temperature, frost and bud burst.Growing season terminates at the cessation of growing degree accumulation (species-specific) and/or by user defined date to account for declining photoperiod length (October 15 in this study).Frost events following termination of the growing season do not cause mortality/ damage except for lethal frost events as recorded in (a).Texture class defines the field capacity and available water storage capacity (AWSC) of a soil in mm/m of soil.Percolation parameter determines if soil losses groundwater when AWSC is exceeded or receives groundwater from upslope sites; for mesic sites, percolation is set to zero.

Fig. 4 .
Fig. 4. Methodological framework for linking TACA model to SORTIE-ND model to investigate the impacts of climate variability and change on species and ecosystem response.

Fig. 5 .
Fig. 5. Summary of measured stand compositions used to generate replicates for each site type.Average stand composition with standard errors of the proportion shown.Sample sizes for each site type: dry-poor (n ¼ 11); mesic-medium (n ¼ 40); moist-rich (n ¼ 20).

Table 1 .
Study area weather station characteristics and average climate conditions.