Topoedaphic and Forest Controls on Post-Fire Vegetation Assemblies Are Modified by Fire History and Burn Severity in the Northwestern Canadian Boreal Forest

Abstract

Wildfires, which constitute the most extensive natural disturbance of the boreal biome, produce a broad range of ecological impacts to vegetation and soils that may influence post-fire vegetation assemblies and seedling recruitment. We inventoried post-fire understory vascular plant communities and tree seedling recruitment in the northwestern Canadian boreal forest and characterized the relative importance of fire effects and fire history, as well as non-fire drivers (i.e., the topoedaphic context and climate), to post-fire vegetation assemblies. Topoedaphic context, pre-fire forest structure and composition, and climate primarily controlled the understory plant communities and shifts in the ranked dominance of tree species (***8% and **13% of variance explained, respectively); however, fire and fire-affected soils were significant secondary drivers of post-fire vegetation. Wildfire had a significant indirect effect on understory vegetation communities through post-fire soil properties (**5%), and fire history and burn severity explained the dominance shifts of tree species (*7%). Fire-related variables were important explanatory variables in classification and regression tree models explaining the dominance shifts of four tree species (R² = 0.43–0.65). The dominance of jack pine (Pinus banksiana Lamb.) and trembling aspen (Populus tremuloides Michx.) increased following fires, whereas that of black spruce (Picea mariana (Mill.) BSP.) and white spruce (Picea glauca (Moench) Voss) declined. The overriding importance of site and climate to post-fire vegetation assemblies may confer some resilience to disturbed forests; however, if projected increases in fire activity in the northwestern boreal forest are borne out, secondary pathways of burn severity, fire frequency, and fire effects on soils are likely to accelerate ongoing climate-driven shifts in species compositions.

1. Introduction

Wildfires are the most extensive stand-initiating disturbance in the northwestern Canadian boreal forest, typically recurring every 50–100 years [1,2]. When wildfires occur, they burn with varying intensities (energy release) in response to fire weather, topography, and fuel type, producing a range of burn severities. Burn severity is defined as the ecological impacts of fire on vegetation and soils [3,4]. Many boreal forest plants have adapted to repeated wildfires though traits such as resprouting or suckering, seed banking, or, in the case of some tree species, serotiny. Serotinous and semi-serotinous conifer tree species have cones that may open in response to and survive some heating, and retain some viable seeds in the canopy following wildfires. Through this mechanism, serotinous species can produce extensive seed rains from aerial seedbanks immediately following fire [1,5]. Wildfire burn severity has important implications for post-fire understory vegetation communities and recruitment of seedlings. Heating and combustion from wildfires kill some trees and may reduce the viability of seeds in aerial seedbanks (including those of serotinous species) beyond a threshold of fire intensity or if the duration of heating is extensive [6,7]. Variable combustion of organic soils provides diverse seedbeds for plants and trees, ranging from thick remnant organic layers to exposed mineral soils, and alters the composition and exposure of post-fire soil seed banks [8,9]. Some burning of organic soils promotes vegetative regeneration, but deep burning may damage roots and rhizomes, negatively affecting the capacity of resprouting species to regenerate following fires [5].

In many ecosystems, burn severity is a dominant and enduring control on post-fire understory vegetation assemblies [10,11,12] and seedling recruitment [11,13,14], influencing the resulting structure and composition of forests. Although burned sites in the boreal forest generally return to a mature forested stand structure within 100 years [15], researchers using remote sensing to examine the post-fire recovery of vegetation following wildfires have found different rates of revegetation amongst burn severity classes. Severely burned sites demonstrated the highest decline in vegetation immediately post-fire [16,17]. In the years following a wildfire, severely burned sites subsequently experienced the largest increases in vegetation, indicating either forest recovery or colonization of these sites by disturbance-favouring plants and trees [16,17]. In North American boreal forests, post-fire understory vegetation communities in black spruce (Picea mariana (Mill.) BSP.) [18], jack pine (Pinus banksiana Lamb.) [19], and mixed broadleaf and coniferous stands [20] are influenced by surface burn severity and depth of burn, in conjunction with the availability of seed sources and vegetative propagules. In these studies, colonizing species such as graminoids and annual forbs established themselves broadly in severely burned areas, whereas slow-growing lichens, evergreen shrubs, and higher species richness were more prominent in low severity and scorched areas [18,19,20]. Lower densities of recruitment of coniferous trees have been observed when sites burned severely and at short intervals [21,22], and increased proportions of early-successional tree species, such as jack pine and trembling aspen (Populus tremuloides Michx.) are associated with high severity burning [23,24]. The relative dominance of different species of trees and the density of post-fire forests are lasting legacies of boreal wildfire severity [21,23,25].

When burn severity is studied at a broader landscape scale, that is, across multiple forest types and wildfires, the effects of burn severity on post-fire vegetation communities and recruitment may be challenging to detect. Burn severity is correlated to pre-fire forest type and stand structure [8,26,27,28], potentially obscuring or explaining observed effects of burn severity on post-fire plants and trees. Studies of burn severity that encompass multiple forest types have identified topoedaphic and pre-fire forest conditions as the primary post-fire drivers of understory plant communities and site suitability for tree species [18,29,30,31,32], leading some researchers to characterize burn severity as a secondary “filtering” effect beneath the dominant landscape and climatological controls.

Ranges of burn severity and the relatively infrequent occurrence of large wildfires (≥200 ha) produce a mosaic of stand ages and patterns on the landscape, in regions with mixed- and high-severity fire regimes [33,34,35]. Wildfires interact with past burns, as previous fires and burn severity determine current fuels. Abnormally short fire frequencies are implicated in the dominance shifts of tree species [36], low stocking in post-fire forests [22], and even near-deforestation [37], with implications for forest resilience [38]. Furthermore, burn severity interacts with fire frequency, potentially reinforcing vegetation type conversions [39]. Wildfires are a weather and, therefore, climate-driven disturbance. Fires are expected to increase in size, frequency, and intensity (and therefore in severity) [40,41,42] in North America as the climate warms and severe fire weather increases [43]. The forests of the Canadian Northwest Territories provide an interesting opportunity to study the effects of extensive free-burning wildfires in an ecosystem with multiple dominant coniferous and broadleaf tree species, across a moisture gradient ranging from hydric to xeric. Given the ecologically important role and actively changing patterns of fire in the boreal forest, studies characterizing the relative importance of fire effects and fire history, and non-fire and climate drivers in determining post-fire vegetation assemblies and species composition shifts will provide insights into the trajectories of future forests.

