Climatic controls of fire activity in the red pine forests of eastern North America

Large-scale modes of climate variability influence forest fire activity and may modulate the future patterns of natural disturbances. We studied the effects of long-term changes in climate upon the fire regime in the red pine forests of eastern North America using (a) a network of sites with dendrochronological reconstructions of fire histories over 1700 – 1900 A.D., (b) reconstructed chronologies of climate indices (1700 – 1900), and (c) 20th century observational records of climate indices, local surface climate and fire (1950s-2021). We hypothesized that (H1) there are states of atmospheric circulation that are consistently associated with increased fire activity, (H2) these states mark periods of increased climatological fire hazard, and (H3) the observed decline in fire activity in the 20th century is associated with a long-term decline in the frequency of fire-prone states. At the annual scale, years with significantly higher fire activity in the reconstructed and modern fire records were consistently associated with the positive phases of the Pacific North American pattern (PNA), either independently or in combination with the positive phase of the El Ni ˜ no-Southern Oscillation index (ENSO). During years with both ENSO and PNA in their positive state, the region experienced positive mid-tropospheric heights and temperature anomalies resulting in drought conditions. The fire-prone climate states identified in the reconstructed records became less frequent in 1850 but re-emerged in the 20th century. While our study did not demonstrate a direct influence of climate on the observed decrease in fire activity in the 20th century, it does reveal a clear climate signal embedded within the fire history reconstruction of the region over the past centuries. This study underscores the importance of considering large-scale climatic patterns in understanding historical fire regimes and highlights their role for future fire dynamics in the region and shaping ecological effects of future fires.


