Fire refugia are robust across Western US forested ecoregions, 1986–2021

In the Western US, area burned and fire size have increased due to the influences of climate change, long-term fire suppression leading to higher fuel loads, and increased ignitions. However, evidence is less conclusive about increases in fire severity within these growing wildfire extents. Fires burn unevenly across landscapes, leaving islands of unburned or less impacted areas, known as fire refugia. Fire refugia may enhance post-fire ecosystem function and biodiversity by providing refuge to species and functioning as seed sources after fires. In this study, we evaluated whether the proportion and pattern of fire refugia within fire events have changed over time and across ecoregions. To do so, we processed all available Landsat 4–9 satellite imagery to identify fire refugia within the boundaries of large wildfires (405 ha+) in 16 forested ecoregions of the Western US. We found a significant change in % refugia from 1986–2021 only in one ecoregion—% refugia increased within fires in the Arizona/New Mexico Mountain ecoregion (AZ/NM). Excluding AZ/NM, we found no significant change in % refugia across the study area. Furthermore, we found no significant change in mean refugia patch size, patch density, or mean distance to refugia. As fire size increased, the amount of refugia increased proportionally. Evidence suggests that fires in AZ/NM had a higher proportion of reburns and, unlike the 15 other ecoregions, fires did not occur at higher elevation or within greener areas. We suggest several possibilities for why, with the exception of AZ/NM, ecoregions did not experience a significant change in the proportion and pattern of refugia. In summary, while area burned has increased over the past four decades, there are substantial and consistent patterns of refugia that could support post-fire recovery dependent on their spatial patterns and ability to function as seeds sources for neighboring burned patches.


Introduction
Extreme weather, rising temperatures, severe drought, and increased ignitions (Westerling 2016) are creating conditions conducive to fires, leading to a global increase in burned area (Abatzoglou et al 2021) and forest loss (Tyukavina et al 2022).These expanding burned areas are not uniform in severity; fire interacts with landscapes to form a complex mosaic of burned and less affected areas.The islands of unburned or minimally burned areas within fire perimeters are commonly referred to as fire refugia (Krawchuk et al 2020).Fire refugia play an important role in enhancing post-fire ecosystem function and biodiversity, offering refuge to species during a fire event and serving as a vital seed source in the aftermath (Kolden et al 2017).As the world continues to warm, refugia could help promote forest persistence and disturbance recovery by maintaining legacy species from previous climate conditions (Dobrowski 2011, Stevens-Rumann et al 2018, Coop et al 2020).Fire behavior can shape the formation of refugia.If climatic or fuel conditions favor higher fire severity and a continuous supply of fuel, fires may burn more completely, resulting in fewer refugia (Whitman et al 2022).Conversely, under moderate fire weather the spatial heterogeneity of factors like fuel moisture, temperature, wind, topography may support the formation of refugia (Collins et al 2019).Fire suppression efforts like backburning may also play a role by reducing fire severity and generating more patchiness (Stephens et al 2009).
Studies have hypothesized that climate change and growing wildfire activity are leading to increases in fire severity and distance to seed sources (Harvey et al 2016) and a decrease in the proportion of refugia within fire perimeters (Kolden et al 2017, Collins et al 2019).There is some evidence of this.For example, studies in the United States have documented increased fire severity in specific regions (e.g.sites in the Northern Rockies, Parks et al 2018), specific vegetation types (e.g. the southern Rocky Mountain lower montane, Rocky Mountain subalpine, and California chaparral vegetation types, Picotte et al 2016), and within patch interiors (Cansler andMcKenzie 2014, Stevens et al 2017).In Alberta, Canada unburned areas within fire perimeters significantly have declined since 1985 (Whitman et al 2022).However, other studies have not found any consistent temporal trend in fire severity or proportion of fire refugia (e.g.Meddens et al 2018b, Chapman et al 2020).Thus, there is no scientific consensus regarding whether the proportion of refugia are changing across a wide swath of ecoregions (Meddens et al 2018b, Parks et al 2020, Buonanduci et al 2023).
One possible reason for the lack of scientific consensus is that refugia studies typically focus on narrow time frames and localized regions.While the characteristics of pixels and patches are important for refugia formation, so are event level factors like fire-driven weather, fuel moisture conditions, and fire suppression activities (such as backburning and slurry drops).Analyzing event-level patterns of fire refugia is important for understanding the spatial and temporal characteristics of fire regimes (Chuvieco et al 2016, Balch et al 2020) and provide valuable insights to guide land managers in planning for post-fire recovery and pre-fire mitigation.Yet, with a few notable exceptions (e.g.Kolden et al 2012, Meddens et al 2018b), few studies analyze refugia at the fire event level across ecoregions.
This study explores how fire refugia within forested areas have changed over time and across ecoregions between 1986 and 2021 across the Western US.We focus the study on fire event-level refugia metrics that reveal the area and patchiness (% refugia, total refugia area, median and maximum refugia patch size, number of refugia patches, refugia patch density, median and maximum distance to refugia).Our hypothesis is that in a preponderance of fires within western US ecoregions, % refugia has decreased over time while the distance to refugia patches within the fire perimeter has increased.This study has several distinctive characteristics that advance our understanding beyond previous work.First, it has the largest spatial and temporal extent of any refugia study (16 ecoregions in the Western US from 1986-2021).Secondly, it classifies fire refugia using a fire severity measure that is spatially consistent across ecoregions (bias corrected composite burn index (CBI), Parks et al 2019) and derived from image composites rather than single imagery.Third, this study focuses on fires that took place in relatively dense forest cover (>50% pre-fire forest), not in lower density savannas or woodlands where management is more prevalent (Chapman et al 2020).If refugia are becoming smaller and more fragmented, it may forecast a different successional pathway for forest vegetation as the planet warms (Blomdahl et al 2019).

