Potential changes in monthly fire risk in the eastern Canadian boreal forest under future climate change

We examined the fireclimate relationship in the Waswanipi area, central Quebec. This territory is about 15000km² and is situated in the western black spruce feather moss bioclimatic subdomain (Robitaille and Saucier 1998) (Fig.1). It is composed of continuous boreal forest, where stands are dominated by black spruce (Picea mariana (Mill.) BSP) growing on thick glaciolacustrine tills originating from the Ojibway proglacial lake (Robitaille and Saucier 1998). Jack pine (Pinus banksiana Lamb.) is abundant on coarse-textured soils, while mixed stands of aspen (Populus tremuloides Michx), white birch (Betula papyrifera Marsh.), white spruce (Picea glauca (Moench) Voss), and balsam fir (Abies balsamea (L.) Mill.) can be found on upland till soils (Rowe 1972). To model the fireclimate relationship we calculated regression analyses for two territories: the Waswanipi area (15000km²) and the western black spruce feather moss bioclimatic subdomain (which encloses the Waswanipi area). We used two different territories, as the regression approach used to establish the historical fireclimate relation usually generates better results for larger territories owing to a smaller number of years with no fire activity.

According to the 19712000 climate normals from the Chapais weather station (4947N, 7451W, altitude 396.20m), the area has 1235 degree-days per year above 5C and around 961mm of precipitation, with one-third falling as snow. February is the coldest month and July the warmest, with a daily mean temperature of16.6 and 16.3C, respectively (Environment Canada 2004).

We used the provincial fire data provided by the ministre des Ressources naturelles et de la Faune du Qubec for 19732007. This database reports all fires from all origins (from lightning as well as from human activities). We selected all fires>10ha for our analyses to eliminate very small fires that were not contributing to the area burned. These data encompass the period during which systematic fire detection in the restricted fire management zone was made by detection planes, which only began in the late 1960s (Blanchet 2003).

We summarized the total area burned and the total number of fires for each day of the fire season considered. For the fireclimate analyses, we set the fire season from 1May to 31 August, since this period accounted for more than 99% of the 19732007 annual area burned in the Waswanipi area.

No transformation allowed normality to be achieved at the monthly step because of the high frequency of months with no fire. Because of these numerous nonfire months, the monthly area burned variable behaves rather like a binomial variable. This property oriented our analyses to logistic regressions to model the monthly fire risk rather than to linear multiple regressions to model the total area burned in a month (see Flannigan et al. 2005). The fire risk is defined here as the monthly probability of having a large (>500ha) or very large (>2000ha) burned area.

The FWI system (Van Wagner 1987) is used across Canada to evaluate the daily fire risk based on local daily weather data. It consists of moisture codes and fire behaviour indices calculated from the daily temperature, 24 h accumulated precipitation, relative humidity, and wind speed (Fig.2 Van Wagner 1987 Wotton 2008). The three moistures codes track moisture in different levels of the forest floor. The Fine Fuel Moisture Code (FFMC) evaluates the moisture content in small, readily consumed fuels on the surface of the forest floor. The FFMC indicates the ease of ignition and flammability of fine fuels. The Duff Moisture Code (DMC) measures the moisture content of the moderate organic layers of the forest floor, where litter begins to decay. It provides an estimate of consumption of these duff layers and of medium-size woody debris. The Drought Code (DC) is an indicator of the moisture content of deep layers of the forest floor and of large down and dead woody debris on the forest floor. The Buildup Index (BUI) is a combination of the DMC and DC. It evaluates the potential fuel available for surface fuel consumption by the passing fire front. The Initial Spread Index (ISI) is a combination of FFMC and wind speed it evaluates the potential rate of spread of a fire. Finally, the FWI is a combination of the BUI and the ISI, and the Daily Severity Rating (DSR) is essentially a logarithmic transformation of the FWI. Here, we evaluated whether one or several of these FWI components could provide good estimates of the regional fire activity in central Quebec, Canada.