This study describes post-fire vegetation communities and seedling recruitment across a broad range of topoedaphic vegetation classes and levels of burn severity, to identify direct and indirect drivers of these assemblies in the northwestern Canadian boreal forest. In support of this goal, our objectives were: 1. To characterize post-fire vegetation assemblies and recruitment of seedlings across burn severity and topoedaphic gradients; 2. To assess the relative importance of climate and pre-fire forests, burn severity and fire history, and post-fire soils to understory vegetation communities and shifts in the dominance of tree species; and 3. To identify direct and indirect effects of fires on post-fire vegetation, as well as drivers of shifts in the dominance of tree species in the post-fire cohort.

2. Methods

2.1. Study Area

Field sites were established in six, large, lightning-caused wildfires (14,000 to 700,000 ha) that burned in 2014 (Figure 1). The year of 2014 was an extreme fire season in the northwestern Canadian boreal forest region, with drought-driven wildfires burning a total area > 3 million ha [44]. The sampled fires burned in the Northwest Territories and Wood Buffalo National Park. The fire regime of this area is one of infrequent stand-initiating wildfires [45,46]. In the Canadian boreal forest, these large wildfires comprise a small fraction of the total number of fires, but they are responsible for the vast majority of the area that was burned [35].

The study area experiences long cold winters and short hot summers. Mean annual temperatures at the field sites ranged from −4.3 °C at the furthest north site to −1.8 °C at the furthest south [48,49]. Topography of the study area is minimal, consisting of level terrain in the southwestern part of the study area, on the boreal plain, and rolling granitic hills on the boreal shield in the northeast [41]. The forests of the study area are dominated by jack pine, black spruce, white spruce (Picea glauca (Moench) Voss), and trembling aspen. Important secondary tree species include eastern larch (Larix laricina (Du Roi) K. Koch), paper birch (Betula papyrifera Marsh.), and balsam poplar (Populus balsamifera L.) [48]. There is also a substantial wetland (chiefly peatlands) component to the region. Peat-forming wetlands may form extensive complexes and cover approximately a third of the total area [50]. Although the study area falls within the discontinuous permafrost zone of Canada [51], no field sites had frozen active layers in the top metre of soil.

2.2. Field Methods

We sampled 51 field sites one year post-fire and resampled 30 sites three years post-fire. The sites were selected using a stratified random sample that was evenly distributed across high-, moderate-, and low-burn severity classes. The mapped burn severity was produced using an initial assessment of a differenced normalized burn ratio (dNBR) image [52], classified with thresholds developed by Hall et al. [53]. Field sites were >100 m and ≤2 km from roads. More isolated sites were also opportunistically accessed by helicopter. The field sites accessed by helicopter were located in order to capture the locally available range of burn severity and topoedaphic vegetation communities (ecosites), ensuring that each sampled site offered a distinct combination of severity and vegetation type. Field sites were positioned in an area of homogenous burn severity, topoedaphic setting (upland or wetland), and dominant vegetation that extended ≥60 m in any direction. The site moisture (from hydric to xeric) and ecosite categories were classified according to Beckingham and Archibald [54]. Ecosites were generalized into the dominant topoedaphic vegetation classes of open wetland, treed wetland, upland spruce, upland mixedwood, and upland jack pine (from wettest to driest). All the sampled wetlands were peat-forming wetlands (peatlands). Plot centres were recorded with a differential GPS unit. The mean distance between the plot centres of all the field sites was 170 km, with a minimum distance of 103 m.

When sampling one year after the fire, the sample plots were 30 × 30 m, with two 30-m transects oriented in the cardinal directions, crossing at the plot centre. A detailed figure of the plot layouts used for field sampling is included in Appendix A: Figure A1. Compositions of tree species, percent overstory mortality due to fire, stem density (stems ha⁻¹), and basal area (m² ha⁻¹) of mature trees in the pre-fire stand were measured at this time for 32 trees ≥ 3 cm diameter at breast height (DBH) using the point-centered quarter method [55,56] at eight evenly-spaced points along the two transects. In very low stem-density areas (i.e., open wetlands), a variable-radius circle plot with a minimum radius of 15 m was used to sample overstory trees. Pre-fire understory stem densities of seedlings and saplings (stems ha⁻¹) were measured using 3-m radius plots at the endpoints of each transect. The number of understory density plots sampled ranged from one to four, depending on the density and evenness of the seedling and saplings.

We collected basal sections from fire-scarred trees to determine the time since the stand origin (TSO) and time since the last fire (TSLF) at each plot. If no scarred trees were identified nearby, a section of a mature dominant tree was sampled. Some open wetlands (fens) had no trees. Samples were sanded and digitally scanned, and annual growth rings and fire scars were dated in CooRecorder [57].

Burn severity was measured in 10 × 10 m subplots at the four corners of each plot. Surface burn severity was measured using the surface Burn Severity Index (BSI) [58]. BSI values range from zero (unburned) to four (ash, mineral soil exposed) using classes defined by Dyrness and Norum [59]. Overstory burn severity was measured using the Canopy Fire Severity Index (CFSI) [60]. CFSI classes range from zero (no tree mortality) to six (no primary branches remaining, pole charring occurred). The percent cover of each BSI and CFSI severity class was estimated within the four subplots, and final values of the two severity metrics were calculated using area-weighted means of each class value, and then averaged for each field site.

We measured the post-fire organic soil depth (cm; up to a maximum of 10 cm) at the inner corners of the same subplots used for estimates of severity and seedling density. The soil cores (13.5 cm in depth, 5.5 cm in diameter) were taken one year post-fire at the plot centre and inner corners of the southwest and northeast subplots, as well as at a complementary set of neighbouring unburned control sites (n = 12) representing unburned examples of all sampled vegetation communities. Cores were inserted to a minimum depth of 8.5 cm and the soil samples were separated into organic and mineral horizons; the three samples from each site were pooled by the horizon. If mineral soil was not present in the top 13.5 cm of the soil profile, it was not collected. Soils were oven-dried and the physicochemical properties of both organic and mineral samples were measured in the lab. These properties were: pH, electrical conductivity (EC; mS cm⁻¹), percent total nitrogen (N), percent total carbon (C) measured by loss on ignition, calcium (Ca; mg kg⁻¹), potassium (K; mg kg⁻¹), magnesium (Mg; mg kg⁻¹), and sodium (Na; mg kg⁻¹). The percentages of sand, silt, and clay in mineral soils were also measured. Measurements from the two pooled horizons from each site were combined using sums weighted by the mean proportion of the core occupied by each horizon.

Estimates of percentage cover of understory vascular plant species were made one year post-fire in five 1 × 1 m plots per field site. Vegetation plots were located at the plot centre and at the inner corners of subplots. Species were identified according to Moss [61] and Cody [62], and the estimated percentage cover for each species was summed across the five plots and scaled to sum to 100%. Carex spp. and Salix spp. were distinguished for counts of species richness but were not identified beyond genus for ordination or indicator species analyses (vegetation analysis explained in detail in Section 2.3).

