Fire severity and plant productivity recovery in a mixed grass prairie wildfire driven by extreme winds

Study area

The study area is located near the provincial border of Saskatchewan and Alberta in the Mixed Prairie of the Northern Great Plains (Coupland 1961). On 17 October 2017, two large fires ignited in southeast Alberta near Acadia Valley (Empress Fire; 51.114 N, −110.121 W) and Hilda (Hilda Fire; 50.535 N, −109.992 W). Both fires crossed the provincial border into southwest Saskatchewan, with approximately 28,000 ha of land burned, 8000 of which was pastureland. High wind speeds and a moderate drought contributed to extreme fire behaviour (Briere 2017). Peak winds on the day of ignition exceeded 120 km h⁻¹ (Glen and Briere 2017); the dominant wind direction up to approximately 5 pm was 224° followed by a shift to approximately 277° (Environment and Climate Change Canada; https://climate.weather.gc.ca/climate_data/hourly_data_e.html) causing a distinct ‘finger’ pattern to the fire footprint.

The Alberta sites are classified as Dry Mixedgrass Natural Subregion, while Saskatchewan sites are classified as Dry Mixed Grass Ecoregion (Coupland 1961; Thorpe 2014). In the Dry Mixed Grass Ecoregion of Saskatchewan, mean annual precipitation and temperature are 250–350 mm and 3.5°C, respectively (Acton et al. 1998). In the Dry Mixedgrass Natural Subregion of Alberta, mean annual precipitation and temperature is 333 mm and 4.2°C, respectively (Natural Regions Committee 2006). Vegetation at the study sites is characterised by xeric graminoid species, such as Bouteloua gracilis, Pascopyrum smithii, Koeleria macrantha, and Hesperostipa comata, forbs such as Selaginella densa, Artemisia ludoviciana, and Phlox hoodii, and shrubs such as Rosa arkansana (Natural Regions Committee 2006; Thorpe 2014). The topography of the study area is composed of knolls, gentle and steep slopes, as well as bordering a steep river valley of the South Saskatchewan River. The elevation of the study area ranges from 590 m to 880 m above sea level.

Satellite images and image processing

Satellite imagery was retrieved from the United States Geological Survey Earth Explorer Database for Sentinel-2 Level 1C. Sentinel-2 Level 1C was chosen for its high spatial resolution and for already being atmospherically and radiometrically corrected. Criteria for image selection included: less than 30% cloud cover with burn scars not covered by cloud, both fires needed to be within one image, image capture date needed to occur as close as possible to peak biomass production (mid-July) and as close as possible to before and after the fires. Five images were chosen for use (Table 1). Images include three pre-fire images. Mid-July is the typical peak in growth in this ecosystem, therefore the 15 July 2017 image captures peak green biomass pre-fire. Senescence is normally well advanced by late September in this ecosystem, therefore 28 September 2017 image captures cured biomass. The immediate post-fire image was collected approximately 11 days after ignition when all active fire had been suppressed.

Upon image selection, bands 2–4, 8, 10, and 11 for each image were transferred into one file in Geomatica (PCI Geomatics 2018). The immediate post-fire (F) image was assessed to determine the bands with the highest visual contrast between burned and unburned areas (RGB: SWIR2, SWIR1, and Red) (Fig. 1). Visual interpretation of these bands was used to draw two irregular polygons outlining the Hilda and Empress fires (Fig. 1). These irregular polygons were combined and used for clipping all five Sentinel-2 images and ancillary data down to manageable file sizes.

Fig. 1.

False-colour (RGB: SWIR2, SWIR1, and Red) composite of both the Hilda and Empress wildfires on 23 October 2017. The original image was cut into an irregular polygon to remove any unnecessary image. The two pink irregular polygons were then used later in the study when the burn scars were undetectable. The northern burn scar (yellow) is the Empress fire, while the southern burn scar is the Hilda fire. Fires originated near the westernmost point of each burn scar.

After each image was cut down to a manageable size the SWIR 1 and SWIR 2 bands were resampled to 10 m by 10 m from their original 20 m by 20 m resolution. This was done to ensure that when spectral indices were calculated, the resolution of all pixels matched. Two spectral indices - the normalised differenced vegetation index (NDVI) and the bare soil index (BSI) - were calculated. NDVI measures the greenness of vegetation through the absorption of red light by the vegetation chlorophyll content while simultaneously reflecting near infrared light off of turgid cell walls (Tucker 1979; Smith et al. 2015). BSI measures the amount of bare ground through the absorption of infrared and red light, and the reflectance of near infrared and blue light (Rikimaru et al. 2002; Diek et al. 2017; Nguyen et al. 2021).

Environmental data

Four types of ancillary environmental data were used. Soil Landscapes of Canada (SLC) version 3.2 polygons of a scale of 1:1,000,000 contained information such as soil name, soil texture components, and soil chemical composition. The Canadian Land Inventory (CLI) maps at a scale of 1:1,000,000 describe the agricultural capability of the landscape based on soil classification, texture, presence of stones, and topography. The Digital Elevation Map (DEM) of Canada at a scale of 1:50,000 described topography. A calculation of solar heat load was made to estimate the amount of heat a slope receives relative to aspect (McCune and Keon 2002). Aspect with respect to the wind (whether a location was on a slope facing into the fire, flanking a fire, or on a reverse slope away from the fire) was calculated based on dominant wind direction. Fire movement prior to 5 pm of the day of fire was assumed to be driven by winds from 224° while the fingers extending to the southwest were assumed to driven by winds from 277°.

Both the SLC polygons and the CLI maps were georeferenced and projected to World Geodetic System (WGS) 1984 (QGIS Development Team 2015). To ensure the CLI maps and the SLC polygons were correctly georeferenced, a preloaded QGIS map of Canada was used to ensure the border and corners of the maps and polygons aligned correctly with Alberta and Saskatchewan. The SLC polygons were clipped to the study area. Upon georeferencing, the CLI maps had individual polygons drawn to reflect the different land use areas (QGIS Development Team 2015). Lastly, the Canada wide DEM was projected to WGS 1984 and clipped to the study area prior to computing aspect and slope.

Both fires covered just over 1 million pixels. To ensure the most efficient exportation of all data, the irregular polygon that outlined each fire was split along the provincial border. The four polygons were auto-populated with points (a maximum of 1 million points 10 m apart) to ensure every part of the polygon was sampled and that each pixel was only sampled once. Following auto population, the points were linked to the three layers of ancillary data using the join function. Lastly, both NDVI and BSI for all five Sentinel-2 images were linked to the points using Extract Multi-Values to Points. Once linking was complete the four attribute tables were exported.

Following data exportation, the fire severity index (FSI) and the recovery index (RI) were calculated. The FSI measured the change in Bare Soil Index (BSI) between fall 2017 (F17) and immediate post-fire (F) images. Before FSI calculation, negative BSI values for both F17 and F images were truncated to zero to ensure that only the effects of the fires were captured. Following the creation of FSI, it was found that some negative values FSI were produced. It was assumed that these negative values were caused by random error due to the large data set, thus a second truncation was performed to change negative FSI values to zero to ensure misinterpretation of FSI was avoided. RI was produced by assessing the overall change in NDVI between post-burn (PF) and pre-burn (PrF17) values. The log response ratio of NDVI (log(PFx + 1)/(PrF17 + 1)), was taken to ensure only the recovery (increased greenness) of the specific PF year being evaluated could be determined. Positive values created from the log response ratio indicate greater biomass than pre-fire conditions. Negative values form the log response ratio indicate a production deficit compared to pre-fire conditions.

Statistical analysis

The role of environment in determining fire severity index (FSI) were assessed using models with slope, aspect with respect to wind, PrB17 NDVI, and F17 NDVI as predictors. The role of environment and fire severity in determining grassland recovery index (RI) were assessed using models with soil order, land capability, solar heat load, and fire severity as predictors. Prior to analysis, correlations between FSI and RI and the environmental predictors were examined using R (Simko 2021; R Development Core Team 2022). Any variables including predictor variables with correlations greater than |0.6| were selected for subsequent analysis. Environmental variables that could not logically affect fire behaviour were excluded from analysis.

Generalised Additive Models (GAM) were used to examine how fire severity and recovery varied with respect to environmental gradients (Wood 2011; R Development Core Team 2022). GAM models were used because they permit complex non-linear relationships between predictor and response variables to be assessed. To test how slope aspect with respect to wind direction (categorical variable with three levels) impacted fire severity, a series of four GAM models were fit. These models included (1) the slope alone model, (2) aspect as a categorical fixed effect, (3) aspect as a categorical fixed effect and the smoothed slope term as separate effects, and (4) the continuous smoothed term with the categorical variables broken down into by = terms (effectively fitting a separate smoothed line for each aspect) (Table 2). Recovery was modelled using a similar approach with separate GAM models using Fire Severity and Heat Load index as continuous predictors, and models combining the continuous and categorical predictors with Soil Order and Land Class as categorical predictors (Table 3). In all cases, function gam from package mgcv was used to fit gam models (Wood 2006). Smoothed relationships with continuous predictor variables were fit using thin-plate regression splines; by = arguments were used when separate smoothed relationships were desired for each level of a categorical variable.

Fire severity

Fire severity, a measure of how much bare ground was exposed, was best predicted by pre-fire (September 2017) NDVI. Fire severity increased rapidly with pre-fire (F17 and PrB17) NDVI (Fig. 2a; Table 2). Peak fire severity occurred between 0 and 0.2 NDVI with a gradual decline to an NDVI of ~0.6 where a secondary peak of fire severity occurred. Similar patterns were observed with July 2017 (PrB17) NDVI, where a rapid increase in fire severity to an NDVI of ~0.35 was followed by a plateau in fire severity at higher NDVI.

Fig. 2.

Relationships between fire severity and (a) NDVI, (b) hill slope, and hill slope divided by aspect with respect to the wind. Fitted lines from GAM models are shown with shaded areas around the lines indicating the standard error of the fitted values.

Aspect had smaller effects on fire severity than NDVI. The relationships between slope and wind direction (aspect) created complex effects on fire severity (Fig. 2b; Table 2). Aspect with respect to wind was divided into three categories: facing upwind into fire ±45° of the dominant wind direction, ±45° downwind away from fire, and flanks. Aspect with respect to wind direction had limited impact on fire severity on slopes less than 10°. On slopes greater than 10° higher fire severity occurred on slopes facing into the wind, and lower severity on slopes facing away from the wind. There was little impact of slope on severity on the flanks.

Recovery

Recovery following fire, a measure of how much green regrowth occurred relative to the amount of pre-fire green biomass, was rapid; Year 1 NDVI was similar to pre-fire greenness (Fig. 3a; Table 3). Recovery was much lower in Year 3. Relationships between recovery and fire severity were complex with greater recovery detected above a fire severity of 0.2. Recovery was similar across heat load values indicating similar recovery rates on both hot south facing slopes and more mesic locations (Fig. 3b; Table 3). More fertile sites (Chernozemic and Vertisolic soils, and Land Class 3) had generally lower recovery than less fertile sites; this trend was most pronounced on sites with lower fire severity and in year 3 (Fig. 4; Table 3).

Fig. 3.

Relationships between (a) fire severity and (b) solar heat load for year one (blue line) and year three (red line) post-fire. Recovery is the log response relationship between ratio of NDVI for the year of recovery divided by pre-burn values; values greater than zero exceed pre-fire greenness and values less than zero have less greenness than pre-fire. Solar heat was calculated using McCune and Keon (2002) method which uses slope and aspect to derive sun exposure. Shaded areas indicate the standard error around the fitted values.

Fig. 4.

Relationships between fire severity and recovery for year one and year three post-fire as a function of site soil order (a, b) and agricultural land class (c, d). Recovery is the log response ratio of NDVI for in the year of recovery divided by pre-fire values greater than zero exceed pre-fire greenness, values less than zero have less greenness than pre-fire. Soil order and land class are proxies for site potential productivity.

Fire severity

This study is one of the first to examine the landscape factors that drive grassland fire severity under extreme (>120 km h⁻¹) wind conditions. Fire severity, or the amount of bare soil exposed by the fire, varied as a function of topographic and environmental conditions including pre-fire NDVI, slope, and aspect with respect to wind with NDVI as the dominant factor. We interpret the peak in fire severity at an NDVI of 0–0.2 and 0.2–0.4, to be driven by fuel loading as these NDVI values in grasslands are associated with areas of high litter content and minimal green vegetation, respectively (Xu et al. 2014). In September, senescence is normally well advanced in these grasslands, thus September NDVI should provide a reliable measure of available cured fuels. A peak in fire severity in areas of high fuel is consistent with other studies (Cruz et al. 2017; Kraaij et al. 2018; Wragg et al. 2018; García-Llamas et al. 2019).

We observed important interactions between slope and slope aspect with respect to wind direction that generated very high fire severity on slopes greater than 15° that faced into the wind (±45° of the dominant wind direction). Slope was unimportant, however, on sites oblique to the wind. We suspect that slopes facing into the wind had high severity due to both the extreme wind directly feeding oxygen to the flames causing a hotter fire (Thomas et al. 2017; Evers et al. 2022), and flames being forced parallel to the ground. When uphill slopes are greater than 24° the tip of the flame to the base of the flame lay parallel to the ground (McArthur 1967). There are no studies assessing the fire severity generated by parallel flames, though it can be hypothesised parallel flames would have a higher fire severity due to longer residence time as the entire flame length (tip to base) would expose a single spot to heat for longer. Parallel flames are typically associated with slopes greater than 24°, but when compounded with the extreme winds, parallel flames could occur at slopes as low as 15°. Many slopes greater than 15° are found in the study area.

The interactions between aspect with respect to wind and slope driving fire severity are novel. Studies that assess aspect and slope are commonly done theoretically or within controlled environments, such as prescribed fire, where fire intensity (energy released) is typically quantified rather than fire severity (McArthur 1967; Weise and Biging 1996; Nelson 2002; Thomas et al. 2017; Hall 2019). Fire severity is used in wildfire studies where direct measures of intensity are challenging to take (Evers et al. 2022). Both measures are typically assessed as independent factors, but when aspect with respect to wind and slope are considered simultaneously, wind either overrides the effects of slope, or is additive to the effects of slope (Weise and Biging 1996; Nelson 2002; Hall 2019; Evers et al. 2022). The importance of aspect with respect to wind in our study demonstrates the importance of studying fire behaviour under both moderate and extreme conditions.

Recovery

Recovery varied as a function of fire severity and of environmental conditions including potential site productivity (soil order and land classification) and heat load. The heat load index showed the highest recovery on the coolest and warmest aspects. Overall recovery was lower in year three than in year one.

The fire severity–recovery relationship was similar for both soil orders and land classes. The four dominant soil orders and land use classifications in the study area vary in potential productivity (Regosols < Chernozems < Vertisols < Gleysols; Class 3 < Class 4 < Class 5 < Class 6) (Bedard-Haughn 2010; VandenBygaart 2010; Brierley et al. 2011; Pennock et al. 2011). Faster recovery might be expected on more fertile sites where soil moisture is less limiting and greater growth rates are possible, but such a pattern was not observed here. The lack of a site potential effect on recovery rate may be because each site type has its own characteristic native plant communities. These communities differ in species composition with more xeric tolerant species on less productive sites; all of these species however have extensive belowground structures that are well adapted to resprouting post-fire, particularly a dormant season fire (Coupland and Johnson 1965; Comfort and Lamb 2023).

Year one of recovery was similar to pre-fire greenness, congruent with the findings of Steiner et al. (2020). Both fires occurred in the late fall, which can cause lower productivity (recovery) the year following when compared to spring fires (Jardine et al. 2017). Recovery could also have been underestimated in year one because the bright green regrowth contrasted with dark soil. NDVI measures the total greenness of vegetation but can be disrupted by large amounts of bare ground, which can create a lower NDVI value (Smith et al. 2015; Steiner et al. 2020). Finally, lower precipitation occurred in June 2018 compared to June 2017. The lower average precipitation for June could have been an additional contributing factor to the lower-than-expected recovery rates for year one. Normal temperature and precipitation trends were seen in 2019, but in 2020 a season-long severe drought drastically reduced productivity (Grover 2024). Biomass production in grassland ecosystems is closely tied to precipitation (Sala et al. 1988; Huxman et al. 2004; Knapp et al. 2015; Kreyling et al. 2017), and while average productivity declined due to drought, the relationships between recovery and fire severity and site heat load did not change substantially. This lack of interactions demonstrates that variation in fire severity does not change the vulnerability of these grasslands to drought. This is a surprising finding given that litter is very important for soil moisture retention (Deutsch et al. 2010; Hilger and Lamb 2017), and increased fire severity implies reduced biomass carryover. In this case, the lack of interaction can likely be explained by the nearly complete consumption of litter by the fire observed at all sites surveyed on the ground (Grover 2024).

Recovery of vegetation was highest in both areas with a moderate heat load and areas with a high heat load. Recovery peaking under a high heat load could be attributed to a community composition of xeric species that have robust belowground perennial organs and can tolerate drought conditions (Natural Regions Committee 2006; Thorpe 2014). Recovery peaking under a moderate heat load could be attributed to the increased soil moisture content in those areas.

Post-fire grazing management decisions can be challenging for agricultural producers. Regulations limiting grazing post-fire (Bureau of Land Management 2007; Government of Alberta 2018) perpetuate the idea that fire within native grasslands ecosystems are not a historically naturally occurring process; yet grasslands evolved under the co-disturbance of fire and grazing (Anderson 2006). Given the landscape scale variation in severity and recovery demonstrated here, one-size-fits-all grazing prescriptions are not appropriate. Grazing animals are attracted to the higher nutritional content of recovering burned vegetation (Fuhlendorf et al. 2009; Vera-Velez et al. 2023), but areas with limited vegetation recovery may be damaged. Understanding the relationships between fire severity slope, and aspect can enable agricultural producers to manage their lands using techniques such as cross fencing or virtual fencing to reduce grazing pressure on areas with higher fire severity while taking advantage of nutritious recovering vegetation in areas of lower severity.