Daily temperature, precipitation, wind speed, and relative humidity from 1 May to 31 August for 19602007 were extracted using BioSim 9 for the centre of the Waswanipi area (Rgnire and Saint-Amant 2008). For each weather variable, BioSim 9 provided the mean daily value from the three nearest weather stations based on an inverse weighted spatial interpolation. Depending on the year considered, the selected weather stations were situated between 26km (1979, weight 75.5%) and 408km (1993, weight 17.6%) with a mean distance of 134 64km to the center of the study area. Values were adjusted for distance and elevation. We calculated annual and monthly mean, minimum and maximum values for each weather variable and FWI components from the daily 19602007 data set. We used these potential predictors (Table1) in the linear regression analysis to model the log-transformed annual area burned and the log-transformed annual number of fires, and in the logistic regression analysis to model the monthly fire risk.

The historical link between annual fire activity (log-transformed area burned and log-transformed number of fires), weather variables, and FWI components was calculated using multiple linear regressions (Flannigan et al. 2005 Bergeron et al. 2006) using SAS version 9.1 (SAS Institute Inc. 2000). The structure of the linear models iswhere Y is the log-transformed annual area burned (or the log-transformed annual number of fires) in the territory investigated (Waswanipi area or western black spruce feather moss subdomain), is a constant, _i is the regression coefficient estimate of the variable x_i selected among the candidate explanatory variables (Table1), and _j is the standard error associated with the model.

Multiple regression models were compared and evaluated using different statistics in addition to the adjusted coefficient of determination and the standard error of the estimates. The F ratio was used to evaluate the predictive capability of the model, considering the number of variables selected. The F ratio is obtained by dividing the explained variance by the unexplained variance. We used the Akaike information criterion corrected for small sample sizes (AICc see Mazerolle 2006) to select the best model among three sets of candidate explanatory variables: weather variables only, FWI components only, and weather variables and FWI components taken together. We also reported the Bayesian information criterion (BIC). The BIC is another information criterion for model selection that we used to compare with the result obtained using the AICc. As underlined by Girardin et al. (2008), the empirical modelling provides evidence of statistical association between forest fire activity and climate, but only suggests possible biological relations (Arbaugh and Peterson 1989). Using a threshold of p< 0.05, there is typically a 5% chance that some of the 24 potential predictor variables may be retained in spite of weak biological justification or because of spurious relationships. Selected variables that had illogical ecological relationships to the model (e.g., positive influence of relative humidity or negative influence of FWI on fire activity) were thus manually removed from the candidate explanatory variables, and models were recalculated until all variables taken by the model were ecologically sound.

The stability of the regression model was tested using a split-sample calibrationverification scheme. Regressions were calculated for the entire time period of 19732007. Variables selected were then entered in a complete regression over two subperiods: 19731990 and 19912007. The regression coefficients estimated for one subperiod (calibration period) were then applied to the selected variables over the other subperiod (verification period) (see Girardin 2007). The strength of the relationship between the regression models and observations was measured using Pearsons correlation coefficients.

When examined at the monthly time step, the area burned produced many months without fire, so this variable was treated as a binomial variable. We used logistic regression analyses to model the monthly fire risk, using monthly weather variables and FWI components as potential explanatory variables. Monthly fire risk is defined here as the probability to have a month with area burned over 500 or 2000ha. The structure of the logistic model iswhere P_FIRE is the monthly fire risk, is a constant, _i is the regression coefficient, and x_i is the variable selected among the candidate explanatory variables (Table1). We reported the same information criteria as for the previous regression analyses. Logistic regressions were calculated using the SAS version 9.1 (SAS Institute Inc. 2000) with a stepwise forward selection procedure. As for the linear regressions, we manually deleted aberrant variables.

To verify the stability of the regression coefficients, we calculated the logistic regression with the variables selected for the 19732007 period for two subperiods: from May 1973 to June 1990, and from July 1990 to August 2007. We evaluated submodels using the percentage of concordance (i.e., the proportion of events and nonevents correctly predicted by the submodels).