The density of seedling recruitment was measured one year post-fire (2015), and subsequently re-measured three years post-fire (2017) in 30 forested sites (excluding open wetlands). Initial measures of seedling density were made in the 10 × 10 m subplots in 2015. In 2017, seedling density was re-measured using a 2-m wide 35-m long belt transect that was oriented north-south, crossing the original plot centre at 17.5 m. Belt-transect length varied by seedling and sapling size classes. Seedlings that were 0–10 cm were counted for the first 10 m of the transect (area 20 m²) and seedlings that were 10–50 cm were counted for the first 20 m (area 40 m²). Seedlings > 50 cm and saplings (live trees > 1.33 m with a DBH < 3 cm) were counted for the entire transect length. In cases of very uneven seedling density, transects of all size classes were extended to better represent the actual composition and density. This set of resampled sites excluded non-forested open wetlands (n = 11) and inaccessible sites (no helicopter or road access, n = 7). A further three sites were abandoned due to subsequent disturbances. The two datasets were combined and the latest available seedling density measurement for each site was used.

We calculated site climatic variables that described the average heat load and moisture stress from 30-year normals (1981–2010) of PRISM climate data [63] downscaled to local elevation [64] using bilinear interpolation and elevation adjustment in ClimateWNA [49]. The climatic moisture deficit (CMD; mm) was calculated as the sum of the monthly difference between Hargrave’s atmospheric evaporative demand and monthly precipitation. Annual heat-to-moisture index (AHM) was calculated as the scaled ratio of mean annual temperature and mean annual precipitation [49].

2.3. Analysis

All statistical analyses were conducted in R [65]. The variance of burn severity explained by the topoedaphic vegetation classes was assessed using a linear mixed-effects model with a random term of the fire name, fitted in the lme4 [66] and lmerTest [67] packages. We examined all model residuals and found them to be normally distributed. The statistical significance of fixed effects was estimated using an analysis of variance (ANOVA) with Type II sums of squares and a Satterthwaite approximation of degrees of freedom [68]. We conducted post-hoc comparisons of least-squares means with a Tukey test for multiple comparisons in the lsmeans package [69]. Species richness and Shannon diversity index of understory vascular plant communities of each site were calculated in the vegan package [70]. Bray-Curtis dissimilarities between understory vegetation communities were ordinated using non-metric multidimensional scaling (NMDS). We fitted vectors of environmental variables to the NMDS axes and assessed the goodness of fit (R²) of these relationships also using vegan (

Table 1. Significant (* p ≤ 0.05) explanatory environmental variables fitted to nonmetric multidimensional scaling (NMDS) of the understory vegetation community data (Figure 2). Unitless variables are identified with a hyphen in the Units column.

Environmental Variable	Abbreviation	Units	Mean	Range
Basal Area	BA	m2 ha−1	10.8	0.00–53.39
Burn Severity Index	BSI	-	2.39	0.54–4.00
Electrical conductivity of soil	EC	mS cm−1	0.73	0.05–3.53
Organic soil depth	OSD	cm	4.7	0–10
Percentage sand in mineral soil	% Sand	%	43.9	0–95
pH	pH	-	6.29	3.21–8.12
Potassium	K	mg kg−1	411.1	74.1–1148.4
Site moisture	Moisture	-	-	Xeric–Hydric
Sodium	Na	mg kg−1	137.8	38.52–494.48
Time since last fire	TSLF	year−1	58	9–151
Total carbon	Total C	% mass	21.9	0.61–52.8
Total nitrogen	Total N	% mass	0.81	0.18–2.66
Total stem density of overstory and understory trees	Density	stems ha−1	5822	0–29,012

). Indicator species for each topoedaphic vegetation class were identified from understory vascular plant assemblies using 1000 permutations of a multi-level pattern analysis in the indicspecies package [71]. We assessed the influence of burn severity on soil properties when controlling for topoedaphic vegetation class (as a proxy for pre-fire site conditions) using an ANOVA of multivariable linear mixed-effects models with a random term of the fire name. Once again, model residuals were examined for normality. We employed Type II sums of squares where interactions between independent variables were not significant. We applied a Type III ANOVA if there were significant interactions between independent variables. The same approach was used to assess the influence of burn severity, TSLF, and topoedaphic vegetation class (pre-fire conditions) on the Shannon diversity index and seedling density. We used comparisons of least-squares means with a Tukey test for multiple comparisons to assess significant differences in species dominance shifts and seedling density between different topoedaphic vegetation classes.

Differences in the pre-fire and post-fire cohorts of trees were examined using compositional data. Overstory basal areas and total (understory and overstory) stem densities of each dominant tree species were converted to proportions relative to the absolute basal area and stem density for each site. Seedling counts were also converted to proportions by species, and these proportions, or compositions, were transformed with a centred log-ratio using the compositions package [72]. We then used paired t tests to identify statistical differences in the pre-fire and post-fire composition of trees, by species. We compared the natural logarithm (log_e) of seedling density in sites that experienced very short fire return intervals to that of sites experiencing more typical fire return intervals with a Wilcoxon signed-rank test. Significant differences in log_e seedling density between topoedaphic vegetation classes were also tested using a Tukey test of least-squares means.

To examine shifts in the relative importance or dominance of tree species we calculated fractional ranks of pre-fire overstory tree species proportions of jack pine, white spruce, black spruce, trembling aspen, and all other tree species combined, by basal area. The most prevalent species received a rank of 1 and the least dominant (or absent) species received a rank of 5. In the case of ties, ranks were split between species, so that total rank values always summed to 15. We chose to use the pre-fire basal area rather than the number of stems as a measure of dominance as we felt that stem density did not adequately capture the potential fecundity and relative importance of less-common but large trees in mixedwood stands (e.g., white spruce). Because the basal area of trees established after fire represents only a small fraction of the pre-fire measure, post-fire tree species proportions were assigned fractional ranks by seedling stem density. The pre- and post-fire fractional rank scores of each species at each site were differenced to characterize shifts in tree species dominance in the post-fire cohort, producing a matrix of shifts in ranked dominance for each species by site. Rank shifts of near-zero indicated minimal change in the species’ prevalence in the post-fire cohort, whereas negative values indicated a decrease and positive values indicated an increase in ranked dominance. Rank shift data had a theoretical range of −4 to 4. Analyses using the shift in rank dominance data were only performed for the four dominant tree species, excluding the combined “other” category. This application of fractional rank shifts characterizes the proportional change of tree species dominance relative to all tree species present in the community, rather than considering a single species at a time (as is the case with ratio data), and offers a normally distributed variable for analysis of dominance changes.

