Wildfires are increasing in importance in many regions of the Canadian boreal forest and are an ongoing risk for forest management activities. We simulated the effects of fires on long-term harvest levels on the 59 forest management units of the province of Québec, Canada, for the 2020–2100 period. Different climate change pathways (stable, RCP 4.5 or 8.5) and salvage logging rates (20% or 70% of mature burned stands) were simulated. Changes in forest flammability due to climate change, species migration, and forest management were also considered. Under stable climatic conditions, the decline in potential harvest levels due to fire, based on 50 simulations per scenario, ranged between 3% and 33% (mean = 11%) when high salvage logging rates were simulated, compared to 6%–45% (mean = 20%) for low salvage rates. Climate change caused increases in burn rates between −3% and 39% for RCP 4.5 and between 33% and 69% for RCP 8.5 at the end of the 21st century, depending on fire zones. However, the effects of these modified burn rates on harvest levels did not differ substantially from those of baseline burn rates, probably because the projected burn rates were highest during the later part of the simulations (2070–2100), when their impacts on harvest level calculations were limited. This study indicates that potential harvest levels calculated without considering wildfires are likely to be non-sustainable.

1. Introduction

Wildfires have an important influence on forest management activities. In the short term, fire suppression is expensive (Stocks and Martell 2016), as are post-fire silvicultural interventions to restore or maintain forest productivity (Cyr et al. 2022a). In the long term, wildfires also cause decreases in sustainable harvest rates (Van Wagner 1983; Boychuk and Martell 1996) and may exert great pressure on other ecosystem services such as carbon sequestration (Colombo et al. 2012). In Canadian boreal forests, negative impacts of fire on forest management and various ecosystem services are expected to increase in the future, because climate change will likely lead to an increase in the frequency and duration of fire-prone weather events (Gauthier et al. 2015a; Hanes et al. 2019).
Maximizing salvage logging of recently burned forests helps decrease wildfire-induced impacts on long-term timber supply (Leduc et al. 2015). However, field implementation of post-fire salvage logging is often impeded by different factors, such as limited road access to remote areas, the presence of burned stands that are too young and have little merchantable value for the transformation industry (Gauthier et al. 2015b), or by the rapid degradation of dead trees by insects or rot (Barrette et al. 2015; Boucher et al. 2020). In some jurisdictions, regulations also require that a proportion of burned stands is left unharvested to create habitat for fire-dependant animal or plant species (Nappi et al. 2004; Saint-Germain and Greene 2009). Overall, due to these constraints, salvage logging is not sufficient to completely offset the negative effects of fires on sustainable harvest rates (Leduc et al. 2015).
Another commonly used approach to minimize fire risks on long-term harvest levels is to periodically re-calculate sustainable harvest levels, e.g., every 5 or 10 years. This “replanning” approach allows managers to adjust management objectives and the corresponding harvest levels to take into account any deviation from the initial projections due to economic, political, social, or ecological factors (Paradis et al. 2013). Specifically, in response to wildfires, replanning allows to gradually adjust harvest levels following individual fire events to reduce the risk of abrupt, unforeseen declines in wood supply that would ultimately occur without periodic replanning (Armstrong 2004). By simulating this replanning process through wood supply calculation models coupled with natural disturbance models, it becomes possible to anticipate declines in harvest levels over longer periods (e.g., 50–100 years), according to different disturbance scenarios (Armstrong 2004; Savage et al. 2010). Knowledge of these potential declines could help managers minimize the associated risks by, for example, pre-emptively reducing harvest levels (Raulier et al. 2014; Leduc et al. 2015) or by planning across broader spatial extents (Armstrong 2004; Cumming and Armstrong 2004).
Weather patterns are widely expected to become more conducive to fire in the future across Canada (Wang et al. 2020), as is already the case in parts of the country (Hanes et al. 2019). However, the extent at which fire-conducive weather will translate into actual increases in the size, frequency, or severity of wildfires is subject to considerable uncertainty (Boulanger et al. 2018). This uncertainty is amplified by other factors such as the effectiveness of fire suppression (Cardil et al. 2019) and complex interactions with changes in the flammability, abundance, and spatial arrangement of vegetation or fuel types that may influence future fire regimes (Girardin et al. 2013; Bernier et al. 2016; Marchal et al. 2020).
The main objective in this study was to develop a modelling approach to evaluate the effects of future wildfires on long-term harvest levels across a large territory (the province of Québec). Specifically, we used a spatially explicit landscape simulation model developed by Bouchard et al. (2019) to integrate the potential influence of different factors on harvest levels: climate-induced changes in fire frequencies, changes in tree species composition (disturbance- and climate-induced), and variation in rates of post-fire salvage logging. This was undertaken by formulating distinct management and climate change scenarios and by determining their effects on potential harvest levels.

2. Methods

2.1. Study area

