Recurrent pattern of extreme fire weather in California

Historical wildfire events in California have shown a tendency to occur every five to seven years with a rapidly increasing tendency in recent decades. This oscillation is evident in multiple historical climate records, some more than a century long, and appears to be continuing. Analysis shows that this 5–7 year oscillation is linked to a sequence of anomalous large-scale climate patterns with an eastward propagation in both the ocean and atmosphere. While warmer temperature emerges from the northern central Pacific to the west coast of California, La Niña pattern develops simultaneously, implying that the lifecycle of the El Niño-Southern Oscillation that takes multiple years to form could be a trigger. The evolving patterns of the Pacific-to-North America atmospheric teleconnection suggest the role of tropical and subtropical forcing embedded in this lifecycle. These results highlight the semi-cyclical hydrological behavior as a climate driver for wildfire variability in California.


Introduction
According to the California Department of Forestry and Fire Protection (CAL FIRE), most of the top 20 fires in California occurred since the year 2000 (bars in figure 1(a)). On top of that, the complex fire which occurred on August 2020 broke records as one of the most severe fires in California (CNN 6 , 6 October 2020), and authorities are now concerned of upcoming exceptional danger of drought and fire in the US west (NOAA 7 , 29 April 2021). These series of catastrophic events urge to identify sequential pattern of fire related climate phenomena to restrict fire damage ahead. Upon examining the burned area and number of fire statistics from the Monitoring Trends in Burning Severity (MTBS, lines in figure 1(a)), it was observed that significant wildfire damages in 6 https://edition.cnn.com/2020/10/06/us/gigafire-californiaaugust-complex-trnd/index.html. 7 www.ncei.noaa.gov/news/us-drought-monitor-update-april-27-2021.
California tend to occur every five to seven years since the late 20th century (grey columns in figure 1(a)), a feature that could be periodic in nature that is studided herein.
Fire is intimately tied to ecosystems through interactions with vegetation and climate. As an ecosystem process, fire regulates spatial distribution and composition of vegetation (McLauchlan et al 2020). Vegetation structure, in terms of fuel source, influences fire regimes by selecting particular traits and species to survive within a given fire regime (Rogers et al 2015) and therefore fire probability and severity. Climate also has a major influence on fire regimes across diverse scales ranging from short-term fire weather to seasonal and decadal variability. Particularly in California, where the wet-dry season is distinct, climate variability plays a crucial role in fire weather condition. Considering this interactive relationship, the apparent oscillation in California raises a possible physical coupling between fire, vegetation and climate. Due to the documented enhancement in fire weather associated with the increasing temperature (Jolly et al 2015), historical observation and climate model experiments have indicated an aggravated risk in widespread fires under climate change Liu et al 2010. Natural climate variability that modulates California's water cycle extremes strongly modulates its fire potential and occurrence . The underlying climate variations in California and their link with the Pacific Ocean have been widely analyzed, and such natural variability, with an important contribution from the evolution cycle of the El Niño-Southern Oscillation (ENSO), affects California's precipitation (Rajagopalan and Lall 1998), vegetation carbon uptake (Keeling et al 1995) and groundwater . Here, we present observation-based climate diagnostics of the puzzling quasi-periodic variation of California's fire occurrences (figure 1(a)) and its connection with the recent drought cycle.