Subsequently, we assessed the relative importance of three groups of variables in the categories “Soils”, “Site”, and “Fire” to understory vegetation community dissimilarities and shifts in the dominance of tree species using variance partitioning. Soils were represented by post-fire soil properties, whereas Site category variables were pre-fire forests, topoedaphic context, and recent spatial climate averages. The Fire category included burn severity and fire history variables (

Table 2. Environmental variables incorporated in the explanatory variance partitioning of understory vegetation community dissimilarities and shifts in tree species dominance. Correlated soil properties collapsed with a principal components analysis for inclusion in variance partitioning are indicated with a ^†. Unitless variables are identified with a hyphen in the Units column.

Environmental Variable	Units	Mean	Range	Variance Partitioning Category
Calcium †	mg kg−1	13,648.6	217.1–60,815.4	Soils
Electrical conductivity †	mS cm−1	0.73	0.05–3.53	Soils
Magnesium †	mg kg−1	1461.9	22.0–5191.3	Soils
Percentage sand in mineral soil	% mass	44	0–95	Soils
Percentage silt in mineral soil	% mass	14	0–51	Soils
pH	-	6.29	3.21–8.12	Soils
Potassium †	mg kg−1	411.1	74.1–1148.4	Soils
Sodium †	mg kg−1	137.8	38.5–494.5	Soils
Total carbon †	% mass	21.9	0.61–52.8	Soils
Total nitrogen †	% mass	0.81	0.02–2.66	Soils
Absolute stem density of overstory and understory trees	stems ha−1	5822	0–29,012	Site
Annual Heat-Moisture Index	-	20.32	16.5–23.6	Site
Black spruce basal area	m2 ha−1	2.37	0–29.07	Site
Climatic Moisture Deficit	mm	191	171–214	Site
Jack pine basal area	m2 ha−1	5.61	0–51.38	Site
Site moisture	-	-	Xeric–Hydric	Site
Trembling aspen basal area	m2 ha−1	0.87	0–6.22	Site
Total overstory basal area	m2 ha−1	10.8	0–53.4	Site
White spruce basal area	m2 ha−1	1.55	0–36.53	Site
Burn Severity Index	-	2.38	0.54–4	Fire
Percentage overstory mortality	%	89	6–100	Fire
Post-fire organic soil depth	cm	4.7	0–10	Fire
Time since last fire	year−1	58	9–151	Fire
Time since stand origin	year−1	104	9–237	Fire

). All measured and downscaled environmental variables were considered for inclusion in variance partitioning models. If the variables were highly correlated (Spearman’s |ρ| ≥ 0.7) one explanatory variable of the pair was selected for inclusion in the model. Several highly correlated soil properties were decomposed using a principal components analysis (PCA; Table 2). Sites with incomplete data were removed (n = 5), and explanatory variables were standardized before variance partitioning. The significance (α = 0.05) of the unique variation explained by each group of environmental drivers (Soils, Site, and Fire) was tested using distance-based redundancy analysis, also in the vegan package.

Finally, we fit explanatory classification and regression trees (CARTs) to shifts in the ranked dominance of each tree species derived from ranked proportions of the pre-fire basal area and post-fire stem density, using the tree package [73]. The regression trees were constrained by requiring a minimum of five field sites per node, and a minimum within-node deviance of 0.05. We intentionally excluded pre-fire basal area and stem density of any tree species as regression tree predictor variables in order to learn about secondary climatic, soil, and burn severity effects on the dominance shifts of tree species. The same suite of environmental, burn severity and fire history, and climate variables were included as potential predictors of shifts in tree species dominance for each species’ CART model (

Table 3. Environmental, climate, burn severity, and fire history variables included in the classification and regression tree models of the dominance shifts of tree species. Unitless variables are identified with a hyphen in the Units column.

Variable	Units	Mean	Range
Annual Heat-Moisture Index	-	20.3	16.5–23.6
Burn Severity Index	-	2.38	0.54–4
Canopy Fire Severity Index	-	2.5	0–6
Climatic Moisture Deficit	mm	191	171–214
Electrical conductivity	mS cm−1	0.73	0.05–3.53
Percentage overstory mortality	%	89	6–100
Percentage sand in mineral soil	%	44	0–95
pH	-	6.29	3.21–8.12
Post-fire organic soil depth	cm	4.7	0–10
Total nitrogen	%	0.81	0.02–2.66
Time since last fire	year−1	58	9–151
Time since stand origin	year−1	104	9–237
Wetland	-	-	Upland or Wetland

).

3. Results

A broad range of burn severity was represented in the field sites (Appendix A: Figure A2). The BSI values of field sites ranged from 0.5 to 4, the CFSI values from 0 to 6, and percent overstory mortality ranged from 6.25% to 100% [28]. The surface (BSI) and overstory (CFSI) burn severity were statistically related to topoedaphic vegetation classes (ANOVA, *** p < 0.001 and * p = 0.02, respectively;

Table 4. Multivariable linear mixed-effects models describing surface Burn Severity Index (BSI), Canopy Fire Severity Index (CFSI), time since last fire (TSLF), topoedaphic vegetation classes (TVC), and ecological outcomes of seedling density and Shannon Diversity Index of understory vascular plant communities. The effect size (F) and significance (p) of terms are tested with Type II sums of squares where there are no significant interactions, and Type III sums of squares in the presence of a significant interaction. Significance of intendent variables is signified as follows: *** p ≤ 0.001, ** p ≤ 0.01, * p ≤ 0.05.

Multivariable Linear Mixed-Effects Model	ANOVA Sums of Squares	Degrees of Freedom	Independent Variable	Sums of Squares	F	p
BSI = TVC + (1|Fire Name)	II	4	TVC	20.39	14.94	*** < 0.001
CFSI = TVC + (1|Fire Name)	II	4	TVC	33.40	3.36	* 0.02
Diversity = TVC × BSI + (1|Fire Name)	III	4	TVC	2.51	3.09	* 0.03
		1	BSI	0.15	0.72	0.40
		4	TVC × BSI	2.64	3.26	* 0.02
Diversity = TVC × TSLF + (1|Fire Name)	III	4	TVC	2.08	2.52	0.06
		1	TSLF	0.00	0.01	0.91
		4	TVC × TSLF	2.89	3.51	* 0.02
Density = TVC + BSI + (1|Fire Name)	II	4	TVC	18.01	3.49	* 0.02
		1	BSI	5.62	4.35	* 0.04
Density = TVC + TSLF + (1|Fire Name)	II	4	TVC	43.94	7.84	*** < 0.001
		1	TSLF	0.59	1.14	0.29

), but overstory mortality was not. Post-hoc comparisons of least-squares means with a Tukey test confirmed some statistical differences in burn severity amongst topoedaphic vegetation classes for BSI (* p ≤ 0.05; Appendix A: Figure A2). Surface burn severity was lowest in open wetlands and highest in jack pine uplands. All other topoedaphic vegetation classes had BSI values that were similar to one of these two groups. The differences in least-squares means of CFSI between topoedaphic vegetation classes were not significant at α = 0.05 (Appendix A: Figure A2).