The study area (740 000 km2, including water bodies) is in the southern portion of the province of Québec, Canada (Fig. 1a). Forests located north of the study area are not dense or productive enough to sustain commercial forest management activities (Jobidon et al. 2015) and were excluded from our study. The southern fringe of the province, which includes the main towns, cities, agricultural areas, as well as smaller private woodlots, was also excluded because this area is highly fragmented by human use and less affected by wildfire (Chabot et al. 2009; Couillard et al. 2022). Thus, our study area includes 59 forest management units (FMUs) located on public lands as of 2020 (Fig. 1d). According to provincial regulations and forest management guidelines, a sustainable, non-declining harvest level is expected from each of these units over the long term, usually defined as a planning horizon of 150 years (Rheault 2013). The study area includes three ecological forest zones, arranged along the north–south climatic gradient: the boreal, mixed, and temperate zones (Fig. 1b). Historical fire regimes in this territory were recently studied by Couillard et al. (2022), based on fire occurrence, ignition type (human or lightning), and fire size data collected for the 1890–2019 period. These fire regimes correspond to 13 fire zones (Fig. 1c).
Fig. 1.
Fig. 1. Maps of the study area, showing (a) location within Canada, (b) bioclimatic zones, (c) fire regime zonation with historical burn rates (% per year), and (d) 59 forest management units (as of 2020). In panel c, historical burn rates were not measured in territories dominated by land uses other than forests such as agriculture or urban areas; wildfires are possible in those territories but were not simulated in this study. In panel d, key characteristics of the numbered forest management units are indicated in Supplementary Table S2. Background map credits: Esri (water bodies) and Québec Government (Bioclimatic zones, fire zones, and forest management units: Données Québec). Projection: NAD83 Quebec Lambert.
All basic forest information, including current forest characteristics such as stand composition and stand age, was obtained from provincial forest geospatial information databases (https://www.foretouverte.gouv.qc.ca/, last accessed 12 December 2020). Forest composition was categorized in nine vegetation types according to the dominant tree species or species groups when no species was sufficiently dominant (Table 1).
Table 1.
Table 1. Relative proportion (%) of the nine different vegetation types among mature forests under initial conditions (2020), per bioclimatic zone.

2.2. QLDM modelling framework

Quebec-LDM (or QLDM, for Quebec Landscape Dynamic Model) is a spatially explicit model that simulates the effect of disturbance processes (anthropogenic and natural) on changes in forest characteristics across large regions (e.g., >10 000 km2) (Bouchard et al. 2019). It allows the simulation of stand-replacing wildfires, clearcuts, partial cuts, and salvage logging (i.e., clearcuts performed immediately after fire). In QLDM, sustainable harvest levels are determined for each FMU using a wood supply model based on the age structure of the forest and stand productivity (i.e., age of maturity of different stand types). Species composition changes are simulated by using disturbance-specific transition matrices. In a climate change context, the transition probabilities are modified according to species migration capabilities and species tolerances to different climatic conditions, as determined by species distribution models (Périé et al. 2014). QLDM also allows to vary the degree of species persistence to reflect uncertainty in the response of different species or populations to future conditions (Bouchard et al. 2019). The main processes and variables used in this model are summarized in Fig. 2. The model uses 2 km × 2 km grid cells and a simulation time step of 5 years. It is programmed in the R stats package (R version 4.0.3; R Core Team 2021), and all scripts are available at https://github.com/QLDM-Organization.
Fig. 2.
Fig. 2. Main processes simulated in the Quebec-LDM modelling framework, with particular emphasis on fire simulation (see detailed descriptions in Methods).

2.3. Fire simulation

2.3.1. Target burned area

Before simulating the ignition and spread of fires, we defined a mean target burned area per 5-year period for each fire zone. This burned area was calculated through the following steps. First, we used historical (1890–2019) information (Couillard et al. 2022) to define a baseline mean periodic burn rate (BRz) calculated for each fire regime zone z (Fig. 1c), corresponding to the proportion of the landscape that burned per 5-year period:
where zone area includes all forest stands that could potentially burn (thus excluding water bodies, agricultural or built areas, etc.). In a second step, projections of changes in the duration of fire-prone weather events were summarized with the spread event day probability index (SEP; eq. 2; Wotton et al. 2017), a statistic derived from the daily Canadian Fire Weather Index (FWI) System (Van Wagner 1987):
Daily values of SEP are summed over a temporal resolution of one fire-year. Yearly changes in the SEP according to climatic projections were obtained for each fire regime zone between 2010 and 2100 and averaged over 30-year periods 2010–2040, 2041–2070, and 2071–2100 (see Section 2.5). These projections were used to calculate burned areas for each 5-year period p until 2100, assuming that the relationship between SEP and the area burned was linear, according to the following equation:
where TBAz,p represents the target burned area in zone z and period p, SEPz,p is the average SEP value in each zone z and period p, and SEPz,p=1 is the average SEP in each zone z during the first period. PAz is the total area that could potentially burn in zone z, a value that doesn't change with time. Even though the relationship between SEP and fire size or area burned corresponds in theory to a power function (Wang et al. 2020), this relationship is noisy and difficult to parameterize when fire frequencies are low and could produce very high values outside the range of observations. In this study, we assumed that this relationship was linear for simplicity, and because a linear relationship is probably less likely to deliver biased results when extrapolating.
Changes in fuel characteristics and in the spatial distribution of different vegetation types are expected to influence fire spread and, ultimately, fire size. In eastern Canada, for instance, a landscape-level increase in young (<40 years old) stands or hardwood-dominated stands due to harvesting is expected to generate conditions that are less fire-prone (Bernier et al. 2016; Marchal et al. 2020). In this study, the different types of dominant vegetation were associated to different spread potentials, depending on their general flammability (Quinby 1987; Pelletier et al. 2009; Terrier et al. 2013; Bernier et al. 2016). We assumed that in each fire zone, a higher proportion of less-flammable vegetation types translated into a lower overall burn probability (Marchal et al. 2020). Hence, a mean vegetation flammability (VFz,p) value was calculated in each zone at every time step according to the relative abundance of the four fuel types (eq. 4):
where is the proportion of cells in zone z that are of fuel type i during period p, and Wi is the spread potential of the fuel type (Table 2). Mean vegetation flammability was then used to modify the target burned areas (TBAz,p; eq. 3), according to eq. 5.
where VFz,p=1 is mean vegetation flammability at the beginning of the initial simulation period, TBAz,p is defined in eq. 3, and TBAfz,p is the modified burned area according to changes in mean vegetation flammability.
Table 2.
Table 2. Vegetation types and corresponding spread potentials (W) (from 0 to 1) used in the fire spread algorithm.

2.3.2. Fire ignition and fire spread

In each zone and during each 5-year period, new fires were simulated until the target burned area (TBAfz,p) was reached. Each fire was simulated according to the following steps. First, a target fire size was drawn from a historical fire size distribution derived from the 1890–2019 baseline data set (Supplementary Table S1). For each individual fire, an ignition point was selected randomly across each fire zone. The fire was allowed to spread from this ignition point until one of two outcomes occurred: (i) the fire reached the pre-defined size or (ii) the fire extinguished as its front reached cells with no fuels or low-flammability fuels. The last fire simulated in a given fire zone for a given period was allowed to burn to its full size even if the target burned area was surpassed, and the excess area was subtracted from the target burned area for the next period. Fires were allowed to cross FMU limits and fire zone limits.
In the fire-spread algorithm, the fire was more likely to spread to neighbouring cells occupied by highly flammable vegetation than to those occupied by less flammable vegetation (Box 1 of Fig. 2). Spread probability was also set slightly (30%) higher to cells downwind from the grid cells that were ignited, based on dominant spring and summer wind directions calculated with Biosim (Régnière et al. 2013) for the 1981–2010 period. We did not simulate long-range fire spotting or the effects of variation in burn intensity on fire spread.

2.4. Forest harvesting simulations

Forest harvesting was simulated in two steps. The first step consisted of calculating a sustainable harvest level (HLu,p) for each forest management unit u and each period p (Box 2 of Fig. 2). This harvest rate was calculated in area (km2 harvested/period). The second step consisted of selecting the cells to be harvested or salvaged at each 5-year period (Box 3 of Fig. 2), according to the harvest level calculated during the first step. These two processes are described below.

2.4.1. Calculation of potential harvest levels during each period

Potential harvest levels were determined at the FMU level. They were calculated by considering the age structure of the forest, the productivity of the different stand types (age of maturity), as well as the type of silvicultural regime through which each stand was treated (clearcutting or partial cutting). Potential harvest levels were recalculated at each 5-year period, to take into consideration changes in forest structure caused by fire events and harvestings. Periodic recalculation of harvest levels usually results in a decline in harvest levels through time because a portion of mature or premature stands are burned at each time step (Armstrong 2004).
Because harvest levels were calculated using an area-based approach, they do not consider the variation in volume associated with species and ecological context. Biomass production is also expected to be influenced by climate change (D'Orangeville et al. 2018), something that is not integrated in the harvest calculations. From an integrated management perspective, several concerns other than forest productivity must also be considered to ensure sustainable harvest levels. Management objectives such as keeping a minimal abundance of old seral stages for biodiversity conservation objectives, social values, operational considerations (e.g., road access or transportation distances to mills), market demand for different wood products, and the potential effects of other disturbance types (e.g., insect outbreaks) are all expected to influence the optimal harvest level. In this study, we did not consider most of these factors and focused on the potential effect of fires. Only stands located in protected areas or those that were not accessible according to government databases were considered non-harvestable during the 2020–2100 interval.
For each FMU, harvest levels were defined separately for stands treated by clearcuts or by partial cuts. Broadly consistent with current practices in Québec, 90% vs. 10% of stands dominated by temperate hardwoods were treated using partial cuts vs. clearcuts, respectively. The proportions were reversed in conifer- and boreal hardwood-dominated stands (90% were treated with clearcuts and 10% with partial cuts). For simplicity, in a stand managed with clearcuts, the age of stand maturity was assumed to be 100 years in the boreal zone, 90 years in the mixed zone, and 80 years in the temperate zone (Fig. 1), irrespective of local forest composition. In the case of stands managed by partial cuts, we considered that trees were harvested in two passes, with a cutting cycle that was half the age of stand maturity as defined above. In the context of our area-based harvest calculations, an area harvested by partial cuts was considered equivalent to half the area harvested by clearcuts.

2.4.2. Selection of harvested cells during each period

The number of cells to be harvested by partial cuts or clearcuts in each FMU was the same as the potential harvest level defined at each period. The cells were randomly selected in the pool of available mature cells in each unit at each period. When fires burned some cells in a FMU during a given 5-year period, the harvest of mature burned cells was prioritized (salvage logging), but the overall harvest level during this period (in burned + unburned cells) could not exceed the potential harvest level calculated beforehand. Moreover, burned trees deteriorate and cannot be used by the transformation industry after a few years. Hence, at the end of each simulation period, the age of burned cells was reset to 0, whether the cell had been salvaged or not, to emulate stand recovery following natural regeneration or plantation. For example, if a 12 km2 fire burned in a FMU where the potential harvest was 8 km2/period, only this 8 km2 was salvaged, while the residual 4 km2 was not harvested and left to regenerate naturally. For simplicity, we assumed there were no regeneration delays or complete regeneration failure in forest stands impacted by fire or harvesting.

2.5. Climatic information

Baseline information about recent (1980–2010) climatic conditions in the study area was taken from Worldclim (Fick and Hijmans 2017). Information about future climatic conditions was produced by the Earth System Model MIROC (Watanabe et al. 2011). Projections of the MIROC model were used in this study because they are intermediate in terms of projected temperature and precipitation compared to those of other Earth System Models (Boulanger and Pascual 2021).We used two representative concentration pathways (RCPs), 4.5 and 8.5, whose radiative forcing rises gradually to reach 4.5 or 8.5 W·m−2 in the year 2100, respectively.
The climatic scenarios were used to calculate two groups of climate variables that drive different ecological processes in the overall model (Fig. 2). First, changes in the mean annual temperature and total annual precipitation were used to track the potential northward displacement of the climatic niche of each tree species (Bouchard et al. 2019). Climatic suitability for the different tree species was inferred from the current climatic envelopes of each tree species (Périé et al. 2014), and the transition probabilities from one dominant species to another were determined by the climatic suitability of both species, as well as soil conditions and the availability of seed sources (see details in Bouchard et al. 2019). Outputs from the climatic model were averaged over a period of 5 years for this first group of variables (i.e., mean annual temperature and total annual precipitation) and updated in the model at each 5-year interval. The second group consisted of climatic variables that were used to predict fire weather (see details in Section 2.3). This second group of variables was averaged over a period of 30 years and updated during simulations at years 2040 and 2070.

2.6. Simulations

The influence of wildfires, salvage logging, and forest dynamics on future harvest levels was examined through the simulation of different scenarios (Table 3). Scenarios 1 and 2 are baseline scenarios isolating the effects of fire and harvesting, respectively. Scenarios 3 and 4 include both harvest and fire, but low rates of salvage logging in mature burned stands following fire were used in scenario 3 (20%), compared to relatively high rates in scenario 4 (70%). These different salvage logging scenarios reflect the uncertainty surrounding the applicability of salvage logging, particularly in the boreal zone (Nappi et al. 2004). For the scenarios that included fire, 50 simulations were realized to capture the variability caused by the random occurrence of fires of various sizes within the different FMUs. Scenario 2 was run only once, given that it did not include the stochastic occurrence of fire. For each scenario, 50 simulations were run for each RCP (stable, RCP 4.5 and RCP 8.5).
Table 3.
Table 3. Parameter settings for the different simulation scenarios.

2.7. Statistical analyses

A statistical analysis was conducted to identify the factors that influence the decline in harvest levels through time, using a linear model with the “lm” base function in R. Since the interpretation of p values can be misleading in the analysis of simulation model results (White et al. 2014), we used tvalues to assess the relative importance of different predictors on harvest levels. The response variable in this analysis, the decline in harvest level (DHL), was calculated according to eqs. 6 and 7:
where min(HLu,sc,s) is the minimum harvest level observed in FMU u, scenario sc, and simulation s, and MinHLu,sc,s is the same quantity expressed as a proportion of the initial harvest level (during period 1). DHLu,sc is the median decline in harvest level in management unit u and scenario sc. The predictors used in the statistical analysis were baseline burn rate at the FMU level, climate change scenario, salvage logging intensity, and log-transformed FMU's size (to obtain a normal distribution).
The potential harvest level varied from one simulation to the next according to variation introduced by the stochastic occurrence of fires. In a long-term planning framework, we wanted to evaluate the probability that a given harvest level could be maintained over the long term. To do this, we first identified the minimal harvest level observed during each simulation s during the 80-year simulation period. Using results from the 50 simulations, we calculated the lowest, median, and highest minimal harvest level for each scenario and each FMU. These threshold values allowed us to define potential harvest levels that would have a low, intermediate, high, and very high probability of being non-sustainable over the 2020–2100 period.

3. Results

3.1. Burned area

Burned areas are reported at the provincial scale by being summed across all FMUs. At this scale, burned areas were higher in the RCP 4.5 and RCP 8.5 climatic scenarios compared to stable climatic conditions, regardless of the management scenarios (Fig. 3). At the end of the simulations, mean area burned increased by between 10% and 11% under the RCP 4.5 pathway compared to the reference scenario (stable climate), whereas mean area burned under the RCP 8.5 pathway increased by 30%–38% (Fig. 3). These results also show that the overall burned area was roughly similar for the different management scenarios but tended to be slightly higher in the “no harvest” baseline scenario. Overall, at the end of the simulations (i.e., year 2100), burn rates in the “no harvest” scenario were 6%–10% higher than in the other scenarios (Fig. 3, upper row).
Fig. 3.
Fig. 3. Variation in total burned area (upper row) and mean vegetation flammability (VF) (bottom row) according to different levels of salvage logging. The curves represent values that were averaged from 50 simulations.
The "no harvest" scenario also resulted in a mean vegetation flammability that was considerably higher than the other scenarios (Fig. 3, lower row). By comparison, the influence of climate change on mean vegetation flammability, when isolated, was relatively low at the provincial scale (Fig. 3, lower row).
We also examined the effects of changes in vegetation flammability and climate on burn rates by fire zone, under the RCP 8.5 pathway. The effects of changes in vegetation flammability on burn rates varied regionally: a marked decrease in mean vegetation flammability was observed in fire zone A, located in the mixed and temperate forest zones, whereas the opposite was observed in fire zone G, which is located in the northern part of the study area, with relatively weak effects elsewhere (Fig. 4). Overall, when climate change and changes in mean forest flammability are considered, mean burned area increased by 33%–69% at the end of the planning horizon in the different fire zones for the RCP 8.5 pathway (Fig. 4) and between −3% and 39% for RCP 4.5 pathway (results not shown).
Fig. 4.
Fig. 4. Effect of climate and vegetation flammability modifiers on the overall burn rates in the different fire zones (A–M). The burn rates are the average of the 50 simulations realized for scenario 4 (Table 3), using the RCP 8.5 climate scenario. Note that the scale of the y-axis is identical in all zones except for zone G.
Fires were allowed to spread preferentially towards more flammable fuels. The effect of this process on model outcomes was assessed by comparing the relative abundance of the three fuel classes (high, intermediate, and low flammability) in burned areas, against their relative abundance over the whole fire zone (Supplementary Fig. S1). As expected, the results show a consistent preference for highly flammable fuel types, whereas the representation of intermediate and low-flammability fuels within fire perimeters was lower than their abundance at the zone level.

3.2. Harvest levels: trends at the province level

Overall, at the scale of the province and at the end of the planning horizon, there was a 10%–18% decrease in harvest levels in 2100 in the scenarios with fire, compared with the no-fire scenarios (Fig. 5). At this scale, the decline in harvest levels was gradual until around 2060, but harvest levels increased afterwards. Despite this increasing trend in the later part of the century, the difference between harvest levels in the no fire baseline scenario and the scenarios with fire remains stable or increases during the same period (Fig. 5, upper row).
Fig. 5.
Fig. 5. Total harvested area for the different salvage logging and climate change scenarios. The upper row shows the area harvested at each period, and the bottom row shows the proportion of the overall harvested area that was cut in recently burned areas (salvage logging). Each curve represents the average of 50 simulations, except for the fictive no-fire baseline scenario for which a single simulation was run.
On average, burned stands that were salvage-logged represented between 2% and 4% of the area harvested at the province level when 20% of mature burned stands were harvested, against 8%–14% when a 70% salvage logging rate was simulated (Fig. 5).

3.3. Harvest levels at the forest management unit scale

In most FMUs, fire caused a decline in harvest level through time. The steepness of this decline and the range of maximal declines observed in the 50 simulations varied among FMUs. Statistical analyses showed that these declines were more pronounced in smaller FMUs, in FMUs for which initial burn rates were relatively high, and in scenarios for which the rates of salvage logging were relatively low (Table 4). Climate change resulted in higher burn rates (Fig. 3, upper row), but these higher burn rates did not translate in markedly lower harvest levels, and the weight of this variable in the statistical results was relatively low (Table 4). On average, under stable climatic conditions, the median decline in harvest levels at the end of the planning horizon (compared to first-period harvest levels) were between 5.7% and 44.7% (mean = 19.8%) when low salvage logging rates were used, and between 3.2% and 33.3% (mean = 11.1%) when high salvage logging rates were used.
Table 4.
Table 4. Results of the statistical model predicting the effects of FMU characteristics (size and mean burn rate), climate change scenario, and the rate of salvage logging on the median decline in harvest levels at the FMU level.
Potential declines in harvest levels relative to initial conditions under the RCP 4.5 scenario are presented in Supplementary Table S2 for all FMUs. Results of the 50 simulations are shown for two FMUs as examples (Fig. 6). These examples show how the ranges of potential declines in harvest levels can vary between FMUs, as well as variation from one stochastic realization to the next.
Fig. 6.
Fig. 6. Variations in potential harvest levels caused by fires in two forest management units. Each black curve represents the result of one simulation. The red curve indicates the baseline harvest level in a fictive no-fire scenario. The horizontal dashed lines indicate, for a given FMU, the minimal, median, and maximum harvest levels across all simulations during the period when the harvest level reached its lowest value (see eqs. 6 and 7). These limits were used to identify harvest levels that would have a very high (yellow), high (light orange), low (dark orange), or very low (red) probability of being sustainable over the 2020–2100 period.

4. Discussion

4.1. Effects of wildfires on harvest rates

In our study, harvest levels were calculated based on forest age structure and age of stand maturity, for each forest management unit (FMU), at each time step. The algorithms used to make these calculations were relatively simple, which allowed us to simulate many replicates of each climate or management scenario, with random fire occurrence, over a large and heterogeneous territory. The results show gradual declines or stagnation in harvest levels caused by wildfires during the 80-year simulations. The steepness of these declines were markedly different among FMUs: in FMUs with relatively high burn rates, the gradual declines in harvest levels that were observed, i.e., between 30% and 50% of the initial harvest levels (Fig. 7; Supplementary Table S2), were similar to those observed by Savage et al. (2010) in Ontario and by Armstrong (2004) in Alberta. These declines were much less important in the FMUs located in the eastern and southern parts of the study area, where burn rates were lower.
A 70% salvage logging rate lessened the negative impacts of wildfires on harvest levels, compared to a 20% salvage logging rate. Salvage logging can reduce fire impacts on future harvest levels by prioritizing the harvest of mature stands that burned recently (<5 years) while re-scheduling the harvest of unburned mature stands to later planning periods (Leduc et al. 2015). We considered the harvesting of 70% of recently burned mature stands to be a relatively optimistic scenario: salvage logging rates higher than 70% would probably be unrealistic considering the road access issues that often impede salvage logging operations, particularly in remote boreal forests (Nappi et al. 2004). As shown in Fig. 7, forcing a 70% salvage logging rate allows to substantially attenuate the impacts of wildfires on harvest levels compared with a 20% rate, except in FMUs located in the northwestern part of the study area, where the burn rates were particularly high. FMU size was also an important predictor of wildfire impacts on harvest levels, with steeper declines observed in smaller FMUs. This happened because some of the smaller units (e.g., <1000 km2; Supplementary Table S2) were proportionately more affected by large fires (e.g., >500 km2). When a large proportion of the FMU is affected by fire during a given period, this results in much sharper declines in harvest levels in subsequent periods. The negative influence of such large fires was dampened in larger management units (e.g., >10 000 km2).
Fig. 7.
Fig. 7. Median declines in harvest levels observed in the different FMUs of the province, as a proportion of the initial harvest level (during the first 5-year period, 2020–2025). Panel (a) corresponds to a 20% salvage logging scenario and panel (b) to a 70% salvage logging scenario. Background map credits: Esri (water bodies) and Québec Government (forest management units: Données Québec). Projection: NAD83 Quebec Lambert.
At the scale of our study area, burned area increased by up to 38% at the end of the simulation under the 8.5 RCP scenario compared to the stable climate scenario (Fig. 3, upper row), increases that were even more important in some fire zones (Fig. 4). However, these temporal changes in burn rates had a relatively weak effect on harvest levels, especially compared to the spatial variations in baseline burn rates (Table 4). This limited effect could be explained partly by the fact that the most important increases in burned area took place at the end of the simulation horizon (after 2070). Conducting simulations over a longer time frame might have increased the weight of climate change on harvest rates in the statistical analyses. Similarly, taking into consideration post-regeneration failure after fire, particularly in northern black spruce forests (Baltzer et al. 2021; Cyr et al. 2022a), could have increased fire impacts on sustainable harvest levels: in this study, we assumed that forest stands are immediately reforested after a fire, naturally or through plantation.

4.2. Changes in burn rates through time

Simulation outcomes showed an increase in area burned during the 2020–2100 period under climate change scenarios. Changes in mean vegetation flammability were also observed during the same period. Even though previous studies already showed that an increase in (less flammable) hardwoods could affect fire regimes (Girardin et al. 2013; Terrier et al. 2013; Marchal et al. 2020), our study indicates that this effect is most likely to take place at the interface between the boreal and temperate forest zones, where the increase in hardwood abundance is likely to occur much faster than in northern boreal forests (Bouchard et al. 2019; Boulanger et al. 2022).
The results concerning the relative importance of fuels and climate on burn rates should be interpreted with caution since the parameterization of fire models is subject to some uncertainty. For example, fire ignition and fire spread are influenced by daily and hourly weather patterns that are impossible to integrate in large-scale forest dynamics simulations, which typically operate on weather data aggregated across multiple years. It was also difficult to use the results from previous studies conducted in regions outside our study area for calibration of the fire simulation module, since fire regimes are a product of nonlinear interactions between landform, forest structure, weather, and anthropogenic influences and are difficult to generalize or extrapolate from one region to another (Allaire et al. 2020).

4.3. Management implications

Ensuring some form of non-declining flow of timber to the wood transformation industry is an essential element of sustainable forest management. However, truly stable harvest levels are difficult to calculate precisely, particularly in a context where unforeseen disturbances such as wildfires may occur (James et al. 2007). Our study indicates that potential harvest levels calculated without considering fire are non-sustainable in many FMUs of the province, for all climate change scenarios (including the baseline). For FMUs in which burn rates are significant, it has been recommended to lower harvest levels as a precautionary approach to increase the stability of the long-term wood supply (Leduc et al. 2015). Even if an a priori decrease in harvest rates is an obvious means to increase long-term stability, determining the appropriate harvest rate can be difficult, as this is not a strictly scientific decision, rather involving a multitude of stakeholders who might have different tolerance levels to different socio-economic or environmental considerations. Moreover, these decisions must be taken at the FMU level because factors such as forest age structure, demand from the transformation industry for different tree species, or ecosystem protection objectives can vary among FMUs. Such analyses should be clearly presented in the management plans, something that is presently poorly developed in Québec and elsewhere in Canada (Gauthier et al. 2015a). We suggest that tolerance to the risk of potential wood supply disruptions or environmental degradation due to the cumulative effects of wildfires and harvestings should be discussed explicitly before deciding of a proper long-term harvest level.
Even if the capacity of fire management agencies to suppress fires is expected to decrease in a climate change context (Podur and Wotton 2011), our study shows that forest managers can mitigate the risks posed by these fires on harvest levels up to a point by maximizing salvage logging and by changing the size or configuration of the management units in which harvest rates are calculated. While the role of salvage logging to reduce fire impacts on long-term harvest levels has been highlighted by previous studies (Leduc et al. 2015), the importance of management unit size has been less emphasized so far (but see Cumming and Armstrong 2004). Ultimately, enlarging or merging some management units might reduce the volatility introduced by large fires on wood supply. However, salvage logging and FMU mergers involve trade-offs that were not examined in detail in this study. Specifically, an increased use of salvage logging will cause an increase in operational costs, and burned trees are less interesting for the transformation industry (Mansuy et al. 2018). As for FMU mergers, they might not be feasible everywhere, as original FMU delineations were not necessarily dictated by operational considerations. For instance, the relatively small FMUs located in fire zone G, to the northwest of the study area, were negotiated as part of the “Paix des Braves” agreement between the Cree Nation and the Québec Government (Cyr et al. 2022b). Changing those delineations could involve complex political negotiations. Hence, even if salvage logging and FMU mergers appear to be solutions that could help mitigate wildfire risks on harvest levels, they might be difficult to implement in practice.
Several factors other than fire that could influence sustainable harvest levels were purposely left aside in this study. For instance, ecosystem management targets, social acceptability of forest management, operational constraints, or intensive silvicultural measures that could increase forest productivity were not incorporated in our modelling framework. Because burn rates are very variable among fire zones, we suggest that if silvicultural investments are to be made to compensate the expected losses due to fire, these investments should be concentrated primarily in fire zones where burn rates are expected to remain relatively low in the future (Reed and Errico 1986; Cyr et al. 2022a). In our study, we also considered that the entire potential harvest level for each FMU was effectively harvested at each time step, something that is rarely achieved in Québec (Paradis et al. 2013). Using an effective harvest level that is somewhat below the potential harvest level would probably help stabilize harvest rates in the long term. Simulation approaches such as the one developed in this paper could be used in the future to further document the trade-offs among different approaches to mitigate fire impacts under different climatic pathways.


We thank the two anonymous reviewers for their insightful comments on the manuscript.


Allaire F., Filippi J.-B., Mallet V. 2020. Generation and evaluation of an ensemble of wildland fire simulations. Int. J. Wildland Fire, 29: 160–173.
Armstrong G.W. 2004. Sustainability of timber supply considering the risk of wildfire. For. Sci. 50: 626–639.
Baltzer J.L., Day N.J., Walker X.J., Greene D., Mack M.C., Alexander H.D., et al. 2021. Increasing fire and the decline of fire adapted black spruce in the boreal forest. Proc. Natl. Acad. Sci. U.S.A. 118: e2024872118.
Barrette J., Thiffault E., Saint-Pierre F., Wetzel S., Duchesne I., Krigstin S. 2015. Dynamics of dead tree degradation and shelf-life following natural disturbances: can salvaged trees from boreal forests ‘fuel’ the forestry and bioenergy sectors? Forestry, 88: 275–290.
Bernier P.Y., Gauthier S., Jean P.-O., Manka F., Boulanger Y., Beaudoin A., Guindon L. 2016. Mapping local effects of forest properties on fire risk across Canada. Forests, 7: 157.
Bouchard M., Aquilué N., Périé C., Lambert M.-C. 2019. Tree species persistence under warming conditions: a key driver of forest response to climate change. For. Ecol. Manage. 442: 96–104.
Boucher J., Hébert C., Bauce E. 2020. A flexible approach for predicting and mapping postfire wood borer attacks in black spruce and jack pine forests using the differenced normalized burn ratio (dNBR). Can. J. For. Res. 50: 880–889.
Boulanger Y., Pascual J. 2021. Boreal forests will be more severely affected by projected anthropogenic climate forcing than mixedwood and northern hardwood forests in eastern Canada. Landscape Ecology 36: 1725–1740.
Boulanger Y., Parisien M.-A., Wang X. 2018. Model-specification uncertainty in future area burned by wildfires in Canada. Int. J. Wildland Fire, 27: 164–175.
Boulanger Y., Pascual J., Bouchard M., D'Orangeville L., Périé C., Girardin M.P. 2022. Multi-model projections of tree species performance in Quebec, Canada under future climate change. Global Change Biol. 28: 1884–1902.
Boychuk D., Martell D.L. 1996. A multistage stochastic programming model for sustainable forest-level timber supply under risk of fire. For. Sci. 42: 10–26.
Cardil A., Lorente M., Boucher D., Boucher J., Gauthier S. 2019. Factors influencing fire suppression success in the province of Quebec (Canada). Can. J. For. Res. 49: 531–542.
Chabot M., Blanchet P., Drapeau P., Fortin J., Gauthier S., Imbeau L. et al., 2009. Le feu en milieu forestier. Manuel de foresterie. 2e éd. Éditions Multimondes, Québec, QC. pp. 1037–1090.
Colombo S.J., Chen J., Ter-Mikaelian M.T., McKechnie J., Elkie P.C., MacLean H.L., Heath L.S. 2012. Forest protection and forest harvest as strategies for ecological sustainability and climate change mitigation. For. Ecol. Manage. 281: 140–151.
Couillard P.-L., Bouchard M., Laflamme J., Hébert F. 2022. Zonage des régimes de feux du québec méridional. Mémoire de recherche forestière no 189. Gouvernement du Québec, Ministère des Forêts, de la Faune et des Parcs, Direction de la recherche forestière, Québec, QC. 23 p.
Cumming S.G., Armstrong G.W. 2004. Divided landbase, overlapping tenures, and fire risk. For. Chron. 80: 478–484.
Cyr D., Splawinski T.B., Puigdevall J.P., Valeria O., Leduc A., Thiffault N., et al. 2022a. Mitigating post-fire regeneration failure in boreal landscapes with reforestation and variable retention harvesting: at what cost? Can. J. For. Res. 52: 568–581.
Cyr F.-X., Wyatt S., Hébert M. 2022b. Power-Sharing between the Cree and Québec Governments in Eeyou Itschee (Québec, Canada): Sovereignties, Complexity, and Equity under the Adapted Forestry Regime of the Paix Des Braves Int. For. Rev. 24, 345–359.
D'Orangeville L., Houle D., Duchesne L., Phillips R.P., Bergeron Y., Kneeshaw D. 2018. Beneficial effects of climate warming on boreal tree growth may be transitory. Nat. Commun. 9: 1–10.
Fick S.E., Hijmans R.J. 2017. WorldClim 2: new 1-km spatial resolution climate surfaces for global land areas. Int. J. Climatol. 37: 4302–4315.
Gauthier S., Bernier P.Y., Boulanger Y., Guo X., Guindon L., Beaudoin A., Boucher D. 2015a. Vulnerability of timber supply to projected changes in fire regime in Canada's managed forests. Can. J. For. Res. 45: 1439–1447.
Gauthier S., Raulier F., Ouzennou H., Saucier J.-P. 2015b. Strategic analysis of forest vulnerability to risk related to fire: an example from the coniferous boreal forest of Quebec. Can. J. For. Res. 45: 553–565.
Girardin M.P., Ali A.A., Carcaillet C., Blarquez O., Hély C., Terrier A., et al. 2013. Vegetation limits the impact of a warm climate on boreal wildfires. New Phytol. 199: 1001–1011.
Hanes C.C., Wang X., Jain P., Parisien M.-A., Little J.M., Flannigan M.D. 2019. Fire-regime changes in Canada over the last half century. Can. J. For. Res. 49: 256–269.
James P.M., Fortin M.-J., Fall A., Kneeshaw D., Messier C. 2007. The effects of spatial legacies following shifting management practices and fire on boreal forest age structure. Ecosystems, 10: 1261–1277.
Jobidon R., Bergeron Y., Robitaille A., Raulier F., Gauthier S., Imbeau L., et al. 2015. A biophysical approach to delineate a northern limit to commercial forestry: the case of Quebec's boreal forest. Can. J. For. Res. 45: 515–528.
Leduc A., Bernier P., Mansuy N., Raulier F., Gauthier S., Bergeron Y. 2015. Using salvage logging and tolerance to risk to reduce the impact of forest fires on timber supply calculations. Can. J. For. Res. 45: 480–486.
Mansuy N., Barrette J., Laganière J., Mabee W., Paré D., Gautam S., et al. 2018. Salvage harvesting for bioenergy in Canada: from sustainable and integrated supply chain to climate change mitigation. WIREs Energy Environ. 7: e298.
Marchal J., Cumming S.G., McIntire E.J. 2020. Turning down the heat: vegetation feedbacks limit fire regime responses to global warming. Ecosystems, 23: 204–216.
Nappi A., Drapeau P., Savard J.-P. 2004. Salvage logging after wildfire in the boreal forest: is it becoming a hot issue for wildlife? For. Chron. 80: 67–74.
Paradis G., LeBel L., D'Amours S., Bouchard M. 2013. On the risk of systematic drift under incoherent hierarchical forest management planning. Can. J. For. Res. 43: 480–492.
Pelletier G., St-Onge J., Bordeleau P., De Rainville P., Bart F., Aubin É., et al. 2009. Classification des peuplements forestiers en tant que combustibles selon la méthode canadienne de prévision du comportement des incendies de forêt. Gouvernement du Québec, Ministère des Ressources naturelles et de la Faune, Québec, QC, 56 p.
Périé C., de Blois S., Lambert M., Casajus N. 2014. Effets anticipés des changements climatiques sur l'habitat des espèces ligneuses au Québec. Mémoire de recherche forestière no. 173. Gouvernement du Québec, Ministère des Ressources naturelles, Direction de la recherche forestière, Québec, QC. 46 p.
Podur J., Wotton B.M. 2011. Defining fire spread event days for fire-growth modelling. Int. J. Wildland Fire, 20: 497–507.
Quinby P. 1987. An index to fire incidence. Can. J. For. Res. 17: 731–734.
R Core Team. 2021. R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. Available from https://www.R-project.org/.
Raulier F., Dhital N., Racine P., Tittler R., Fall A. 2014. Increasing resilience of timber supply: how a variable buffer stock of timber can efficiently reduce exposure to shortfalls caused by wildfires. For. Policy Econ. 46: 47–55.
Reed W.J., Errico D. 1986. Optimal harvest scheduling at the forest level in the presence of the risk of fire. Can. J. For. Res. 16: 266–278.
Régnière J., St-Amant R., Béchard A. 2013. BioSIM 10–user's manual. Canadian Forest Service, Québec, Canada.
Rheault H. 2013. Rendement soutenu. In Bureau du forestier en chef: Manuel de détermination des possibilités forestières 2013–2018. Gouvernement du Québec, Roberval, QC. pp. 19–21.
Saint-Germain M., Greene D.F. 2009. Salvage logging in the boreal and cordilleran forests of Canada: integrating industrial and ecological concerns in management plans. For. Chron. 85: 120–134.
Savage D.W., Martell D.L., Wotton B.M. 2010. Evaluation of two risk mitigation strategies for dealing with fire-related uncertainty in timber supply modelling. Can. J. For. Res. 40: 1136–1154.
Stocks B., Martell D.L. 2016. Forest fire management expenditures in Canada: 1970–2013. For. Chron. 92: 298–306.
Terrier A., Girardin M.P., Périé C., Legendre P., Bergeron Y. 2013. Potential changes in forest composition could reduce impacts of climate change on boreal wildfires. Ecol. Appl. 23: 21–35.
Van Wagner C. 1983. Simulating the effect of forest fire on long-term annual timber supply. Can. J. For. Res. 13: 451–457.
Van Wagner C. 1987. Development and structure of the Canadian forest fireweather index system. In Forestry Technical Report 35. For. Tech. Rep. Can. For. Serv. 35p.
Wang X., Studens K., Parisien M.-A., Taylor S.W., Candau J.-N., Boulanger Y., Flannigan M.D. 2020. Projected changes in fire size from daily spread potential in Canada over the 21st century. Environ. Res. Lett. 15: 104048.
Watanabe S., Hajima T., Sudo K., Nagashima T., Takemura T., Okajima H., et al. 2011. MIROC-ESM 2010: model description and basic results of CMIP5-20c3m experiments. Geosci. Model Dev. 4: 845–872.
White J.W., Rassweiler A., Samhouri J.F., Stier A.C., White C. 2014. Ecologists should not use statistical significance tests to interpret simulation model results. Oikos, 123: 385–388.
Wotton B., Flannigan M., Marshall G. 2017. Potential climate change impacts on fire intensity and key wildfire suppression thresholds in Canada. Environ. Res. Lett. 12: 095003.

Supplementary material

Supplementary Material 1 (DOCX / 14.7 KB).
Supplementary Material 2 (DOCX / 113 KB).
Supplementary Material 3 (DOCX / 23.2 KB).

Information & Authors


Published In

cover image Canadian Journal of Forest Research
Canadian Journal of Forest Research
Volume 53Number 2February 2023
Pages: 119 - 132


Received: 21 July 2022
Accepted: 7 October 2022
Accepted manuscript online: 21 October 2022
Version of record online: 17 January 2023

Data Availability Statement

The R code used to generate the data analyzed in this study is available at https://github.com/QLDM-Organization. Primary data comes from government repositories mentioned in the manuscript.

Key Words

  1. forest landscape modelling
  2. harvest level
  3. wildfire
  4. climate change
  5. salvage logging



Département des sciences du bois et de la forêt, Pavillon Abitibi-Price, Université Laval, 2405 rue de la Terrasse, Québec, QC G1V 0A6, Canada
Author Contributions: Conceptualization, Formal analysis, Methodology, Writing – original draft, and Writing – review & editing.
Département de Science et technologie, TÉLUQ, 5800, rue Saint-Denis, Montréal, QC H2S 3L5, Canada
Forest Science and Technology Centre of Catalonia (CTFC), Solsona, Lleida, 25280, Spain
Author Contributions: Conceptualization, Formal analysis, Methodology, Validation, Visualization, Writing – original draft, and Writing – review & editing.
Département de Science et technologie, TÉLUQ, 5800, rue Saint-Denis, Montréal, QC H2S 3L5, Canada
Author Contributions: Formal analysis, Funding acquisition, Methodology, Writing – original draft, and Writing – review & editing.
Jonathan Boucher
Laurentian Forestry Centre, Canadian Forest Service, Natural Resources Canada, Québec, QC, G1V 4C7, Canada
Author Contributions: Formal analysis, Methodology, Writing – original draft, and Writing – review & editing.
Natural Resources Canada, Canadian Forest Service, Northern Forestry Centre, Edmonton, AB, T6H 3S5, Canada
Author Contributions: Formal analysis, Methodology, Writing – original draft, and Writing – review & editing.
Marc-André Parisien served as an Associate Editor at the time of manuscript review and acceptance; peer review and editorial decisions regarding this manuscript were handled by Yu Wei.

Author Contributions

Conceptualization: MB, NA, EF
Formal analysis: MB, NA, JB, EF, MAP
Funding acquisition: MB, EF
Methodology: MB, NA, JB, EF, MAP
Validation: NA
Visualization: NA
Writing – original draft: MB, NA, JB, EF, MAP
Writing – review & editing: MB, NA, JB, EF, MAP

Competing Interests

The authors declare there are no competing interests.

Funding Information

Funding was provided by the ministère de la Forêt, de la Faune et des Parcs du Québec as well as Québec's Fonds vert (Plan d'action 2013–2020 sur les changements climatiques).

Metrics & Citations


Other Metrics


Cite As

Export Citations

If you have the appropriate software installed, you can download article citation data to the citation manager of your choice. Simply select your manager software from the list below and click Download.

There are no citations for this item

View Options

View options


View PDF

Get Access

Login options

Check if you access through your login credentials or your institution to get full access on this article.


Click on the button below to subscribe to Canadian Journal of Forest Research

Purchase options

Purchase this article to get full access to it.

Restore your content access

Enter your email address to restore your content access:

Note: This functionality works only for purchases done as a guest. If you already have an account, log in to access the content to which you are entitled.





Share Options


Share the article link

Share on social media

Cookies Notification

We use cookies to improve your website experience. To learn about our use of cookies and how you can manage your cookie settings, please see our Cookie Policy.