Meteorological data
The Vegetation Health Index (VHI) (Kogan 1995) is a combined estimation of moisture and thermal condition using the Normalized Difference Vegetation Index (NDVI) and the brightness temperatures provided from Blended Vegetation Health Products. The live fuel moisture (LFM) (Pollet and Brown 2007), which is a ratio of weight between the fresh and dry vegetation, represents the amount of water content and the flammability of vegetation (Dimitrakopoulos and Papaioannou 2001). The LFM database is in-situ measured and maintained by the U.S Forest Service-Wildland Fire Assessment System. For the risks of wildfire and drought, we use the Keetch-Byram Drought Index (KBDI) (Keetch and Byram 1968) provided from ERA5 (Hersbach et al 2020(Hersbach et al ) (0.25 • × 0.25 • , 1979 and self-calibrating Palmer Drought Severity Index (scPDSI) (Wells et al 2004) from NCAR (Dai et al 2004)  . To access longterm drought data, summer averaged Palmer Modified Drought Index (PMDI) (Palmer 1965, Karl 1986, Heddinghaus and Sabol 1991, reconstructed by treering, from Living Blended Drought Atlas (Cook et al 2010) is also used for 1500-2017. Climatic Research Unit (CRU) also provides long-term land domain temperature and precipitation (0.5 • × 0.5 • ) for time since year 1901 (Harris et al 2020). Also, we use a precipitation dataset from the Japanese 55 year Reanalysis (Kobayashi et 1854-present), and geopotential height is from JRA55. The soil moisture content at 100 cm depth is compared with three different land surface models from the North American Land Data Assimilation System (NLDAS-2): Mosaic (Koster and Suarez 1992), Noah (Chen et al 1996) and VIC (Liang et al 1994).

Multiple-taper spectrum estimation method with singular value decomposition (MTM-SVD) analysis
An important part of the proposed analysis is to identify the shared frequency in the coupling between fire, vegetation, and climate. The frequency domain SVD provides coherent climate structure across a multivariate dataset. Before a set of climate variable series are decomposed into orthogonal modes, each component is standardized by removing climatological seasonal cycle and dividing by the standard deviation. The constituent series are transformed the spectrum domain using multitaper spectral analysis (Thomson 1982, Park et al 1987, Mann and Lees 1996. With these series x, the orthogonal sequences of K data tapers are calculated at each frequency f, where ∆t is the sampling interval (1 month), M is the number of grid points and N is the time dimen- is the kth member in orthogonal sequence of K data tapers. Each taper is based on a frequency band of half-bandwidth of pf R about a given frequency f, where f R = (N∆t) −1 is the Rayleigh frequency. Here, we choose p = 2, K = 3, which provides reasonable frequency resolution and sufficient spectral degrees of freedom Park 1994, 1996).
The K eigenspectras at the M grid points are organized as M × N matrix A(f ), where the w is latitude weightings in proportion to grid point area. The matrix A(f ) is decomposed by means of a complex SVD as into orthogonal modes, where λ k represent the relative fraction of total variance explained by the kth mode, left-eigenvectors u k represents empirical orthogonal functions in the spatial domain, and right-eigenvectors v k represents in the spectral domain.

Signs of recurring pattern over the recent 40 years
The recurring pattern is shown in various wildfirerelated observations over California, including vegetation indices. The VHI in figure 1(b) (black line) appears to fluctuate similarly yet opposite with burned areas in figure 1(a). The trend of decline in VHI for the last 20 years is consistent with the increase of major fire damages (r = −0.45, since year 1987). These results describe that a drier state of vegetation provides more combustible fuel loads and a higher potential for large wildfires. Given that VHI includes greenness, it may merely reflect losses of vegetation by the wildfire occurrences rather than its flammability. However, we found consistent results from the LFM (r = 0.67 with VHI and −0.33 with burned area, since year 1987), which was manually collected from hundreds of sites in California to measure the water content in vegetation. As shown in figure 1(b), LFM supports that the aridity is closely associated with the danger of extreme wildfires and both variables depict a 5-7 year frequency (see the spectral analysis next). The implication from this correspondence between the burned area data and the vegetation growth factors is twofold: (a) fires tend to occur at low-moisture/drier years and (b) large fires trail high-growth years by 1-2 years with more fuel. A similar oscillating pattern also appears in three drought indices (figure 1(c)): KBDI, scPDSI and PMDI. These 'drought cycles' are in agreement with the previous indices in that more (less) vegetation grows in less (more) dry/drier years, which corresponds to the trailing enhancement (suppression) of fires, and these appear to occur in a cyclic manner.
To examine the correspondence among fire, vegetation and climate, we further conducted the cross-spectral analysis. Figure 2 shows that the crossspectral power of VHI with all other indices of environmental conditions, including soil moisture (figure S1 (available online at stacks.iop.org/ERL/16/094031/ mmedia)), reveals common periodicities that are consistently within the 5-7 year frequency band. This finding supports the inference made from figure 1 that the apparent oscillatory feature of fire intensity (burned areas) is potentially driven by a predominant intra-decadal climate oscillation. One of the prominent climate drivers of this oscillatory feature in California can be linked to the ENSO lifecycle, including its precursor and decay phase, that appears to have amplified since the late 20th century in the warming climate .
The power spectral analysis of the aforementioned variables (figures S2(a)-(d) shows a common spectral peak within the 5-7 year frequency band. The 5-7 year spectral power is significant at p < 0.05 among all the variables regardless of their period of record, suggesting that a 5-7 year 'climate oscillation' is a persistent feature. Among the climate factors that affect fire weather in California, precipitation exhibits a more energized 5-7 year oscillation than surface air temperature (figures S2(e) and (f)), echoing a previous observation (Los et al 2001) that the NDVI has a closer link with precipitation than temperature. In terms of soil moisture (figures S2(g)-(i)), the power in the 5-7 year frequency band is even more pronounced, ecohing the previous observation (Wang et al 2016) that the land surface processes damp atmospheric signals of higher frequencies, and the longer-term variability of precipitation remains subsequently 2015. Additional frequency analysis of the longer-term data is shown next, in section 3.2.

Validation of the frequency in the long term records
To examine whether the modern-era intra-decadal climate oscillation is similarly robust before the 21st century, we conducted the wavelet spectral analysis (Torrence and Compo 1998) for the available datasets that have longer than 100 years of record. As shown in figures 3(a) and (b), both precipitation and surface temperature show significant oscillatory characteristics in the 5-7 year frequency band throughout the 20th century. Despite their episodic occurrence, each period with an active 5-7 year oscillation in precipitation lasts around 15 years and recurs every 30 years, corresponding to the post-2000 condition. The episodes of increased 5-7 year power are even more pronounced in scPDSI in its 150 years of data (figure 3(c)), with a recurrence interval of ∼30 years.
To put the apparent oscillation into the preinstrumental perspective, we further examined a drought proxy established from tree rings (Keen et al 2020) of more than 500 years, the PMDI ( figure 3(d)). The wavelet spectrum of the PMDI reveals repeated occurrences of the amplified 5-7 year powers from the year 1500-2000. Nonetheless, the common significant results in the 5-7 year band suggest that the recurring wet-dry pattern is not a new phenomenon due to climate change, but rather is a part of natural variability due to a clear signal of 5-7 year power in precipitation rather than temperature (figure S3). However, climate change and human factors have apparently amplified this climate oscillation's impact on droughts (Diffenbaugh et al 2015) and wildfires , which can be clearly seen in increasing temperature and associated fire risk , Son et al 2021 indicated by KBDI ( figure 1(c)).

Climate forcing of the oscillation
To examine possible sources of the 5-7 year recursive fire weather conditions, the first line of analysis was conducted using regression (figures S4-S6). We find that the recurring aridity in California is associated with atmospheric-ocean coupled patterns over the subtropical North Pacific. Anomalously high SST appears in the Central North Pacific and propagates during about two-year period over the mid-latitudes adjacent to Western North America (figures S4(e) and (f)). This oceanic transport is closely matched by the track of land-approaching precipitation deficits from the Central North Pacific to California (figures S5(d)-(f)), along with anomalously high pressures in 850 hPa (figures S6(d)-(f)). These results reflect the relationship in SST forcing and atmospheric variability that induced the previous extreme droughts in California (Wang andSchubert 2014, Seager et al 2015).
To provide a more concise depiction of these results, we perform the MTM-SVD on major climate components on a global scale. The regression and the MTM spectral analysis are calculated between VHI in California and all other variables, while the illustration of the MTM analysis (figure 4) is focused on the 4-8 year frequency range to account for truncation. In figure 4, the amplitude and phase of the MTM coherence are visualized as the length and direction of vectors, respectively, following previous studies Park 1994, Wang et al 2012). The MTM-SVD spectral coherence summarizes the lead-lag relationship as phase vectors in figure 4(a). The 0 • vector indicating the north represents concurrent phase with the trough of VHI in California. The degree of angle describes phase difference with VHI, for instance, the 180 • vector directing the south means the opposite phase with VHI. In the coast of California, the vectors show about 90 • directing the east (a quarter-phase shift) implying warmer SST appears roughly 1.5 years earlier than the aridity in California, while the middle of the Pacific at mid-latitude is 270 • (opposite phase relative to coastal California). All data are considered only for California region. Linear trend is removed and 12 months moving average is applied in advance. Wavelet spectrum is based on Morlet parameter-6 approach and contour black line is for 95% confidence level (CI, p < 0.05). The red noise in 95% CI are hatched.
These two opposite phases are smoothly connected with anticlockwise rotation from the middle of the Pacific to the coast of California, accompanying a pair of negative-positive SST anomalies. These results portray that the SST anomalies propagate eastward during the few years leading up to the low VHI or drought year in California (figures S4(b) and (e)).
Similar with the SST, the precipitation shows negative anomalies in the Western North America ( figure 4(b)) accompanied by a La Niña-like response in the equatorial Pacific. The developing La Niña pattern in SST and precipitation (figures S5 and S6(d)-(f)) suggests that the ENSO cycle may serve as the driver of this multi-year propagation (Wang S-Y et al 2015). However, the 0 • vectors around California are continued rotating not only to the tropical warm pool region, but also to the North Western Pacific, which is discussed as an ENSO precursor (Wei et al 2016). This result, in turn, echoes the previous observation that the extreme phasing of ENSO is not the only factor responsible for the climate and vegetation variability in California , Capotondi et al 2019.
In terms of the atmospheric teleconnection, we observe a consistent variation in the 850 hPa geopotential height regression, which shows a significantly amplified pattern (figure 4(c)) representing the North American dipole with a broadened ridge in the west and a deepened trough in the east (Wang et al 2015 . Regression and MTM-SVD spectral coherence based on VHI in California. Regression (shading) and maximum coherence (vector, amplitude is for length and phase is for direction) between VHI in California and anomalies of (a) SST, (b) precipitation, (c) geopotential height 850mb (zonal mean is removed). Only statistically significant (p < 0.01) regression and coherences are visualized. Linear trend is removed and 12 months moving average is applied in advance. and Singh et al 2016). The spectral coherence displays 180 • shifts of phase between California (directing north) and the north central Pacific (directing south), suggesting a wave train across the North Pacific (figures S6(b) and (e)). Two years prior to the lowest VHI, the geopotential height pattern reaching western North America (figure S6(e)) resembles the Pacific North American (PNA) (Renwick and Wallace 1996) that is known to be associated with California's drought (Lin et al 2017) and wildfire (Trouet et al 2009) . The noticeable tropical phase vectors pointing 180 • away in between the Dateline (270 • to the west and 90 • to the east) suggest an ENSO pattern about two years before the driest state of vegetation in California.

Putting all together
To examine California's VHI connection with prominent climate modes, we provide the correlation matrix in figure S7. VHI has the correlation coefficients of above 0.3 with all ENSO-related indices, such as Niño 4 of the previous year, the Pacific Meridional Mode (PMM) SST index (Chiang and Vimont 2004) and the Trans-Niño Index (Trenberth and Stepaniak 2000) without lead (figures S7(a)-(c)). We additionally compare other known atmosphere-based indices, which have significant climatological effects in North America, including the North Atlantic Oscillation (NAO) and 'ridge-trough dipole index'  (figures S7(d)-(f)). Among them, the PNA shows the highest correlation (0.4), three months prior to the vegetation peak in California. It has been known that the positive PNA in winter leads to a warmer climate in California, a decline of the snowpack and less soil moisture in the summer, resulting in more aridity and wildfire risks (Abatzoglou 2011). However, our results are contradictory, showing a significant positive correlation (0.56, p < 0.01) between spring PNA and summer VHI (table S1). When the spring PNA is regressed on the geopotential height in 850mb, the positive phase of PNA pattern appears on the atmosphere ( figure S8(a)). The ocean reflects El Niño pattern in the previous summer in relation to the spring PNA ( figure S8(b)), being more pronounced through the comparison of the time-series of Niño 4 and time lagged PNA (figure S9). Considering these results, even though there have been arguments that ENSO and PNA Figure 5. Diagram of the teleconnection hypothesis for hydrological oscillation in California. The lower/upper bound of wildfire risk in California is divided into blue/red. The premature state (∼2 years before) for each phase is represented by light blue/red, the early warning state (3∼6 months before) is by blue/red and the peak state is by dark blue/red. Each progress is marked with a colored inverted triangle on the cycle line and climatological features are summarized on the tables. are independent of each other (Straus and Shukla 2002), we infer that the positive PNA pattern in spring, which plays a role of precursor of vigorous stage of vegetation in California summer, is associated with ENSO dynamics. Table S1 summaries statistically significant relationship with time differences between ENSO related atmospheric-ocean circulation and hydrological components, such as precipitation and snowpack. Only the spring PNA and the snow depth, however, are not significantly correlated. This is probably obvious because high pressure is enhanced over the western US in the positive PNA phase, inducing warmer and drier state in California (Ault et al 2011). Instead, the negative phase of NAO, which is shown in figure S7(e), is known for the spring snowmelt to be delayed in the negative NAO phase (Myoung et al 2015(Myoung et al , 2017. A combination with the North Pacific Oscillation is also suggested to explain the drought process in California (Lin et al 2017).
Even though there is no consensus yet on what climate modes govern the variability of vegetation phenology in California, several combinations between oceanic and atmospheric oscillations are suggested for the better depiction of the climate variability (Lin et al 2017, Liu et al 2018. In our study, the results of the MTM-SVD showed that all the variables have similar phase differences with roughly 180 • rotation from the common starting point in the tropical warm pool to California and to the Peruvian coast ( figure 4). Also, the pattern of the eastward propagating wave train originated from the western Pacific ( figure S6(b)) is also of interest and deserves further exploration. Based on these results, we assume that the PNA-associated patterns occurring in the midlatitude may share the same source with the full lifecycle of ENSO, which with emphasis in the western Pacific, which has been suggested as a trigger for the PNA (Straus andShukla 2002, Soulard et al 2019).

Conclusions
Wildfire has been one of the major natural disasters damaging California. Unusually large wildfires in the early wildfire season with the driest winter has raised concerns on this year's fire season. In recent years, the episodic catastrophes lead us to detect the 5-7 year frequency in agreement with multiple fire-related climate and vegetation components. Furthermore, this frequency is attested with several long-term records of more than 100 years, proving that they are not a coincidence, nor will they stop occurring in the future. We have shown a potential linkage between the frequency and teleconnections, such as ENSO and PNA. The ENSO pattern and its lifecycle are divided into two types: zonal and meridional modes. The meridional mode has a greater influence on the growth stage of El Niño (Di Lorenzo et al 2015). The PMM, which is more likely to trigger the Central Pacific type of ENSO, forms positive SST anomalies in the central-eastern North Pacific, 1.5 years prior to the maturing of La Niña. The similarities between the contrasting SST pair with a twoyear lead (figure S4(e)) and the meridional mode strongly suggest a connection with the developing phase of ENSO and California's vegetation growth, i.e. one that includes the transition or precursor of ENSO. Figure 5 summarizes the time flow for the 5-7 year oscillation of fire-related climate anomalies in California. In the early stages of La Niña (El Niño), atmospheric high (low) pressure propagates from the northern central Pacific. Once the high (low) pressure reaches the west coast, 3-6 months prior to the dry season in California, there are less (high) winter precipitation and snow. This is followed by negative (positive) PNA in the next spring. By the upcoming summer, the high (low) pressure will mark more (less) distress the vegetation with increasing (decreasing) temperature, eventually providing more (less) vulnerability hydroclimate environment to the extreme wildfires in California. The dynamical linkage of these empirically revealed processes requires further modeling studies to verify.

Data availability statement
The data that support the findings of this study are available upon reasonable request from the authors.