The richness of the understory vascular plant species sampled at the field sites ranged from three to 20. Both understory vegetation community diversity and seedling density were statistically related to topoedaphic vegetation classes (Table 4; Appendix A: Figure A2). Interactions between TSLF and topoedaphic vegetation classes, and BSI and topoedaphic vegetation classes significantly explained the variability in understory plant diversity (Type III ANOVA, * p ≤ 0.02; Table 4). The density of seedlings was significantly explained by both topoedaphic vegetation classes and BSI (Type II ANOVA, p < 0.04), but not by TSLF or CFSI, or by the interactions between topoedaphic vegetation classes and these two metrics (Table 4; Appendix A: Figure A2). All sites with zero seedling establishment were open wetlands (n = 7). Of those sites that experienced some regeneration, seedlings ha⁻¹ ranged from 25 to >75,000. The natural logarithm of the density of seedlings was significantly lower in open wetlands, treed wetlands, and upland spruce, and higher in upland jack pine and upland mixedwood topoedaphic vegetation classes (comparison of least-squares means, with a Tukey test, α = 0.05). The post-fire seedling density was statistically greater in sites that experienced >17 years between fires (Wilcoxon signed-rank test, * p = 0.02).

The two-dimensional NMDS of understory vascular plant communities had a stress of 0.20 (Figure 2). Similarity of understory species communities was primarily related to the physicochemical properties of the soil; however, pre-fire forest structural characteristics of basal area and absolute stem density of overstory and understory trees were also influential (Figure 2). BSI was also statistically related to understory species community dissimilarity, and TSLF was nearly significant (p = 0.052, 999 permutations). Although soil properties were explained by topoedaphic vegetation classes; the organic soil depth, total nitrogen, total carbon, potassium, calcium, and magnesium were also statistically (Type II ANOVA; α = 0.05) related to BSI when controlling for the effect of topoedaphic vegetation class. Therefore, some soil properties were affected by fire (Appendix A:

Table A1. Descriptive multivariable linear mixed-effects models explaining post-fire organic soil depth and soil chemical properties as a function of topoedaphic vegetation classes (TVC) and observed surface burn severity, represented by the Burn Severity Index (BSI). The statistical significance of independent predictor variables in explaining soil properties was determined using an ANOVA considering Type II sums of squares. Significance of independent variables to soil properties is signified as follows: *** p ≤ 0.001, ** p ≤ 0.01, * p ≤ 0.05. Some interactions are nearly significant (α = 0.1).

Soil Property	Multivariable Linear Mixed-Effects Model	Independent Variable	Degrees of Freedom	Sums of Squares	F	p
Organic Soil Depth (cm)	OSD = TVC + BSI + (TVC × BSI)|Fire Name	TVC	4	91.50	10.10	*** < 0.001
		BSI	1	30.68	13.54	*** < 0.001
		TVC × BSI	4	8.80	0.97	0.43
Total Nitrogen (%)	TN = TVC + BSI + (TVC × BSI)|Fire Name	TVC	4	2.34	8.50	*** < 0.001
		BSI	1	0.81	11.77	** 0.001
		TVC × BSI	4	0.26	0.96	0.44
pH	pH = TVC + BSI + (TVC × BSI)|Fire Name	TVC	4	0.21	0.06	0.99
		BSI	1	002	0.02	0.88
		TVC × BSI	4	0.85	0.25	0.91
Electrical Conductivity	EC = TVC + BSI + (TVC × BSI)|Fire Name	TVC	4	3.16	2.08	0.10
		BSI	1	0.32	0.84	0.36
		TVC × BSI	4	0.59	0.39	0.82
Total Carbon (%)	TC = TVC + BSI + (TVC × BSI)|Fire Name	TVC	4	1697.20	13.20	*** < 0.001
		BSI	1	382.75	11.90	** 0.001
		TVC × BSI	4	91.91	0.71	0.59
Sodium	Na = TVC + BSI + (TVC × BSI)|Fire Name	TVC	4	121,344	6.20	*** < 0.001
		BSI	1	26,453	5.41	* 0.02
		TVC × BSI	4	41,222	2.11	0.09
Potassium	K = TVC + BSI + (TVC × BSI)|Fire Name	TVC	4	59,853	0.60	0.67
		BSI	1	35,436	1.41	0.24
		TVC × BSI	4	226,446	2.26	0.08
Calcium †	Ca = TVC + BSI + (TVC × BSI)|Fire Name	TVC	4	2,000,913,013	4.00	** 0.007
		BSI	1	672,020,248	5.37	* 0.02
		TVC × BSI	4	656,747,510	1.31	0.28
Magnesium	Mg = TVC + BSI + (TVC × BSI)|Fire Name	TVC	4	23,579,801	5.54	*** 0.001
		BSI	1	5,452,454	5.13	* 0.03
		TVC × BSI	4	6,042,409	1.42	0.25

^† Calculated using Kenward–Roger approximation of degrees of freedom, due to the mathematical failure of Satterthwaite’s approximation.

). Topoedaphic vegetation classes tended to occupy characteristic areas of ordination space, but there was some overlap between the normal confidence ellipses of classes. Upland mixedwood and upland jack pine groups were especially intermingled (Figure 2a), and mixedwood communities occurred in a sub-region of the broader environmental space occupied by jack pine. Similar patterns are identifiable in the post-fire understory indicator species of each topoedaphic vegetation class (

Table 5. Significant (* p ≤ 0.05) indicator species identified using multi-level pattern analysis within six topoedaphic vegetation classes. The indicator species uniquely associated with one group are indicated with a ^†.

Vegetation Group
Open Wetland	Treed Wetland	Upland Mixedwood	Upland Jack Pine	Upland Spruce
Potentilla palustris (L.) Scop. †	Rubus chamaemorus L. †	Viburnum edule (Michx.) Raf. †	Cornus canadensis L.	Vaccinium uliginosum L. †
Betula glandulosa Michx. †	Vaccinium caespitosum Michx. †	Cornus canadensis L.	Geranium bicknellii Britt.	Geocaulon lividum (Richards.) Fern. †
Epilobium palustre L. †	Vaccinium oxycoccos L. †	Geranium bicknellii Britt.	Rosa acicularis Lindl.	Ledum groenlandicum Oeder
Myrica gale L. †	Rubus arcticus L.	Rosa acicularis Lindl.	Linnaea borealis L.	Equisetum scirpoides Michx.
Rubus arcticus L.	Ledum groenlandicum Oeder	Linnaea borealis L.	Elymus innovatus Beal	Arctostaphylos rubra (Rehder & Wils.)
Carex L. spp.	Equisetum scirpoides Michx.	Elymus innovatus Beal		Rosa acicularis Lindl.
Salix L. spp	Arctostaphylos rubra (Rehder & Wils.) Fern.			Linnaea borealis L.
	Carex L. spp.			Elymus innovatus Beal
	Salix L. spp.			Carex L. spp.
				Salix L. spp.

). All topoedaphic vegetation classes had unique significant indicator species, with the exception of jack pine uplands, which shared all significant indicator species with the upland mixedwood group, and some with the upland spruce group (Table 5). Potentilla palustris (L.) Scop., Betula glandulosa Michx., Epilobium palustre L., and Myrica gale L. were unique indicator species of open wetlands. Treed wetlands had unique indicator species of Rubus chamaemorus L., Vaccinium caespitosum Michx., and Vaccinium oxycoccos L. Viburnum edule (Michx.) Raf. was a unique indicator species of upland mixedwood sites. Furthermore, upland mixedwood sites shared significant indicator species of Cornus canadensis L., Geranium bicknellii Britt., Rosa acicularis Lindl., Linnaea borealis L., and Elymus innovatus Beal with upland jack pine sites. Vaccinium uliginosum L. and Geocaulon lividum (Richards.) Fern. were unique indicator species in upland spruce communities.