To anticipate climate conditions under climate change, we used daily 19612100 outputs from the CRCM. The CRCM contains a numerical description of the physical processes within the climate system. It is a limited-area regional climate model that uses boundary conditions provided by the Canadian GCM 2 (CGCM2) simulation (Plummer et al. 2006 Laprise 2008). While the Intergovernmental Panel on Climate Change (IPCC) recommends using a multimodel and multiscenario approach to determine an envelope of possible future impacts resulting from climate change (Bernstein et al. 2007), only one regional climate model and one climate change scenario (A2) with two realizations (simulations with different boundary conditions) were available for this study. The A2 scenario corresponds to a status quo situation (Nakicenovic et al. 2007): greenhouse gas emissions are continuing to rise at the current rate in a very heterogeneous world with a rapid population growth and regional-oriented economic development. It is generally considered the most pessimistic climate change scenario (Bernstein et al. 2007). Moreover, the examination of recent climate observations suggested that the A2 scenario used by the IPCC tends to underestimate changes in atmospheric CO₂ concentration and temperature (Rahmstorf et al. 2007), so the A2 scenario may actually provide conservative estimates of future climate conditions. We have chosen to use the CRCM rather than several GCMs with several scenarios (see Drever et al. 2009), as we were interested in a relatively small territory and wanted more accurate estimates of future climate conditions at the regional scale, which are more appropriate for forest management.

Continuous 19612100 daily weather data were obtained for the two available A2 realizations for the 21 CRCM cells covering the Waswanipi area (Fig.1). Daily minimum and maximum temperatures were also extracted to calculate daily relative humidity at noon using the GoffGratch equation (Goff and Gratch 1946). We calculated the mean value of each meteorological variable for the 21 CRCM cells of each realization. Then we calculated the mean value of these two realizations. Monthly and annual mean and maximum were calculated using these mean daily values. The 19612007 median of temperature, wind speed, and relative humidity were adjusted according to the median of the corresponding historical series (by subtracting the difference between both medians). The rain series did not necessitate adjustments (the medians were close to each other). Then we used the adjusted daily meteorological variables to calculate the daily values of FWI components.

When compared to other months of the fire season, June accounted for the highest proportion of the annual area burned between 1973 and 2002 (42% for the Waswanipi area and 55% for the western black spruce feather moss bioclimatic subdomain). For this reason, we examined the CRCM meteorological variables and FWI components for this month in particular. Annual area burned and annual number of fires were log-transformed to achieve normality according to a one-sample KolmogorovSmirnov test. Potential temporal trends were tested using a linear regression with time as predictor. To examine future trends in weather variables and FWI components, we extracted the mean and maximum June values, as this month accounted for the largest part of the annual area burned. We tested linear temporal trends in meteorological variables and FWI components using a simple linear regression with time as a predictor in Sigmaplot 11 (Systat Software, Inc. 2006). Normality of mean and maximum series was tested using a one-sample KolmogorovSmirnov test. Maximum rain and mean DSR failed this normality test.

We substituted historical weather variables and FWI components with those from CRCM in the best regression models previously identified to estimate 19612100 fire activity. We tested linear temporal trends in 19612100 log-transformed annual area burned and annual number of fires using linear regression with time as a predictor. Then we calculated mean values for 19752005 (1 CO₂), 20302060 (2 CO₂), and 20702100 (3 CO₂) to calculate the rates of change in fire activity (2 CO₂ / 1 CO₂ and 3 CO₂ / 1 CO₂). We calculated rates of change in the monthly fire risk between the three reference periods for all months taken together, as well as for each month taken separately to verify whether future monthly distribution of fire risk will change under climate change.

Linear regression models of log-transformed annual area burned and log-transformed number of fires had adjusted R² from 0.20 to 0.45. The adjusted R² values were generally higher for the number of fires than for the area burned (Table2). They were equivalent for the two territories considered. According to the lowest AICc, the best models always corresponded to the combination of weather variables and FWI components (Table2). Verification analyses indicated good correlations between observed data and submodels (Pearson correlation coefficients>0.5, Table3).

Monthly fire risk for the western black spruce feather moss subdomain was best predicted by a combination of temperature and a FWI component (BUI or FFMC) (Table4). Conversely, the Waswanipi fire risk was best explained by BUI only. For all models, concordance between modeled and observed data was around 80%. The verification displayed percentages of concordance>71% for all submodels, with percentages being more stable for the Waswanipi area (Table5).