Introduction
Fire is a fundamental Earth system process influencing most terrestrial ecosystems (He and Lamont 2018;McLauchlan et al., 2020), as well as global processes (e.g., surface albedo and the global carbon balance) (Liu et al. 2019).Changes in fire disturbance severity and/or frequency have a strong impact on the Earth's surface radiative budget and climate, especially at high latitudes (Liu et al. 2019).A positive trend in wildfire activity in high latitudes has been particularly pronounced, with increases in greenhouse gases and atmosphere aerosols providing positive feedback that further accelerates climate change (Zhao et al. 2021).The 2023 fire season in Canada exemplified this impact, burning a record-breaking 18 401 197 hectares and causing thousands to be evacuated, major cities to suffer poor air quality, and necessitating The blue lines delimit the study period.
extensive firefighting efforts with assistance from various countries (Natural Resources Canada 2023).
Fire has been an essential natural disturbance in the mixed-pine forests of the boreal-temperate interface in eastern North America (Drever et al. 2006;Nyamai et al. 2014;Meunier et al. 2019;Frelich et al. 2021).Historically, frequent low-to moderate-severity surface fires maintained the dominance of fire-dependent conifers (i.e., red pine [P.resinosa Ait.] and white pine [Pinus strobus L.]) in these forests (Guyette and Dey 1995;Drobyshev et al. 2008a;Stambaugh et al. 2018).This fire regime created the necessary conditions for pines' regeneration by removing or reducing the organic layer of the soil (Nyamai et al. 2014;Stambaugh et al. 2019), augmenting light conditions at the forest floor by reducing canopy cover (McRae et al. 1994), and eliminating shade-tolerant competitors (Nyamai et al. 2014;Stambaugh et al. 2019).In this way, frequent low-severity fires with occasional mixed-severity fires, led to complex, multi-cohort stands (Drobyshev et al. 2008b;Fraver and Palik 2012;Meunier et al. 2019).
Studies of the historical fire regimes of red pine forests have suggested a significant role of Native American land-use practices in shaping the frequent-surface fire regime.This assertion has gained support from multiple lines of evidence, including the documented widespread use of fire by Native Americans and comparisons with the modern rate of lightning-caused fires, which alone fail to explain the historical high fire frequency (Loope and Anderton 1998;Johnson and Kipfmueller 2016).Nonetheless, even with dense pre-European population in some regions in eastern North America, there were regions with forests that experienced a limited anthropogenic impact (Oswald et al. 2020).Indigenous fire management practices were often localized (Roos 2020), indicating that at the regional scale, climate remained a significant determinant of the historical fire regime (Bergeron 1991;Drobyshev et al. 2012;Oswald et al. 2020).
The 1930s witnessed a sharp decline in area burned in the eastern North American forests, including the region corresponding to the natural distribution range of red pine (Drobyshev et al. 2008a) (Fig. 1).These changes in fire activity had a profound negative impact on the regeneration of fire-adapted species, such as red pine and mixed pine forests, at the regional scale.Fire exclusion (Nowacki and Abrams 2008) and logging (Friedman and Reich 2005) gave competitive advantage to fire-intolerant and mesophytic hardwoods (e.g., maples [Acer spp.] and birch [Betula spp.]).The higher foliar moisture content of these deciduous trees and, therefore, the lower flammability of these fuels as compared to coniferous fuels, promoted moist and generally less fire-prone conditions (Nowacki and Abrams 2008).The decreased flammability in forests historically dominated by pines (often in combination with oaks) contributed to the conversion of mixed pine forests to northern hardwoods stands with abundant aspen (Populus spp.), birch, maples, and oaks (Quercus spp.).This pattern has been reported throughout eastern North America over most of the 20th century (Friedman and Reich 2005;Schulte et al. 2007).In contrast, in the northern fringes of the red pine distribution, where conifers (Abies and Picea) are present in the subcanopy, the lack of fire increases the amounts of easily burnable ladder fuels.Such fuels favor crown and lethal fires, decreasing the relative proportion of surface fires in the total pool of fires.Dominance of crown fires in the boreal forest has been reported as the factor explaining the northern limit of red pine (Flannigan and Bergeron 1998).These changes have led to the disappearance of these forests, with red and white pine forests now covering only 0.6 % of their pre-settlement range (Frelich 1995).The well-documented effect of European colonization on fire regimes and the subsequent changes in the composition of mixed pine forests led to broad acceptance of human land-use as the main driver of observed changes (e.g.Guyette et al. 2016;Kipfmueller et al. 2017).However, at the northern limit of the distribution range of red and white pines, where human population density is low and colonization began later, climate variability has been suggested as the main driver of the decline in fire activity in the mid-1800s (Bergeron and Archambault 1993).
Atmospheric circulation anomalies can create fire-conducive conditions and synchronize fire activity at the regional and sub-continental scales (Fauria and Johnson 2006;Fauria and Johnson 2008).This observation warrants a study of broad-scale climate forcing upon regional fire regimes, including the red pine distribution range.Anomalies in atmospheric circulation and local winter climate of the red pine distribution area (RPDA) have been linked to various climate oscillation indices, including El Niño-Southern Oscillation (ENSO), Pacific North American pattern (PNA), and Pacific Decadal Oscillation (PDO) (Rodionov and Assel 2000;Rodionov and Assel 2003;Bai and Wang 2012;Wang et al. 2018).
The PNA pattern plays a crucial role in shaping the region's local climate by directly affecting the strength and position of tropospheric anomalies over North America, the Polar Front Jet Stream (PFJ), and ultimately the regional climate.The PNA manifests as a distinct Rossby wave train with four action centers extending from the North Pacific across the North American continent (Wallace and Gutzler 1981).The most prominent action center of the PNA is a low-pressure system (cyclone) in the North Pacific around the Aleutian Islands known as the Aleutian Low (AL).South of the AL is a high-pressure system (anticyclone) around Hawaii.Downstream, the PNA features a mid-tropospheric ridge-trough system over North America, characterized by an above-normal geopotential height anomaly over western Canada and a below-normal height over the southeastern United States (Wallace and Gutzler 1981).This configuration strongly correlates with positive temperature anomalies over western Canada and negative ones over southeastern U.S. (Wallace and Gutzler 1981), ultimately linking these dynamics to forest fuel conditions at the regional scale (Fauria and Johnson 2008).The ridge-trough pattern significantly influences the position of the jet stream and associated storm tracks over North America.The PNA results in a meridional flow of the jet stream, which brings fewer Pacific-derived winter storm tracks and dry anomalies, particularly in the West coast of North America (Rodysill et al. 2018).However, variability in the intensity and position of the AL can alter the ridge-trough pattern, giving rise to different local climate patterns (Rodionov and Assel 2003;Wang et al. 2012;Lin et al. 2023).During a negative PNA phase, a more zonal and northerly jet stream over eastern North America allows for increased moisture from the Gulf of Mexico to penetrate in this region (Rodysill et al. 2018).
ENSO is another important factor driving regional climate in eastern North America as it plays a key role in triggering the PNA pattern and shaping drought conditions over the continent.Warm SST anomalies in the Tropical Pacific occur during El Niño phase of ENSO and cooling anomalies during La Niña events.El Niño events are generally associated with the deepening of the AL and the PNA pattern (Taschetto et al. 2020).However, the influence of ENSO on the AL and its broader extratropical effects varies significantly among individual events.This variability is attributed to several factors collectively termed "ENSO diversity": the asymmetry between El Niño and La Niña events, the magnitude of ENSO events, and the spatial distribution of the ENSO-related SST anomalies (Capotondi et al. 2015).Previous studies have highlighted the asymmetrical response of the winter climate in the Great Lakes region to ENSO, with strong El Niño events associated with milder winters, while the link between cold winters and La Niña events are notably weaker (Bai et al. 2012).The response of the winter climate in the region to ENSO is also nonlinear, where strong El Niño events result in milder winters, yet the effect of moderate El Niño events are less predictable (Rodionov and Assel 2003).
The PDO signal in the eastern North American winter climate is realized through its association with the AL and the PNA.The PDO is the dominant low-frequency climate variability in the North Pacific and is significantly influenced by anomalies in the AL (Mantua et al. 1997).The PDO`s positive phase is characterized by warm SST anomalies along the west coast of the U.S. and cold SST anomalies from the western North Pacific, and the negative phase exhibits the opposite pattern (Mantua et al. 1997).The positive phase of PDO is associated with a deeper AL, whereas the negative phase is associated with a shallow AL (Newman et al. 2016;Larson et al. 2022).Although ENSO plays an important influence on the PDO (Zhao et al. 2021;Larson et al. 2022), the PDO results from complex interactions beyond that of tropical variability (Newman et al. 2016).The "atmospheric bridge" connecting tropical Pacific to the North Pacific is bi-directional and the PDO can also influence the tropical Pacific (Alexander et al. 2002;Newman et al. 2016).Studies have reported constructive and destructive interference between the ENSO and PDO according to their relative phases, with stronger and more coherent ENSO effects when the ENSO and PDO are in phase (e.g., El Niño and positive PDO) (Gershunov and Barnett 1998;Hu and Huang 2009).
While circulation patterns represented by these indices typically exert their strongest effects on the local climate during winter, their impacts may extend beyond the cold season.For example, winter ENSO and PDO anomalies have been associated with variations in summer moisture availability in Canada (Shabbar 2006).Mild winters can result in a thinner snowpack that melts earlier, accelerating fuel drying and widening the period for effective ignitions (Westerling et al. 2006).Consequently, the connection between winter climate indices and subsequent fire activity has been widely reported across North America (Heyerdahl et al. 2002;Duffy et al. 2005;Brown 2006;Kitzberger et al. 2007).Previous research has demonstrated links between climate indices and fire activity in eastern North American boreal forests, located north of the red pine distribution range (Girardin et al. 2004;Girardin et al. 2006a;Goff et al. 2007).Within the RPDA itself, historical associations between fire activity and local climate conditions have been documented (Drobyshev et al. 2012), despite abundant human fire ignitions likely masking the climate signal during the colonization waves (Stambaugh et al. 2018).However, despite these insights, the association between historical fire activity across RPDA and modes of climate variability, such as ENSO, PNA, and PDO remains poorly documented.
The current study addresses the gap in understanding the influence of large-scale climate variability, such as ENSO, PNA, and PDO on fire regimes in eastern North American mixed-pine forests with red pine.We compiled a network of sites with dendrochronologically reconstructed fire histories and focused specifically on the High Fire Activity Years (HFAYs), which are more likely to be climate-driven (Farris et al. 2010;Falk et al. 2011;Drobyshev et al. 2014).By broadening our analyses to a regional scale, we aimed to capture climate signals influencing the variability of fire occurrence, which might have been obscured by unique characteristics of single sites and the possible impact of human-related ignitions (Stambaugh et al. 2018).
Geographically, the network covered the present range of red pine and stretched over a 200-year period, 1700-1900 A.D. The period partially covered the Little Ice Age (LIA, ~1300-1800s) and the subsequent warming period, and ended before the onset of the industrial era and active fire suppression (1920s and 1930s).The selected period featured abundant fires and lack of fire suppression policies that could affect climate signal in fire chronologies.We focused on red pine forests, which we defined as forests with red pine in the overstory, since these stands are surface fire-dominated.High fire frequencies facilitated the detection of climate-driven changes in fire regimes, as compared to forests that burn infrequently and are often dominated by deciduous trees or coniferous shade-tolerant and fire-intolerant species (e.g., spruce and fir, Picea and Abies, respectively).
We explored the association between fire activity and the climate indices ENSO, PNA, and PDO.We tested three hypotheses: (H1) there are states of atmospheric circulation captured by climate indices that are consistently associated with increased fire activity over the distribution of red pine forests in eastern North America; (H2) these states mark periods of increased climatological fire hazard over the study area, as suggested by instrumental observations of the local surface climate; and (H3) the observed decline in fire activity in the 20th century is due to a long-term decline in the frequency of fire-prone states.We tested H1 on both dendrochronological and modern (1959-2021) fire datasets using a contingency analysis framework.We tested H2 by contrasting the fireprone states we identified in the reconstructed records with observational high-resolution climate over the modern record .In doing so, we borrowed heavily from the climate analogue technique applied in a study of wildfire synchrony over western North America (Kitzberger et al. 2007).We tested H3 by conducting a regime shift analysis for each identified fire-prone climate state over the 1675-2021 period.Through testing these hypotheses, we re-approached the discussion on the contribution of the climate in addition to human factors, which have already been established at more localized scales, in shaping fire regimes in eastern North America.In addition to building a knowledge base on climate controls of disturbance regimes in mixed pine forests, our results contribute to quantifying the range of natural variability in the levels of fire activity and establishing reference points for conservation-oriented management of these ecosystems.The study helps quantify large-scale climate controls of historical, current, and future dynamics of fire activity in an ecoregion that faces a high risk of negative climate change impacts (Montour et al. 2020).
A humid-continental climate with cold winters and warm to hot summers is prevalent over the region (Peel et al. 2007).The study area encompasses a considerable gradient in average annual precipitation, average minimum temperature for January, and average maximum temperature for July: 889.8 mm, − 24.6 • C, and 23.5 • C at the northernmost site, Duparquet, QC, for 1971-2000(Environment Canada 2011); 866.9 mm, − 5.8 • C, and 28.8 • C at the southernmost site, Pike Knob Preserve, WV, for 1906-2008(Young et al. 2018); 520 mm, − 15.6 • C, and 25.7 • C at the easternmost site, Fourth Machias Lake, ME, for 1948-1998(Young et al. 2018); and 644.4 mm, − 20.9 • C, and 27 • C at the westernmost site, Itasca State Park, MN, for 1911-2021(Young et al. 2018).
The study area roughly corresponds to the distribution of red pine, which extends eastward from Manitoba and Minnesota to Newfoundland, and southward from Central Quebec to West Virginia (Roberts and Mallik 1994).Our network of sites had high replication in the area where red pine attains its maximum abundance, i.e., along a northwest-southeast line close to the southern limit of its range, the Great Lakes-St.Lawrence-Acadian forest (Flannigan and Woodward 1994).The network contained 28 sites within the Great Lakes-St.Lawrence forest region (Rowe 1972), and covered an area of approximately 2.73•10 6 km 2 .Sites were located in seven U.S. states (Maine, Michigan, Minnesota, New York, Pennsylvania, Virginia, and Wisconsin) and two Canadian provinces (Quebec and Ontario).
To select the sites for this study, we focused on fire histories developed in forests with a red pine component.Our site selection criterium was solely based on the availability of fire history data in the forests where red pine was present.A review of the metadata from our sites and the other published sites indicated that red pine was an important but not always the dominant component in the forest canopies.Red pine is an effective recorder of fire history because (1) it is adapted to survive low-(surface fires) to moderate-severity fires (mixture of surface and stand replacing fires); and ( 2) it has a relatively slow decomposition rate, especially in xeric sites (Johnson and Kipfmueller 2016).Sites were, at the time of sampling, comprised of (a) pure red pine stands; (b) conifer stands or mixed-wood stands where red pine was a co-dominant species in the canopy; or (c) hardwood-dominated stands (e.g., dominated by red oak, Quercus rubra) that previously consisted of red pine, as suggested by dendrochronological reconstructions.

Fire history database
We compiled a network of 28 sites (Table Supplementary Information SI 1) comprised primarily of existing datasets via the authors or published literature.We also analyzed limited unpublished datasets collected and dated by the authors covering southwestern Quebec and Ontario.Our database contained coordinates, fire years by site, and the period covered by site reconstructions.
We ensured that all site chronologies represented spatially independent observations by establishing a minimum distance of 40 km between neighboring sites.The rationale for this threshold was based on the spatial reconstruction of fires at Seney National Wildlife Refuge (SNWR), MI, covering the 1700-2007 period (Drobyshev et al. 2008).The reconstruction revealed several years during the 1700s and 1800s, when fires covered most of the SNWR landscape (386 km 2 ).A circle covering this area would have a diameter of approximately 20 km.We doubled this value to arrive at the selected threshold of 40 km.We combined sites in proximity to each other, within 40 km radius, so that they became one site with a composite fire chronology.The implementation of this threshold resulted in 25 % of all sites (n = 7) being combined in this manner.By adopting 40 km threshold, we addressed spatial autocorrelation, which could potentially inflate estimations of High Fire Activity Years (HFAYs) frequencies.
An exception to the 40 km threshold was made for the Seney site within the SNWR given that sampling at Seney was more exhaustive as it was designed to reconstruct annual area burned, and the number of fires recorded there was higher.Thus, this site was not pooled with nearby sites in the Hiawatha National forest despite their close proximity.The SNWR study (Drobyshev et al. 2008) stands as the sole example of a spatially explicit reconstruction of fire activity in a landscape characterized by (a) abundant red pine, particularly as the dominating mixed pine forest species; and (b) a lack of permanent waterbodies that would make "perfect" firebreaks, in contrast to the conditions described in Kipfmuller et al.'s (2017) study.
We acknowledge that any solution to group original sites involves a trade-off between keeping the number of areas sufficiently high to allow for a more statistically stable identification of HFAYs and removing autocorrelation.Such autocorrelation may, in fact, originate from both regional climate forcing and annual synchrony in the implementation of land use policies within parts of the RPDA.We were also aware of the growing subjectivity in data treatment that arises when introducing sitespecific clustering algorithms.We thoroughly explored various options of site clustering in the early stages of our analyses, seeking to establish a minimal number of straightforward rules that could accommodate the considerations presented above.
Our statistical framework was based entirely on the analysis of fire occurrence data because the spatial estimates (i.e., fire sizes) were not available for most of the sites.We selected the historical (reconstructed) period based on two criteria: (1) the time period covered ≥ 95 % of all sites, and (2) it did not include periods with declining fire activity, whether due to climatic reasons (Drobyshev et al. 2017) or due to fire suppression.The inclusion of the second criterium was due to our concern that a declining number of fire events would compromise our ability to quantify climate-fire links in the dataset and likely introduce non-climatic trends in the chronology.To objectively define the period of decreased fire activity, we conducted regime shift analyses for each site on smoothed decadal fire occurrence using Rodionov's sequential t-test algorithm (Rodionov 2004).We explored different choices for the L parameter, or moving timeframe, including both 5 and 10 years.We chose an L of 5, which allowed us to make algorithm more sensitive to the changes in fire regime.The Hubert weight parameter was set to 1, and α was set at 0.05 to ensure statistical robustness setting the L parameter, or moving timeframe, to five years; the Hubert weight parameter to 1, and α = 0.05.
We analyzed the reconstructed fire activity during 1700-1900.The start year, 1700, marked the onset of the regime shift towards a higher fire activity at the oldest part of the reconstructions (Fig. SI 2).By selecting the start of the record in this way, we address the "fading record problem" that refers to the decline in replication at the oldest part of the reconstruction, thereby compromising its reliability (Swetnam et al. 1999).The end year, 1900, marked the onset of a lower fire activity period that broadly coincided with the introduction of fire suppression policies.
We opted not to apply a minimum threshold for the number of scarred trees required to identify a fire year.This decision was based on several considerations.First, including years where fires were recorded by only one tree allowed us to capture the full extent of fire occurrence.The focus of our study was on the annual synchrony of the fire years across the distribution range of red pine, rather than the evolution of fire frequencies at individual sites.This made the size of individual fires at the site level of little relevance.Second, the design of the contingency analysis accounted for eventual effect of "fire replication" at a site level and amount of "recording sites" (see the following sub-section).Third, we ensured the robustness of our analyses by using only climate-driven HFAYs for later analyses.Finally, applying a consistent threshold across all sites was not feasible as we did not have access to raw data from all sites.

Identifying climate-driven fires and assessing their association with local climate
We used synchrony in annual fire occurrence across the network of sites as a proxy of climate forcing upon regional fire activity (Farris et al. 2010;Falk et al. 2011;Drobyshev et al. 2014).We identified HFAYs as those years when the synchrony of burned sites was higher than expected under the assumption of a random distribution of fire years across sites and years (Drobyshev et al. 2015).To test for non-randomness in the synchrony level of burned sites, we used contingency analysis along with non-parametric bootstrapping (1 000 runs with replacement).Specifically, we compared the observed synchrony levels of burned sites to a null distribution generated through non-parametric bootstrapping.To define HFAYs, we assumed a binominal distribution of the fire occurrence across sites and calculated expected frequencies using the following formula: where N is the total number of recording sites in the analysis of a specific period; X is the number of burned sites in a single year; p is the probability of a site burning in any year, and q is the inverse of this probability.Calculating p(X) in this way accounted for variability in both site replication and fire frequencies at each site, i.e., number of fires recorded in a site.To further account for the changing trends in the number of burned sites over the studied period, revealed by a fitted segmented linear regression on the number of burned sites versus time (Fig. SI 1), we partitioned our reconstructed fire data into four subperiods (1700-1749, 1750-1799, 1800-1849, 1850-1900) and calculated p(X) for each using the period-specific "process density", which was the number of fire years at each site.
We then examined if the reconstructed HFAYs occurred concomitantly with fire-prone local climate.Specifically, we tested if summer drought, temperature, and mid-troposphere pressure at 500 hPa geopotential height (Z500) differed significantly from their means during HFAYs for the reconstructed record.We conducted spatial composite analyses and produced maps for the region 35 • to 50 • N and − 95 • to − 65 • E. The temporal resolution and coverage differed for the different gridded climate fields because we used different sources.We examined drought during historical HFAYs at the annual scale for the whole reconstructed record .On the other hand, for temperature, and 500 hPa geopotential height, we conducted the analyses at a higher temporal resolution, i.e., seasonal, looking at two seasonal periods: (1) winter-early spring (January-April) and ( 2) late spring-summer (May-August).These seasonal analyses, however, covered only part of the reconstructed record .For the drought field, we used the North American Drought Atlas (NADA), a millennial gridded (2.5 • x 2.5 • ) annual tree-ring based reconstruction of the summer (June-August) Palmer Drought Severity Index (PDSI) (Cook et al. 2004).For the temperature and mid-tropospheric pressure fields, we used the Twentieth Century Reanalysis version 3 (20CRv3) datasets, which comprises grided (1.0 • x 1.0 • ) three-hourly estimates of meteorological variations extending back to 1836 (Slivinski et al. 2019).We used the online tool ClimateExplorer (Trouet and Van Oldenborgh 2013) to conduct these analyses.

Climate indices
We examined three climate indices: ENSO, PNA, and PDO as previous research has demonstrated their influence on local climate conditions within our study area (Rodionov and Assel 2000;Rodionov and Assel 2003;Bai and Wang 2012).For analyses over the historical record (1700-1900), we relied on reconstructions of these indices, while for analyses spanning the modern record (1950s-2021), we used observational data.
In the case of ENSO, we focused on the canonical Eastern Pacific (EP) ENSO because the EP El Niño has a more significant impact on winter climate conditions over the RPDA compared to the Central Pacific (CP) El Niño (Yu et al. 2012).This effect arises from the EP El Niño's eastward displacement of the AL, which brings the PNA's ridge over the North American continent closer to the RPDA (Beniche et al. 2024).We used Freund et al.'s (2019) coral-based reconstruction of the winter (December-February) EP ENSO.This reconstruction used Empirical Orthogonal Function (EOF) analysis to separate the two ENSO flavors, which is the most common method to separate the two ENSO flavors.However, we acknowledge that EP and CP events can spill into each other's zones, and the identification of EP and CP ENSO may vary depending on the methodology used (Kao and Yu, 2009).
For the PNA pattern, we used Liu et al.'s (2017) millennium-long winter (December-March) PNA reconstruction developed with the help of tree-ring records from different PNA-sensitive regions.We chose this reconstruction because it is the only PNA reconstruction that covers the entire historical period we studied and because it incorporated diverse geographic sources providing potentially a more accurate depiction of the PNA.
For the PDO, we used Macdonald and Case's ( 2005) annual (January-December) tree-ring-based reconstruction based on midlatitude western North American trees.While this reconstruction aligns with the prevailing methodological approach of using tree-ring data to reconstruct this index, we acknowledge the potential impact of reconstruction selection on our findings.As highlighted by Kipfmuller et al. (2012), different reconstructions of the PDO have led to varied conclusions in previous studies regarding its association with regional fire activity.Hence, interpretations of historical PDO effects may entail a degree of uncertainty.Nonetheless, this study provides a valuable starting point for understanding these complexities.
For analyses focusing on the modern record, we used monthly observational data of the climate indices EP ENSO, PNA, and PDO.To represent the canonical EP ENSO, we employed a qualitative approach, considering the region where the action center for the EP occurs.Despite potential overlap between EP and CP events, our selection of the Niño 3 region (5 S-5 N and 150-90 W), aligns closely with the eastern equatorial Pacific, making it an appropriate index to measure EP ENSO (Kug et al. 2009).To ensure consistency with reconstructed indices, we averaged monthly observational data over the corresponding months.Specifically, we averaged December-February for ENSO, December-March for PNA, and January-December for PDO.Monthly observational data for all climate indices was obtained from the National Oceanic and Atmospheric Administration (NOAA)-Climate Prediction Center (CPC) website (https://www.cpc.ncep.noaa.gov/data/indices/).

Historical fireclimate associations
To test H1, we analyzed the associations between HFAYs within the RPDA and climate index states or combinations of states using contingency analysis with non-parametric bootstrapping.To this end, we categorized each climate index reconstruction chronology into positive or negative "states" (Bai et al. 2012) based on whether the annual climate index value exceeded or fell below its mean during the study period.While this binary classification of climate states did not differentiate between a strong index value and a moderate one, it facilitated standardized comparisons across all climate indices.It also enabled us to extend this analysis to the modern fire record (see next section), spanning only 60 years, as it ensured sufficient frequency of occurrences in all climate state combinations.
We examined the associations between HFAYs and states of single climate indices as well as pairwise combinations (ENSO-PNA, ENSO-PDO, and PDO-PNA).We limited the number of indices to two in a single analysis to maintain meaningful frequencies in the contingency table, n ≥ 5, for robust bootstrap routines.We ran the analyses for the whole historical period 1700-1900 and for two subperiods, 1700-1800 and 1800-1900.We summarized the occurrence of HFAYs and non-HFAYs for each climate state and pairwise combination of climate states to identify historical fire-prone climate states where the occurrence of HFAYs exceeded chance expectations.

Fire-prone climate states in the modern record
To validate the fire-prone climate states identified in the reconstructed record, we extended the analysis to the instrumental record, employing the same methodology.For the fire data, we used a multisourced burned area record from the Canadian National Fire Database (CNFDB).This dataset covered only a portion of our study area that corresponds to Quebec and Ontario, but provided a relatively long temporal coverage (1959-2021) (Canadian Forest Service 2019).We defined HFAYs as years within the upper 20 % of the distribution of the annually burned area.We chose this threshold to capture the most extreme fire activity years while ensuring we had enough occurrences for meaningful statistical analysis.For the climate data, we focused on modern observations of the canonical EP ENSO, PNA, and PDO indices.
We examined associations between the fire-prone climate states and fire-conducive conditions across the study area (H2).The fire-prone climate states studied in this analysis were those identified in the reconstructed record and validated on the modern record.We ran spatial composite analyses for each fire-prone climate state and the 500 hPa geopotential height, drought, temperature, and precipitation fields for the winter-early spring and for the late spring-summer season over 1950-2021.We obtained the 500 hPa geopotential height from the gridded (0.25 • resolution) fifth generation European Centre for Medium-Range Weather Forecast (ECMWF) Re-Analysis (ERA5, Hersbach et al. 2020), the temperature and monthly precipitation data from the Climatic Research Unit gridded (0.5 • ) Time Series (CRU TS) Version 4.04 (Harris et al. 2020), and the drought data from the global land map of monthly self-calibrated Palmer's Drought Severity Index (scPDSI) (van der Schrier et al. 2013;Dunn et al. 2021).

Long-term trends of fire activity
To test whether the observed decline in fire activity since the end of the LIA is associated with a long-term decline in the frequency of fireprone states (H3), we conducted regime shift analyses on the occurrence of fire-prone climatic states.We identified these states in the reconstructed record and validated them with modern instrumental record of area burned and local climate conditions.We employed Rodionov's Regime Shift algorithm ( 2004) independently for each fire-prone climate state over the period jointly covered by reconstructions and modern data.We attempted to separate the spring (March to May) from the summer (June to August) HFAYs, but discovered that they are the same, so we used only one set of HFAYs.We used an L parameter of 10 and binned the data at 5-year intervals.The Hubert weight parameter was set to 1, and the α = 0.05.This analysis produced binary annual chronologies for all fire-prone states according to whether or not (1/0) the states occurred for the analyzed time period.

Results
We identified 32 HFAYs for the reconstructed historical record .The threshold for detecting significant synchrony in the annual number of sites recording fire varied across the different subperiods: ≥6 for 1700-1749, and ≥8 for 1750-1799, 1800-1849, and 1850-1900 (Table 1).

Historical fire association with the local climate
During historical HFAYs, the winter-early springs and the late springsummer season were associated with mid-troposphere ridging in the northern part of the study area (p < 0.1) (Fig. 2A).The late springsummer season was also associated with increased temperatures (Fig. 2A).Historical HFAYs showed lower PDSI levels for a modest region in the western parts of our study area (Fig. 2B).

Fire-climate associations in the modern record
Contingency analysis revealed five fire-prone climate states for the modern period : ENSO+, PDO-, PNA+, ENSO+/PNA+, and PDO-/PNA+, (Table 2 and Table SI 5).Two of these states, ENSO+/ PNA+ (Fig. 3) and PNA+ (Fig. 4), also appeared consistently in both subperiods and the entire reconstructed record.Consequently, we examined local climate conditions during years characterized by these two states from 1950 to 2021.
Both ENSO+/PNA+ (Fig. 5) and PNA+ (Fig. 6) were associated with high fire hazard weather.Specifically, during ENSO+/PNA+ years, midtroposphere ridging prevailed over northern areas with troughing in the south during January-April, transitioning to widespread ridging during May-August.Positive temperature anomalies were mainly observed in western areas for January-April but expanded across the region during May-August.Negative precipitation anomalies exhibited a more coherent pattern during January-April than during May-August.Negative PDSI anomalies were concentrated in the northwest areas during both periods, January-April and May-August (Fig. 5).
State PNA+ exhibited a similar Z500 field pattern, with ridging in the north and troughing in the south for January-April, transitioning to widespread ridging during May-August (Fig. 6).However, the coherence of positive temperature anomalies, negative precipitation and PDSI anomalies was lower for this state than for ENSO+/PNA+.Notably, the negative precipitation anomaly during PNA+ years was coherent during January-April but patchier during May-August (Fig. 6).

Long-term trends of fire-prone states
The fire-prone climate state ENSO+/PNA+ was frequent from 1700 until around the mid-1700s and from around 1985 to around 2020.However, from 1750 to 1985, its frequency was low (Fig. 7).The climate state PNA+ was frequent from the 1750s until just before the 1850s, and then again from after the 1900s until 2021 (Fig. 7).Around 1990, PNA+ became even more frequent.Consequently, from the 1850s to the 1900s, the two fire-prone climate states were at a low frequency.In general, from the oldest part of the chronology until just before 1850 and then from after 1900s to 2021, at least one of these two fire-prone climate states was in a high-frequency mode (Fig. 7).

Discussion
In this study, we explored the effect of large-scale modes of atmosphere-ocean variability over the Pacific domain upon fire activity in the red pine forests of eastern North America.We detected PNA+ and PNA+/ENSO+ as fire-prone climate states (H1).These climate states were associated with a fire-conducive surface climate (H2).Our results do not show a decline in frequency of fire-prone climate towards the end of the LIA.However, such states became more frequent in the 1900s, and further increased around 1990.This pattern prevented us from drawing

Table 1
Synchrony thresholds in the annual fire activity for four sub-periods and for the entire period with reconstructed data.We used thresholds at the 0.999 quantile to identify high fire activity years (HFAYs).Synchrony levels refer to number of sites recording fire in a single year.a definite answer to H3.Since the climate-driven changes in fire activity likely affect other ecosystems of this part of the North American continent, our results may have implications beyond mixed pine forests.

Reconstructed fire activity and local climate
The reconstructed HFAYs were associated with fire-conducive weather (mid-troposphere ridging and warm temperatures) during the 1836-1900 period (Fig. 2).The association with the summer PDSI over the whole reconstructed period (1700-1900) was modest and limited to the western parts of the study area.Inherent limitations of the reconstructed PDSI record might explain this result.First, the Midwest and northern New England were the weakest calibration areas of the PDSI grid of the contiguous U.S. due, possibly in part, to weak precipitation signal in tree-ring data in these regions (Cook et al. 1999).Second, the NADA's PDSI in our region has a weak to absent winter precipitation signature (St. George et al. 2010).Thus, important processes for the fire season that occur during the cold season, i.e., the dynamics of recharge of soil water reserves, are possibly not accounted for by this reconstruction the RPDA.

Fire-climate associations
At the annual scale, climate states characterized by two indices had a statistically stronger effect upon fire activity of the RPDA than states described by single climate indices.For example, while ENSO+ was not   and the 500 hPa geopotential height (Z500, m2/ s2) and temperature (K) from the Twentieth Century Reanalysis version 3 (20CRv3) dataset.Colored contours represent the significant anomalies (p < 0.10) from the 1836-1900 averages.The results are presented for winter-early spring (January-April), and late spring-summer (May-August).(B) Annual spatial composite analysis of HFAYs (1700-1900) and the PDSI reconstruction from the North American Drought Atlas (NADA).Note that color bars are not always centered on 0. identified as fire-prone in the reconstructed record, it was identified as a fire-prone climate state when in combination with other climate indices (ENSO+/PNA+).These results further highlighted the value of looking at composite climate states, i.e., those represented by a combination of teleconnections (Brown 2006;Sibold and Veblen 2006).In our region, a study of teleconnections associated with the dynamics of Great Lakes ice cover found that even when the single indices, that included ENSO and PDO, showed no significant associations with ice cover, their effect in combination was significant (Wang et al. 2018).
Positive states of ENSO and PNA showed consistent association with HFAYs both in the historical and modern records (Table 2).During EP El Niño events, the AL deepens and is also displaced eastwards to the Gulf of Alaska (Rodionov and Assel 2003).This leads to a PNA configuration with an eastwards-displaced ridge closer to the RPDA, resulting in milder winters and springs in the region (Rodionov and Assel 2003).In line with these dynamics, our analysis identified the combination ENSO+/PNA+ in the reconstructed and the modern record (Fig. 3) as a fire-prone climate state and the state that had the most coherent patterns with fire-conducive local climate in the region (Fig. 5).However, in the subperiod 1800-1900 and the entire period (1700-1900), the combination ENSO-/PDO+ was also associated with HFAYs (Table 2).This result might be explained by the ENSO diversity phenomenon, which causes variability on the extratropical effects among individual ENSO events (Capotondi et al. 2015).
The consistent identification of the positive phase PNA as a fireprone state, whether alone or in combination with other climate indices (Table 2), can also be explained by its impact on the midlatitude jet stream, the PFJ.The PNA pattern forces a meridional jet stream flow, which results in fewer Pacific-derived winter storm tracks and dry anomalies primarily in western North America (Rodysill et al. 2018).This effect may extend to our study area and bring drought conditions.

Table 2
Results of contingency analyses of years with high fire activity (HFAY) and combinations of circulation indices for the entire reconstructed period , two sub-periods of the reconstructed period (1700-1800 and 1800-1900), and the modern times .Only fire-prone climate states significant at 0.90 are shown.The PNA pattern is also out-of-phase with the Tropical Northern Hemisphere (TNH) pattern, which features a deep trough over Hudson Bay associated with cold and moist conditions over eastern North America (Mo and Livezey 1986;Soulard et al. 2019).Thus, the positive PNA entails a negative TNH, preventing formation of the trough and associated wetter conditions in our study area.
In the reconstructed record, the fire-prone climate states involving the PDO predominantly featured the PDO in its positive state (PDO+, ENSO-/PDO+, and PDO+/PNA+) except for one (PDO-/PNA+).The positive PDO phase is associated with a deeper AL (Mantua et al. 1997), and thus with the PNA pattern, which brings milder winter conditions to the RPDA.However, in the modern record, the identified fire-prone climate states only featured the PDO in its negative state (PDO-and PDO-/PNA+).This discrepancy might be due to the relatively short length of the modern fire record, spanning only 60 years, potentially insufficient to capture PDO cycles.Another factor that might contribute to the lack of consistency between results obtained on the reconstructed the modern records could be a less significant PDO effect at the annual scale as this index is known for its decadal periodicities (Mantua et al. 1997).

Fire-prone climate states and local climate
Years with fire-prone climate states were characterized by a fire conducive local climate in the modern record .The fire prone climate states ENSO+/PNA+ and PNA+ were associated with high fire hazard, as indicated by positive mid-troposphere height and temperature anomalies, and low precipitation and PDSI anomalies for the winter-early spring and the late spring-summer season periods (Fig. 5 and 6).However, mid-tropospheric ridging especially in the late spring-summer season was the most widespread over the area and thus the strongest pattern for all three climate-state combinations.A mild winter-early spring season as well as a warm and dry late spring-summer season resulting from mid-troposphere ridging.Atmospheric ridges block zonal air flow, generating a warmer and drier climate.These conditions are favorable for the drying of fuels (Fauria and Johnson 2008).A warm and dry late-spring-summer season featuring blocking ridges, high temperatures, and low precipitation also augment the fire risk.The formation as well as the breakdown of high-pressure blocking systems augment fire risk during the fire season.The formation of blocking ridges create warm, dry conditions (see above), and the breakdown of these systems is linked to enhanced lightning, strong surface winds, and limited rain to moisten fuels (Flannigan and Wotton 2001;Skinner et al. 2002).Lightning provides the natural source of ignition in high latitude forests that is responsible for most of the area burned therein (Skinner et al. 2002;Veraverbeke et al. 2017), while strong winds increase fire spread (Flannigan et al. 2000).

Low frequency trends in occurrence of fire-prone states
Regime shift analysis spanning from 1700 to 2021 revealed a decline in the frequency of fire-prone climate states at the end of the LIA, around 1850 (Fig. 7).This result supported H3 and the earlier interpretation of the decrease in fire activity in the northern periphery of the RPDA in eastern Canada starting in 1850 as a climate-driven phenomenon (Bergeron et al. 2001;Bergeron et al. 2004).However, we observed a resurgence of high-frequency regimes for PNA+ after the 1900s, with a notable further increase around the 1990s, which contradicts H3.Contrary to this pattern, the decrease in fire activity in eastern Canada continued and became more pronounced after 1920 (Bergeron et al. 2001).Fire suppression measures have been proposed as factors contributing towards a decline in fire activity during the 20th century (Nowacki and Abrams 2008).However, fire suppression has been shown to be of limited importance at least at the northern fridges of RPDA, where climate remains the dominant control of the fire regime (Bergeron et al. 2001;Danneyrolles et al. 2021).A decrease in traditional subsistence practices has been proposed to explain the decrease in fire activity in the 20th century for some areas in the Great Lakes (Kipfmueller et al., 2021).More research is clearly needed to explain why the observed low levels of fire activity in eastern Canada (Bergeron et al. 2001) and in the rest of the RPDA (Drobyshev et al. 2008a) over the 20th century do not reflect the long-term patterns of the frequencies of fire-prone climate state (PNA+) reported here.
An atmospheric circulation pattern associated with drought conditions in the area became more frequent since 1985 (Girardin et al., 2006b), which coincides with the shift to a higher frequency of ENSO+/PNA+ state around the 1980s (Fig. 7).This atmospheric pattern reported by Girardin et al. (2006b), resembles the eastward displaced  .Fields are the 500 hPa geopotential height (Z500, m 2 /s 2 ), temperature ( • C), precipitation (mm/month), and PDSI (unitless).Significant anomalies (p < 0.10) from the 1950-2021 averages are shown as colored contours.Results are presented for early spring (January-April), and late spring-summer (May-August).
PNA pattern, which is forced by ENSO EP, the ENSO configuration we focused on in this study.Another atmospheric circulation pattern increasing its frequency since 1850 was a TNH-like pattern with a trough around the Hudson Bay, bringing wetness and a decline in fire activity in eastern and central Canada (Girardin et al. 2006a).This configuration became more common during 1850s -1990 (Girardin et al. 2006a), which does not agree with the evolution of PNA+ reported here (Fig. 7).
For instance, when considering ENSO, PDO, and NAO at the interdecadal scale, only PDO explained extreme precipitation variability over Canada (Tan et al. 2016), and large wildfire

Study limitations
We acknowledge several limitations of this work.We used frequency data as opposed to spatially explicit reconstructions, which may potentially compromise our view of large fire years, defined here as years with large burned areas.Although we assumed that the synchrony of fire occurrence is a legitimate proxy of fire size and, by extension, climate forcing upon fire activity (Farris et al. 2010;Falk et al. 2011;Drobyshev et al. 2014), frequency data may not fully capture the dynamics of the burned areas, which is a better proxy of landscape and Fig. 6.Spatial composite analyses of local climate conditions during the fire-prone state PNA+ years .The fields are the 500 hPa geopotential height (Z500, m 2 /s 2 ), temperature ( • C), precipitation (mm/month), and PDSI (unitless).Significant anomalies (p < 0.10) from the 1950-2021 averages are shown as colored contours.Results are presented for winter-early spring (January-April), and late spring-summer (May-August).region level fire activity. Spaial fire reconstruction on dendrochronological data, however, are largely missing over eastern North America.
Limiting analyses to combinations of two climate indices might oversimplify and potentially miss important interactions within the climate system.It is difficult to analyze more indices in a single analysis due to the low expected frequencies of three-state combinations in the contingency tables, effectively precluding such analyses with the available record.To overcome this obstacle, we would need to considerably expand the length of the regional fire chronology by at least 200-300 years.This is, however, problematic due to a limited residency time of dendrochronologically datable remnant wood on the forest floor of increasingly mesic forests.Finally, the use of the modern fire data for the Canadian part of RPDA might introduce another source of uncertainty.Although the choice of Canadian data was justified by the larger proportion of forested areas and longer fire records, we realized that the modern fire dataset might not be fully representative of the red pine range.

Dynamics of fire activity in the red pine forests and its ecological implications
The lower frequency of ENSO+/PNA+ during the 20th century could have facilitated successional development away from xeric mixed pine forests towards mesic forests dominated by deciduous trees, reported earlier across the American Midwest (Nowacki and Abrams 2008).Fire suppression and infrequent use of prescribed fire have been viewed as the primary factors behind this conversion (Frelich et al. 2021).Our results, however, show that climate change likely also contributed to decreased fire activity starting in the 1850s, leading to mesophication of historically fire-maintained ecosystems and supporting the view of climate-driven disturbance regimes as having lasting legacy on regional species composition (Pederson et al. 2014).
The increase in the frequency of fire-prone climate states since the 1980s may foreshadow the anticipated increase in climatological fire hazard projected for eastern Canada in the future (Wang et al. 2017).Such an increase will likely interact with the legacies of the less fire-prone period, such as fuel load accumulation, further increasing probabilities for potentially large and high severity fires.Fires developing under the conditions of higher fuel loading, including ladder fuels, would likely increasingly be of higher severity and potentially stand-replacing.A similar effect of a "fire deficit" has been recently proposed for the modern Canadian boreal forest (Parisien et al. 2020).For the red pine, high-severity fires generally remove on-site seed sources and limits post-fire recovery (Van Wagner 1971).Red pine has a limited seed dispersing capacity and, therefore, is highly dependent on the presence of scattered surviving individuals for recolonization after fires (Horton and Bedell 1960).Scarcity of red pine seed banks and high fire severity would likely favor regeneration of more easily dispersed species.In the boreal and sub-boreal section of red pine range, jack pine will likely replace red pine (Nyamai et al. 2014).In fact, it is a change from surface to stand-replacing fire as the dominant fire type at the northern fringe of red pine distribution that has been attributed to the disappearance of red pine from the forest canopy (Bergeron and Brisson 1990;Flannigan and Bergeron 1998).To this end, prescribed fires could be a critical element in the conservation policies aiming to maintaining red pine forests during the transition to more fire-prone climate conditions of the present and future (Montour et al., 2020).
It is likely that the legacies of fire-free periods in shaping future disturbance regimes are not limited to red pine ecosystems.Humanrelated or climatically driven modification of successional pathways, leading to changes in fuel loads, may become increasingly important for the control of future disturbance regimes, and specificallyof the fire behavior.We speculate that the feedbacks in the climate-vegetation-fire system likely affect the strength of association between large-scale climate drivers and the properties of current and future disturbance regimes.

Declaration of competing interest
The authors declare the following financial interests/personal relationships which may be considered as potential competing interests: Igor Drobyshev reports financial support was provided by University of Quebec at Abitibi Temiscamingue.Fig. 7. Regime shift analysis for each fire-prone climate state for the period 1700-2021.Red circles represent years with fire-prone states of climate indices, and the blue circles are the years with non-fire-prone conditions.The thick black lines marks climate regimes, which are defined based on the proportion of fire-prone vs. non-fire prone states for a specific circulation index.A higher position of the black line for a specific index indicates a regime with a higher frequency of fire-prone states.Note that the vertical position of black lines in relation to circles is chosen with the only purpose to facilitate visibility of graphical elements.The black line is the observed area burned in Canada for the period 1959-2021, based on the data from the Canadian National Fire Database (CNFDB).

Fig. 1 .
Fig. 1. (A) Map of the studied region with distribution range of red pine shown in pink, and sites with dendrochronological reconstructions of forest fire history, red dots.(B) Fire chronologies from sites included in the analyses.The black dots denote fire years.The vertical red lines represent the High Fire Activity Years (HFAYs).The blue lines delimit the study period.

Fig. 2 .
Fig. 2. (A)Seasonal spatial composite analyses of the reconstructed High Fire Activity Years (HFAYs) and the 500 hPa geopotential height (Z500, m2/ s2) and temperature (K) from the Twentieth Century Reanalysis version 3 (20CRv3) dataset.Colored contours represent the significant anomalies (p < 0.10) from the 1836-1900 averages.The results are presented for winter-early spring (January-April), and late spring-summer (May-August).(B) Annual spatial composite analysis ofHFAYs (1700HFAYs ( -1900) )  and the PDSI reconstruction from the North American Drought Atlas (NADA).Note that color bars are not always centered on 0.

Fig. 3 .
Fig. 3. Contingency analysis of the climate states, classified into positive and negative phases, of ENSO and PNA and high fire activity years (HFAYs) for (A) 1959-2021 and (B) 1700-1900.The maroon vertical lines represent the HFAYs.The scatter blue points represent the ENSO values and the red points represent the PNA values.The blue horizontal lines represent the mean of the ENSO chronology, and the red horizontal lines represent the mean of the PNA chronology.The color bars code the combinations ENSO-/PNA-, ENSO-/PNA+, ENSO+/PNA-, and ENSO+/PNA+.

Fig. 4 .
Fig. 4. Contingency analysis of the climate states, classified into positive and negative phases, of PNA during high fire activity years (HFAYs) for (A) 1959-2021 and (B) 1700-1900.The maroon vertical lines represent the HFAYs, the scatter points represent the PNA values, with the horizontal line representing the mean of the PNA chronology.The color bars code the states PNA-and PNA+.