The dominance of pre-fire and post-fire tree species (represented by log-ratios of compositions of basal area and stem density) were significantly different for jack pine, black spruce, and white spruce (Paired t tests, Bonferroni-corrected * p ≤ 0.04). When pre-fire dominance was characterized using total stem density (Appendix A: Figure A3), post-fire tree species compositions were significantly different for black spruce and trembling aspen (Paired t tests, Bonferroni-corrected * p ≤ 0.02). Having confirmed significant differences between pre-fire and post-fire tree species compositions, we examined the rank shifts in dominance of tree species in order to capture directionality of species-specific changes. Jack pine both increased and decreased in dominance in the post-fire cohort, but the slim majority of sites were neutral (−0.5 to 0.5 shift in rank; 39% of sites). Furthermore, not all plots burned with completely stand-initiating lethal wildfires (Figure 3). Of the sites that experienced declines in the ranked dominance of jack pine, 41% had some surviving jack pine basal area post-fire (Figure 3). Aspen dominance increased in the post-fire cohort, with 54% of sites gaining 1 or more ranks of dominance post-fire, and no sites declining by <−0.5 of a rank (Figure 3). Both varieties of spruce primarily demonstrated no change or declines, in the post-fire cohort (57% of sites were neutral and 37% showed a decrease for black spruce; 69% of sites were neutral and 29% showed a decrease for white spruce). Of those sites with declines in the dominance of black spruce, only 9% (n = 2) had incomplete mortality of black spruce trees (Figure 3). No sites demonstrating declines in trembling aspen or white spruce had live individuals of these species post-fire (Figure 3). Declines and increases in tree species dominance were significantly related to topoedaphic vegetation classes. Increases in jack pine dominance were especially associated with wetlands and spruce uplands, whereas jack pine declines occurred in upland communities where jack pine was already established, especially mixedwood stands where the suckering of trembling aspen was prevalent (Appendix A: Figure A4; Tukey test of least-squares means, * p ≤ 0.05). Although increases in aspen dominance were more pronounced in uplands (Wilcoxon signed-rank test, ** p = 0.009), there were no statistical differences in aspen dominance shifts between topoedaphic vegetation classes (Tukey test of least-squares means; Figure A4). Decreases in black spruce dominance were the most pronounced in upland spruce sites, whereas black spruce dominance was largely stable in treed wetlands and other vegetation classes (Tukey test of least-squares means, * p ≤ 0.05). There were no significant differences in post-fire changes in white spruce dominance between topoedaphic vegetation classes (Figure A4).

When representing the variance in understory vegetation and tree species dominance shifts explained by soils, we found that many soil properties were highly correlated. To address this, we decomposed the correlated soil physicochemical properties using PCA, and included only the first principal component (PC1) as an explanatory variable in the variance partitioning (Table 2). We chose to retain percent sand and PC1, and organic soil depth and site moisture despite high correlations (ρ = 0.8) between these two pairs, as these variables characterized important elements of the three environmental driver groups. Post-fire soils (Soils); pre-fire forests, topoedaphic context, and climate (Site); and burn severity and fire history (Fire) together explained 28% of the variance in understory vegetation communities, and 33% of the variance in the dominance shifts of tree species (environmental variables included in model reported in Table 2; Figure 4). There was a substantial shared variance explained between Soils, Site, and Fire. Overall, Site explained the largest portion of the variance in post-fire vegetation communities (8%) and tree species dominance shifts (13%; Figure 4). Soils significantly explained 5% of the variance in understory vegetation but did not significantly explain tree species dominance shifts. Conversely, Fire was of substantial importance to tree species dominance shifts (7% of variance) but did not significantly explain post-fire vegetation communities (p = 0.08; Figure 4).

Classification and regression trees of the dominance shifts of tree species had R² values ranging from 0.65 to 0.43. Jack pine dominance increased in the post-fire cohort where the total soil N was ≥0.48% and decreased in the post-fire cohort in stands that experienced partial mortality (Figure 5). These sites were typically mixedwood stands, and often had some remaining live basal area of jack pine trees, suggesting that these declines in the dominance of the post-fire cohort do not necessarily indicate persistent shifts away from jack pine dominance, although aspen suckering outpaced the establishment of pine seedlings (Appendix A: Figure A3). Trembling aspen increased in dominance in nearly all plots, but increases were somewhat limited in lightly burned plots and plots with higher N availability, both of which tend to be characteristic of wetlands (Figure 5). Changes in black spruce dominance were neutral in young stands where black spruce was essentially absent pre-fire and in wetlands (TSO < 80.5). Declines in black spruce dominance were augmented in moderate aged (TSO < 103.5) uplands and stands experiencing severe canopy burning (CFSI ≥ 4; Figure 5). White spruce dominance increased slightly in sites with low N availability (Total N < 0.11%) and lower moisture deficits (CMD < 191.5). White spruce dominance declined in historically drier sites, especially in those sites that experienced some canopy involvement in the fire (CFSI ≥ 2.65; Figure 5).

4. Discussion

4.1. Post-Fire Vegetation Communities of Vascular Plants

Post-fire understory vegetation communities were primarily explained by site conditions, but burn severity and fire history had significant secondary effects. Topoedaphic vegetation classes occupied distinct areas of the environmental space and had characteristic indicator species, with the exception of jack pine uplands and mixedwood stands. These indicator species were identifiable for topoedaphic vegetation classes despite recent disturbances, underscoring the importance of non-fire drivers to post-fire vascular vegetation. Understory vegetation communities of mixedwood forests appear to predominantly occur in a subset of environmental conditions that are also characterized by jack pine forests. Site moisture, climate, and pre-fire forest structure and composition, all of which were dominant drivers of understory vegetation assemblies, are independent of fire effects. Therefore, these communities may be somewhat robust to disturbance from fire, as they are primarily controlled by non-fire drivers.