A positive linear trend was identified for the maximum and mean temperature (p < 0.001) and for maximum and mean relative humidity (p < 0.001, both) (Fig.3). These trends were, however, less evident in the FWI components, as only maximum FWI and maximum DSR displayed a positive linear trend with time (p = 0.030 and p = 0.017, respectively) (Fig.4).

A linear trend was detected for the 19612100 log-transformed annual area burned (R² = 0.13, p < 0.001), while no trend was detected for the annual number of fires (p = 0.47) (Fig.5). The increasing trend in 19612100 log-transformed annual area burned was observable on the 10year moving average, but the standard error associated with the model indicated a high interannual variation. When examining the rates of change between future and current periods, we observed a general increase in Waswanipi fire activity. This increase is more pronounced for the 3 CO₂ period (20702100) than for the 2 CO₂ period (20302060). By 2100, the annual area burned could increase by 7%, and the monthly fire risk could increase by 30%, with the highest increases observed in July (70%) and August (100%) (Table6, Fig.6). However, the annual number of fires did not exhibit changes, and the May fire risk displayed a decrease (20%) (Table6).

The monthly fire risk increased for June, July, and August, while the May fire risk slightly decreased between current and future time periods. These results suggest that spring fires may be less frequent in the future, while the fire risk would increase considerably in July and August in the future. This contrasts with the study of Wotton and Flannigan (1993), who suggested an earlier start of the fire season under climate change. Our study did not consider April, September, and October, because the area burned and the number of fire events during these months were too low to allow statistical analysis. However, these months may play an important role in the future fire season (Wotton and Flannigan 1993 Nitschke and Innes 2008b). The change in the fire risk distribution across the fire season represents a major challenge for the fire management strategy in Quebec. While June is currently the main fire month, the future fire season would present a prolonged high fire risk from June to August. This suggests a prolonged effort of fire management, implying an increasing investment in fire control efforts. The fire management agency may more often face situations where its fire suppression capacity may be overwhelmed. The decrease observed in May could be linked to the increase in winter precipitation associated with climate change (Christensen et al. 2007), as the longer snow melt could delay the drying of forest fuels in spring.

According to our results, the increase in fire-weather conditions (maximum FWI) as well as in fire activity (area burned) is relatively modest when compared to the results of other studies that used GCM data over the same territory (Flannigan et al. 2005 Bergeron et al. 2006). Future climate change would trigger fire-weather conditions more favourable to forest fires, and annual area burned will continue to increase when compared to the 19752005 reference period. Using the same linear regression approach, Flannigan et al. (2005) anticipated an increase of about 50%100% in the area burned observed in the reference period 19751995 for the 20802100 period in the eastern part of the Boreal Shield ecozone. Bergeron et al. (2006) estimated that fire rate (annual proportion of area burned) would increase by 30% for 20802100 when compared to the rate for 19402003 in the Waswanipi area. The standard error associated with the fire activity estimates along with the use of a single climate change scenario (A2) constitute the main limitations of the interpretation of an increase in fire activity under climate change in the Waswanipi area. The main challenge in developing an ecosystem-based forest management plan is less this slight increase in fire activity and more the interannual variation in fire activity estimates. Yet our results suggest that fire will remain an important constraint for forest management in the context of climate change.

Wind is a key variable controlling fire activity, notably area burned and fire behaviour. This is confirmed by the selection of the maximum wind speed in the annual area burned models for the western black spruce feather moss subdomain and for the Waswanipi area. However, no linear trend was detected in the 19612100 CRCM wind speed. Interestingly, the clear increase observed in temperature data over the 19612100 period is less perceptible in the FWI components, as only maximum FWI displayed a positive linear trend with time. Maximal DSR also displayed the same positive linear trend, as this index is essentially a log-transformed FWI. As speculated in previous work (Bergeron and Flannigan 1995), our results suggest that the increase in temperature alone is not a sufficient condition to lead to weather conditions favourable to forest fires. We postulate that in a certain measure, the increase in temperature would be compensated by the increase in relative humidity. While this might appear contradictory to the increasing linear trend observed in maximum FWI, this could be reconciled by considering the seasonality aspect. The variance explained by our linear models (R² from 0.36 to 0.48, using one to three predictors) is comparable to that of Flannigan et al. (2005 R² from 0.36 to 0.64, using one to four predictors).