Methods
To answer the research questions, we summarized refugia characteristics of 2968 fires that met criteria for inclusion in the study.We used the following workflow to answer the research questions (figure 1).The following sections follow each step of the workflow.

Step 1: Select wildfires within the 16 forested ecoregions in the Western US
This study focuses on wildfires that occurred in forested ecoregions in the Western US (figure 2).The dataset we used to identify wildfires is the Monitoring Trends in Burn Severity (MTBS) fire boundaries, which comprise all large wildfires (405+ ha in the Western US) starting in 1984 (Eidenshink et al 2007).MTBS fire boundaries are manually digitized based on standardized analysis of Landsat 4-9 imagery and incident perimeters where available.While the MTBS dataset omits smaller fires, it does not show a temporal bias (Picotte et al 2020, Iglesias et al 2022).
We selected fires included in the MTBS database (accessed November 2022) using the following criteria:  (1) Fires must be from the years 1986-2021 so that all fires have adequate before and after imagery.

Step 2: Mask out non-forest and edge pixels within each event
For each of the selected wildfires, we masked out pixels that were classified as non-forest one year before the fire (LCMS land cover product, Housman 2021).This step removes savannas and open woodlands from the analysis.We also masked out all pixels that were less than 100 m from the edge of the fire perimeter.

Step 4: Summarize event level refugia characteristics
After classifying refugia, we summarized the refugia characteristics within each fire in two ways.First we classified the composition of each fire (e.g. the area of each fire that was classified as forest refugia, forest non-refugia, non-forest vegetation, unvegetated, and no data).Secondly, we calculated fire refugia metrics for each fire event: • % refugia (proportion of pre-fire forest within fire boundary that is classified as refugia (i.e.CBI < 1.25) • total refugia area • median and maximum refugia patch size (patches are groupings of adjacent refugia pixels, connected by edge or diagonally).• number of refugia patches • refugia patch density • median and maximum distance to refugia (Euclidean distance from non-refugia pixel to refugia pixel, up to a maximum of 1.5 km).
For each of the metrics, we then used two estimators of the slope of change: OLS and Theil-Sen.The OLS estimator is based on the weighted mean of slopes between data pairs, while the Theil-Sen Slope is based on the median slope between data pairs (Sen 1968).
Thus OLS slope is sensitive to extreme values while Theil-Sen is robust, making it well suited for identifying trends in climate and weather data (Fernandes and Leblanc 2005).We calculate slope using both estimators since we are interested in the trends in refugia metrics both with and without extreme values.For % refugia and slope of % refugia, we also summarized the result within ecoregions.

Step 5: Explore results for the Arizona/New Mexico mountain ecoregion
We assessed two possible reasons for why AZ/NM experienced a significant increase in % refugia, while the other 15 ecoregions did not.A second possible explanation is that fires in AZ/NM could be increasingly occurring in forests with characteristics that are conducive to refugia formation (e.g.lower elevation, less green, or sparser forests).Lower elevation forests tend to have fire regimes with frequent, low severity fires associated with refugia (Covington and Moore 1994).To evaluate this, we calculated slope of change in the following within fire events: mean elevation, mean NDVI in pre-fire forest, and % of fire boundary that comprised forest prior to the fire (henceforth: % pre-fire forest).Note that % pre-fire forest ranges from 50% to 100%, since fires with <50% pre-fire forest cover were excluded from the analysis.

Number, size, and composition of fire events
The MTBS database containd 2968 fires that meet criteria for inclusion in the study.Over the full time series, the median fire size was 1398 ha, mean size was 5721 ha, and interquartile range was 3368 ha.The number of fires have increased over time (figure 4(a)) and so has the mean size of fires (figure 4(b)), though size continues to be highly variable.As the total area of MTBS perimeters increased, so did the area of both refugia and non-refugia (figure 4(c)).Over the full time series, the MTBS fire boundaries contained 27% forested refugia, 66% forested non-refugia, and 7% other (non-forest, non-vegetated, or no data).These percentages remained fairly stable over time (figure 4(d)).In summary, both fire perimeters and the refugia within increased over time, but there has been no obvious trend in % refugia.

Slope of change of fire refugia metrics
The slope of % refugia is significantly positive for both OLS and Theil-Sen estimators, indicating that % refugia has increased during the time period (table 1).
However, this change was driven entirely by the increase in % refugia in the AZ/NM.When excluding AZ/NM, the slope of % refugia was not significantly different from zero.For some refugia metrics (total refugia area, maximum refugia patch size, number of refugia patches, and maximum distance to refugia), results were equivocal.The OLS slope of change was positive and significant but the Theil-Sen slope was not significantly different from zero.The interpretation is that the median slope shows 'no change' while the mean slope shows 'positive increase' driven by fires with unusually high values.For the remaining refugia metrics (refugia patch size, refugia patch density, or mean distance to refugia), slope of change was not significantly different from zero regardless of the estimator (table 1).

Ecoregion-level results
Within the forested areas of each fire boundary, % refugia varied considerably by ecoregion.Fires in the Rockies, Cascades, and Sierra Nevada ecoregions were characterized by median refugia of 30% or less, while fires in AZ/NM and the Coast Range ecoregions had much higher median refugia (figure 5(a)).While there were big differences in % refugia between ecoregions, there was no significant change over time in most ecoregions (figure 5(b)).The Coast Range ecoregion had a decrease (negative slope) in % refugia, and the Cascades and Blue Mountain ecoregions had an increase (positive slope) in % refugia.However, these changes were not significant with the number of fires in the dataset.In AZ/NM, however, the change in % refugia was positive and significant (figure 5(b)).In the following section, we explore reasons why the trajectory of refugia was different in AZ/NM.

Exploring the increase in refugia in the AZ/NM ecoregion
To help understand why AZ/NM experienced a significant increase in % refugia, we evaluated two possible explanations.The first explanation is that re-burns and other forms of forest loss (e.g.thinning and harvesting) could have increased in AZ/NM more than elsewhere, thus promoting refugia.The data offers some support to the idea that AZ/NM is experiencing a greater increase in reburning than the other 15 ecoregions.For both AZ/NM and the 15 other ecoregions, % reburn increased significantly with the OLS estimator (mean slope) but not the Theil-Sen estimator (median slope) (table 2).The AZ/NM ecoregion had a much higher OLS slope than the 15 other ecoregions, suggesting that the increase in % reburn was higher there, driven by extreme values (table 2).The change in % non-fire forest loss was not significant with either estimator.
The second possible explanation is that fires in AZ/NM could be increasingly occurring in forests with characteristics that are conducive to refugia  formation (e.g.lower elevation, less green, or sparser forests).For fires in AZ/NM, elevation and pre-fire NDVI did not change significantly (table 2).In contrast, for fires in the other 15 ecoregions elevation and pre-fire NDVI increased significantly.

Discussion
We expected to observe a decrease in % refugia in a preponderance of Western US forested ecoregions, but instead found a significant increase in 1  An notable exception is that there was an increase in the maximum distance to refugia, driven by extreme values (i.e.significant with OLS slope, not Theil-Sen slope).This suggests that very large fires may pose a risk for recovery.
While our results show no evidence of widespread changes in % refugia across the West, AZ/NM was an exception.That AZ/NM had a high % refugia is not surprising, as the ecoregion's fire regime is historically characterized by frequent low severity fires that produce little mortality in the structurally dominant forest vegetation (Covington and Moore 1994).Fires in AZ/NM experienced a greater increase in reburning than elsewhere, and did not experience a change in non-fire forest loss (i.e.caused by thinning, harvesting, windthrow, or other non-fire disturbances).This raises the possibility that the increase in fire refugia could be related to increased reburning and consequent reduction in fire severity.Unlike the other 15 ecoregions, fires in AZ/NM did not experience a shift to higher elevations or a decrease in pre-fire NDVI.A shift to higher or greener forests, had it occurred in AZ/NM, might have muted the increase in % refugia.Additional research would be needed to determine whether there is a causal link between these factors.
One limitation to the exploration of the AZ/NM ecoregion investigation is that we do not have spatially explicit data on forest treatments during the 1980s-1990s, so we cannot distinguish natural disturbances from those that were part of a management plan.Another limitation to the exploration of AZ/NM is that we did not investigate changes within the forest land cover type.For example, Gambel oak (Quercus gambelii) resprouting and replacing conifer forest (Coop et al 2016) could look like refugia on imagery.Finally, in both the main analysis and AZ/NM exploration, we did not control for the increase in familywise error rate in statistical tests.Due to the relatively large sample size and modest number of tests, we believe this to be a reasonable decision but would encourage replication by other studies.
The study shows that dense forest (pixels >50% cover) show a consistent survivorship in wildfires, with the exception of very large fires where 'maximum distance to refugia' has increased.Overall, the results of this study imply that there are substantial refugia areas that could be used in planning for recovery or restoration efforts.Because post-fire climate conditions may restrict regeneration (Davis et al 2019), it is important to place remaining refugia in areas with better microclimates and lower fire severity, facilitating enhanced regeneration (Davis et al 2023).Future efforts should explore the persistence of refugia, whether finer-scale resolution data (e.g., at the tree level) reveal different trends in refugia, and whether changes within the forest land cover class can explain some of the patterns.For example, Gambel oak resprouting and replacing conifer forest could be associated with the observed changes in refugia in AZ/NM, but this hypothesis remains to be tested at the resolution of species.As the trends in area burned have been concerning, the prevalence of refugia offers some hope that ecosystems continue to have some of the raw material to begin the process of recovery.

Conclusion
We did not find evidence that supports the proposition that refugia are becoming rarer in Western US forests.The total area of refugia and number of refugia have increased proportionally to the total fire boundary area.While % refugia has not changed significantly overall, it has changed in specific ecoregions most notably a significant increase in AZ/NM.Future studies, over a longer time period, may reveal stronger trends.There is also a need to compare fire severity and % refugia in 'newly expanded burns' (e.g.shoulder seasons, or at night) with 'traditional burns' , which may display diverging trends in refugia.This study provides compelling evidence that refugia continue to make a considerable contribution to post-fire landscapes, even in the face of an overall increase in burned area across the Western US.In other words, refugia are not disappearing despite the growing extent of wildfires.Refugia studies like this, which characterize broad event-level trends, are important for identifying wide trends, understanding the ecological contexts and landscape mosaics in which refugia exist.

References
For example, previous studies in the United States have focused on the Colorado Front Range (Chapman et al 2020), the central Sierra Nevada (Blomdahl et al 2019), the Pacific Northwest (Krawchuk et al 2016, Meddens et al 2018b, Meigs et al 2020) and elsewhere.Furthermore, most refugia studies use individual pixels or patches as units of analysis (e.g.Krawchuk et al 2016, Meddens et al 2016, Meigs and Krawchuk 2018, Blomdahl et al 2019, Collins et al 2019, Chapman et al 2020, Downing et al 2021, Meigs et al 2020, Mackey et al 2021, Talucci et al 2022).
(2) The fire incident type must be MTBS type 'wildfire' (excluding the MTBS types 'prescribed' , 'wildland fire use' , and 'unknown').(3) Fire boundaries must overlap forested level 3 ecoregions in the Western US (Omernick and Griffith 2014).(4) Fire boundaries must have 50%+ of the area classified as forest cover in the year before the fire (Harvey et al 2016, Meigs et al 2018).To identify forest cover we used the USGS landscape change monitoring system (LCMS) data (Healey et al 2018), which classifies pixels as trees if most of the pixel's area comprises trees and no other land cover covers more than 10%.

Figure 3 .
Figure 3. Example refugia classification: the Waldo Canyon Fire near Colorado Springs, CO, USA in 2012.

Figure 4 .
Figure 4. Change over time of (a) number of fires, (b) fire size, (c) area of fire perimeters in 5 classes, and (d) percent of perimeters in 5 classes.
M A, Haire S L, Coop J, Parisien M A, Whitman E, Chong G and Miller C 2016 Topographic and fire weather controls of fire refugia in forested ecosystems of northwestern North America Ecosphere 7 e01632 Krawchuk M A, Meigs G W, Cartwright J M, Coop J D, Davis R, Holz A and Meddens A J 2020 Disturbance refugia within mosaics of forest fire, drought, and insect outbreaks Front.Ecol.Environ.18 235-44 Mackey B, Lindenmayer D, Norman P, Taylor C and Gould S 2021 Are fire refugia less predictable due to climate change?Environ.Res.Lett.16 114028 Meddens A J, Kolden C A and Lutz J A 2016 Detecting unburned areas within wildfire perimeters using Landsat and ancillary data across the northwestern United States Remote Sens. Environ.186 275-85 Meddens A J, Kolden C A, Lutz J A, Abatzoglou J T and Hudak A T 2018b Spatiotemporal patterns of unburned areas within fire perimeters in the northwestern United States from 1984 to 2014 Ecosphere 9 e02029 Meddens A, Kolden C, Lutz J, Smith A, Cansler C, Abatzoglou J, Meigs G, Downing W and Krawchuk M 2018a Fire refugia: what are they and why do they matter for global change?Bioscience 68 944-54 Meigs G W, Dunn C J, Parks S A and Krawchuk M A 2020 Influence of topography and fuels on fire refugia probability under varying fire weather conditions in forests of the Pacific Northwest, USA Can.J. For.Res.50 636-47 Meigs G W and Krawchuk M A 2018 Composition and structure of forest fire refugia: what are the ecosystem legacies across burned landscapes?Forests 9 243 Miller J and Thode A 2007 Quantifying burn severity in a heterogeneous landscape with a relative version of the delta normalized burn ratio (dNBR) Remote Sens. Environ.109 66-80 Nagy R C, Fusco E, Bradley B, Abatzoglou J T and Balch J 2018 Human-related ignitions increase the number of large wildfires across US ecoregions Fire 1 4 Omernik J M and Griffith G E 2014 Ecoregions of the conterminous United States: evolution of a hierarchical spatial framework Environ.Manage.54 1249-66 Parks S A and Abatzoglou J T 2020 Warmer and drier fire seasons contribute to increases in area burned at high severity in western US forests from 1985 to 2017 Geophys.Res.Lett.47 e2020GL089858 Parks S A, Holsinger L M, Koontz M J, Collins L, Whitman E, Parisien M A, Loehman R A, Barnes J L and Bourdon J F 2019 Giving ecological meaning to satellite-derived fire

2.3. Step 3: Classify refugia pixels based on CBI within each event
This excluded the approximate areas where conifer seeds could readily spread from outside the fire boundary (Steel et al 2018) and reduced edge effects and digitizing errors (Parks et al 2018, Stevens-Rumann et al 2018).Finally, pixels were masked out if there were no valid pre-or post-fire pixels (rare, but evident in a few of the earliest fires).
the Parks et al (2019) method, which uses random forests regression to predict CBI based primarily on One possible explanation is that re-burns and other forms of forest loss (e.g.thinning and harvesting) could have increased in AZ/NM, thus promoting refugia.Forests with lower stand density that burn frequently tend to experience lower severity fires(Covington et al 1994), and are associated with higher % refugia (Meddens et al

Table 1 .
Slope of change in event-level fire refugia characteristics 1986-2021 (full results in supplementary materials).

Table 2 .
Mean To evaluate characteristics within the last 10 years of fire, date range is 1995-2021.that can act as seed sources.State changes (e.g.forest to grassland) may be more likely caused by climate-driven recruitment failure (Davis et al 2023) than seed availability from lack of refugia.
environmental characteristics of fire perimeters in AZ/NM versus 15 other ecoregions.(fullresults in supplementary materials).AZ/NM ecoregion 15 other ecoregions OLS slope (std err) Theil-Sen Slope (std err) OLS slope (std err) Theil-Sen slope (std err)a Significant at p < .05. b Significant at p < .001.c % of the fire boundary that had experienced a fire in the previous 10 years.Date range is 1995-2021.d % of the fire boundary that had experienced forest loss caused by a disturbance other than fire in the previous 10 years.Non-fire forest loss definition: pixels from the LCMS that were classified as 'forest' 10 years prior to the fire and then classified as 'fast loss' in subsequent years up to the year before fire.'Fast loss' can be caused by harvest, mechanical, debris (e.g.landslides), and fire.Pixels within previous MTBS fire boundaries were excluded, so the forest loss excludes all but the smallest fires (Housman et al 2021).ated patches of surviving trees are disproportionately important for tree regeneration (Coop et al 2019) so missing such refugia could potentially be important.Regardless of the reason, the assumption that refugia are becoming rarer is not borne out by this analysis.Fires continue to create heterogeneous landscapes that include refugia Balch J K and Travis W R 2022 US fires became larger, more frequent, and more widespread in the 2000s Sci.Adv. 8 eabc0020 Kolden C, Bleeker T, Smith A, Poulos H and Camp A 2017 Fire effects on historical wildfire refugia in contemporary wildfires Forests 8 f8100400 Kolden C, Lutz J A, Key C H, Kane J T and Van Wagtendonk J W 2012 Mapped versus actual burned area within wildfire perimeters: characterizing the unburned For.Ecol.Manage.286 38-47 Krawchuk