Although the climate of the boreal forest is changing, changes to topography through background rates of uplift (1.3 cm year⁻¹) [74] and erosion (~0.005–0.0005 cm year⁻¹) [75] reported in parts of the biome are slow, relative to the velocity of climate change (~1 km/year⁻¹ in the boreal forest) [76]. Hydrological feedbacks may also reinforce the persistence of some features, such as peatlands, in the face of drying and warming [77,78]. The persistence of these topoedaphic drivers over time, in the face of ongoing climate change, should encourage the re-establishment of understory vegetation communities following fires. As understory vegetation communities have a substantial influence on below-canopy light availability and nutrient cycling, this vegetation layer’s persistence may, in turn, reinforce the similarity of post-fire communities to pre-fire conditions, with implications for seedling establishment [79,80,81]. Conversely, regions of the boreal forest experiencing rapid topoedaphic changes due to permafrost thaw and thermokarst formation may be more susceptible to shifts in vegetation communities as these dominant drivers undergo substantial short-term changes [82].

Fire affects post-fire vegetation assemblies, both directly through surface burn severity and time since last fire, and indirectly via fire-mediated changes to soil properties. Despite the overarching importance of fire-independent site characteristics, burn severity and fire history were associated with understory vegetation community dissimilarity. Additionally, these variables also had significant effects on diversity, when controlling for topoedaphic vegetation class. Furthermore, post-fire soils alone explained a substantial portion of the variance in understory vegetation communities and several soil properties that were affected by wildfire, as expected in boreal soils [83]. Despite the importance of non-fire site drivers to understory vegetation communities, changes in burn severity or fire return intervals will likely translate to shifts in understory vegetation communities through these secondary pathways. Increases in burn severity, in particular, may lead to lasting, directional compositional changes in understory species assemblies [31]. Where understory vegetation communities exist at the boundaries of their ideal environmental space (or in areas of overlap with other topoedaphic vegetation classes), the effect of burn severity and fire history may be more apparent, and potentially override fire-independent controls.

4.2. Post-Fire Shifts in Tree Species Dominance

Topoedaphic vegetation classes significantly explained total seedling density, but surface burn severity significantly interacted with these conditions, likely indicating an effect of seedbed availability on recruitment. Climate, pre-fire forest composition and structure, and site moisture were important variables in explaining post-fire shifts in the dominance of tree species, and in which sites such shifts occurred. Post-fire seedling density was significantly lower in wetlands and spruce-dominated sites and highest in jack pine and mixedwood uplands. Although the post-fire recruitment was lowest in open wetlands, some forested sites that reburned with very short intervals between stand-initiating fires (≤16 years) also experienced near-failures in the recruitment of all tree species, including those that increased in overall dominance post-fire, and had significantly lower seedling densities than all other sites. Although we observed declines in the post-fire dominance of both spruce species, black spruce maintained its dominance in treed wetlands (peatlands) and white spruce dominance was stable in sites with a lower climatic moisture stress. Furthermore, increases in the dominance of aspen in the post-fire cohort were least prevalent in wetlands. The drivers of species persistence and types of sites where spruce species retained dominance, despite the broader neutral or declining trend across sites, reflect the importance of topoedaphic and climatological drivers to post-fire tree species shifts. Local variability in site moisture may offer refugia from climate change for both of these species, in a landscape with limited topography [78,84].

Despite the importance of climate and pre-fire forests to seedling recruitment, burn severity and fire history had detectable and important effects on post-fire shifts in the dominance of tree species. Canopy fire severity was implicated in both positive and negative shifts in tree species dominance for all conifer species. Black spruce dominance decreased in uplands where stand-initiating wildfires occurred at a frequency of fewer than ~100 years between fires. The dominance of both spruce species was reduced in sites that experienced high-severity crown fire. Surface burn severity was the primary driver of post-fire increases in dominance of trembling aspen, and likely had further indirect effects on the post-fire dominance of white spruce, aspen, and jack pine through nutrient availability due to organic soil combustion and heating. Increases in burn severity and combustion in sites that tended to protect tree species that were susceptible to declines in dominance (i.e., peatlands) or in sites where species declines were particularly pronounced (e.g., upland spruce sites) may have important implications for future tree species compositions [85].

Jack pine and trembling aspen made substantial gains in dominance in the post-fire cohort. These two species are shade-intolerant and require canopy openings from disturbances such as fire to regenerate and are, therefore, successful post-fire species. Spruce species establish shortly after a fire, but appeared to be atypically uncommon in sampled fires from this severe drought-driven fire season, compared to previous studies of the mixedwood boreal zone [86,87] or the northern boreal forest [23,24]. Although spruce trees can persist as suppressed individuals, if seedlings fail to establish following fire they are unlikely to go on to become stand dominants through succession, as the cohort of seedlings established immediately post-fire (1–20 year⁻¹) in boreal forests goes on to make up the future forest [25]. In light of this, jack pine and trembling aspen appear to have gained, at the cost of longer-lived, “late-successional” spruce tree species [88,89]. Additionally, some jack pine stems regularly remained alive post-fire, whereas this was less common for the other three tree species—wildfires killed almost all individuals in burned patches. The successional pathways identified here suggest that increases in burn severity and fire frequency would continue to promote a growing component of jack pine and trembling aspen in northwestern boreal forests, despite topoedaphic, climate, and forest structure controls on post-fire dominance shifts.

In boreal forests, tree species adaptations to wildfire tend to promote “direct regeneration”, where post-fire stands return to pre-disturbance compositions over time. Black spruce is a semi-serotinous species that has demonstrated stand self-replacement following fires in the northern boreal forest [90,91]. The post-fire decreases in the ranked proportional dominance of black spruce that we observed may suggest that increasing fire frequencies and severity may surpass the capacity of this species to re-establish following fires at the proportions previously expected, especially in drier uplands, if there is substantial combustion in both the overstory and understory [24,38,85]. Additionally, in the topoedaphic vegetation classes where the dominance of black spruce was stable (wetlands), and where spruce was previously dominant (spruce uplands), the seedling density was significantly lower than that measured in jack pine and mixedwood uplands. Declines in black spruce dominance relative to early-successional tree species, or through deforestation following severe fires, were observed in Alaska [92,93], the Yukon Territory [37], and in the eastern Canadian boreal forest [24], and this research, provides additional evidence for the potential occurrence of this phenomenon in northwestern Canadian forests.

4.3. Implications for Northwestern Boreal Forests

Wildfire is the stand-initiating disturbance with the largest extent in the northwestern Canadian boreal forest [94]. Therefore, drivers of post-fire vegetation assemblies are an important determinant of future forest composition in this region. The post-fire understory vegetation communities, seedling density, and shifts in the dominance of tree species were primarily attributable to pre-fire forests, climate, and topoedaphic context, suggesting that there is substantial capacity for forests and understory vegetation communities to regenerate post-fire. Although some variability in post-fire communities was attributed to burn severity, burn severity in this region is also associated with pre-fire forest structure and composition [28], further reinforcing the importance of pre-fire drivers to observed vegetation assemblies.

Despite this resilience, long-term shifts in tree species compositions are ongoing in parts of the western Canadian boreal forest, with proportions of early-seral shade-intolerant species such as jack pine and trembling aspen demonstrating increasing prevalence, driven by climate change [88]. Simultaneously, droughts appear to have caused decreases in forest productivity, altered seedling establishment and caused large-scale die-offs of mature trees in northern forests [95,96,97]. Although strong non-fire controls on understory vegetation and seedling establishment offer some resilience to change, the secondary direct and indirect effects of fire will likely serve to accelerate these ongoing changes if fire size, frequency, and severity increase as projected [40,41,98].

An increasing broadleaf component in northwestern boreal forests, such as that observed in this study, may reduce fire severity and flammability of boreal forests [28,99,100] and raise the surface albedo [101], potentially offering a negative feedback to shifts driven by climate change and impeding increases in fire activity [16,102]. This effect would be transient if increases of the proportion of trembling aspen in boreal forests do not persist. Increases in the frequency and severity of droughts may lead to a subsequent decline in this drought-sensitive species [97,103]. Drought stress would also likely further exacerbate black spruce declines and potentially favour more drought-tolerant upland conifers such as jack pine [104,105,106]. Furthermore, droughts increase the susceptibility of fuel-limited young forests to reburning [33,107], which could yet again reinforce reductions in black spruce dominance through the reduced availability of viable seeds. Observed post-fire seedling density was highly variable, substantially different from pre-fire species compositions, and several sites experienced near regeneration-failures when severely burned at short fire frequencies. This research contributes to the growing body of literature indicating that changes to forests of this region are ongoing, despite the overarching resistance to such shifts conferred by regeneration mechanisms and topoedaphic controls [83,92,93].

4.4. Limitations and Future Research

Due to the opportunistic nature of this study’s sampling design, we were unable to measure changes in the composition of the understory vegetation assemblies from pre-fire to post-fire communities. Although we partitioned the variance in post-fire vegetation communities to identify some role of wildfire on their determination, studies where prescribed burns are planned or existing plots are burned over in natural fires are better positioned to measure shifts in the dominance of understory vegetation species from pre-fire to post-fire conditions. Soils are important to post-fire understory vegetation communities and they are also relevant to seedling recruitment through the provision of seedbeds [108,109,110]. Such studies would allow researchers to measure the changes in soils as a result of fires, including changes in the organic layer depth.

We conducted our field sampling one year and three years post-fire. Studies spanning a longer time period can provide additional insights into post-fire vegetation recovery for both understory plants and trees e.g., [18,26], but this was beyond the scope of this work. An assessment of whether a forest has recovered to a state similar to pre-fire conditions would require an extensive period of time, reflecting the growing conditions at high latitudes and local disturbance regimes (e.g., stand ages at the time of burning ranged from 9 to 237 years in this study). Although the post-fire recruitment pulse for some tree species may not be complete three years post-fire, there is ample evidence that the recruitment occurring within the first few years post-fire largely determines the future species composition and structure of the stand in boreal forests [25,108]. Just over one-third of our plots did not have repeated measures of seedling recruitment data sampled three years post-fire; however, the majority of these sites that were not revisited were non-forested wetlands (fens), with no trees pre-fire. This data gap affected closer to a quarter of the forested sites. Spruce trees are slower to establish, and this may have biased our results; however, we did observe some spruce seedling recruitment in most plots with a pre-fire spruce presence. By converting our measurements of post-fire seedling density to compositional log-ratio data, and calculating shifts in the ranked dominance of species, we captured changes in post-fire tree species composition and normalized the highly skewed seedling density data. This method does not permit us to assess structural changes that may have occurred; for example, whether post-fire forest density increased or decreased and whether these outcomes vary by species. Although we did characterize some variability in seedling density by topoedaphic vegetation communities and fire frequencies, future research could combine ranked dominance shift data with seedling and pre-fire stem densities to directly capture regeneration failures and structural changes, in addition to the shifts in proportional dominance measured here.

At the time of burning there was an ongoing multi-year drought in the study area, which continued into 2015, and may have affected the post-fire recruitment of seedlings, as well as their growth [95,111]. The drought conditions may also have affected burn severity of the fires, as fire weather is significantly related to overstory and understory combustion in boreal forests [28,112]. The identified impacts of fire on understory vegetation and seedling recruitment may have been influenced by these pre- and post-fire environmental conditions and therefore, the observed vegetation assemblies and ecological outcomes may be most representative of severe fire years. Sampling in wildfires that occurred in different years would capture a wider range of pre- and post-fire climates.

5. Conclusions

In this study, the primary determinants of post-fire outcomes for boreal forest vegetation communities and shifts in tree species dominance were pre-fire forests, topoedaphic context and climate. Burn severity, fire history, and post-fire soils were significant secondary drivers. Burn severity and fire history did not significantly explain the variability in understory vegetation communities; however, post-fire soils were related to understory vegetation community dissimilarities. Furthermore, burn severity was significantly related to understory vascular plant diversity. Severely burned vegetation communities tended to have lower understory species richness and diversity, as did very wet sites, which typically burned at low levels of severity. Post-fire shifts in tree species dominance, as characterized by differences in ranked proportional compositions, were significantly related to fire history and burn severity, but this effect was less important than pre-fire and climatological conditions. The overriding control of fire-independent drivers on post-fire vegetation may provide some resilience to forests in the face of climate change, as they are less susceptible to fire-mediated type conversions due to site moisture and pre-fire forest drivers. Despite this potential for resilience, changes to forest vegetation community compositions due to altered climates are occurring, and burn severity and fire history were important explanatory variables in our models of shifts in tree species dominance. In a forest with potentially increasing frequency, size, and severity of fires, the long-term resistance to change conferred by topoedaphic and forest controls may be overwhelmed by the direct and indirect effects of wildfires, which offer pathways to change. Burned sites will also experience altered post-fire climates, with potential increases in moisture stress and droughts which would exert additional pressures on initial post-fire vegetation. Ongoing shifts in the dominance of tree species are the result of both climate and fire. If these disturbances continue to increase, the observed shifts towards early-seral species such as jack pine and trembling aspen could produce large-scale changes in vegetation dominance that may lead to substantial—and perhaps unanticipated—ecological changes.