Summer and Fall Extreme Fire Weather Projected to Occur More Often and Affect a Growing Portion of California throughout the 21st Century

: Annual burned area has increased in California over the past three decades as a result of rising temperatures and a greater atmospheric demand for moisture, a trend that is projected to continue throughout the 21st century as a result of climate change. Here, we implement a bias-correction and statistical downscaling technique to obtain high resolution, daily meteorological conditions for input into two ﬁre weather indices: vapor pressure deﬁcit (VPD) and the Canadian Fire Weather Index System (FWI). We focus our analysis on 10 ecoregions that together account for the diverse range of climates, ecosystems, topographies, and vegetation types found across the state of California. Our results provide evidence that ﬁre weather conditions will become more extreme and extend into the spring and fall seasons in most areas of California by 2100, extending the amount of time vegetation is exposed to increased atmospheric demand for moisture, and heightening the overall risk for the ignition and spread of large wildﬁre. The ecoregion-level spatial scale adopted for this study increases the spatial speciﬁcity of ﬁre weather information, as well as the resolution with which ﬁre and land managers can implement strategies and counter-measures when addressing issues related to climate change.


Introduction
Numerous research studies indicate that wildfire activity in the western United States has increased over the past few decades due primarily to increases in temperature, decreases in precipitation, and increases in atmospheric aridity [1][2][3][4][5]. Furthermore, by the end of the 21st century, anthropogenic global warming will contribute to increased risk of extreme fire weather through the effect of rising temperatures [6]. However, the impacts of climate on fire regime vary dramatically by ecosystem, predominant vegetation types, topography, and human activity. These complex relationships between climate and wildfire are apparent in California, where climate, vegetation, and topography vary significantly from east to west and north to south.
For example, large parts of southern California experience near continuous high fire danger in the summer and fall, as high temperatures, seasonal drought, and the onset of the Santa Ana wind season combine to increase the risk of severe wildfire. In these areas, wildfire is rarely limited by weather conditions governing fire spread, but by availability of fuel to burn. On the other hand, in forested ecosystems, where deep rooted vegetation is less susceptible to short-term variations in moisture, the limiting factor may be fire spread potential or dried fuel available to burn [7,8]. Therefore, the effects of climate change on wildfire in California, specifically projected increases in temperature and aridity, will not impact all ecosystems equally, but will be strongly dependent on vegetation distribution and will be more likely to place flammability-limited ecosystems (i.e., forested) at a greater risk for catastrophic wildfire [9,10]. Furthermore, understanding the broad range of outcomes regarding the impacts of climate change on extreme fire weather risk is critical because climate change may lead to the loss of existing climate regimes and the emergence of novel ecosystem states, impacting the distribution of biomass (i.e., continuity and flammability), fire spread potential, and ignition [11,12].
One of the most valuable tools for assessing historical and projected climate variability, especially the meteorological variables that govern wildfire on a regional scale, are global climate models, or GCMs [13]. The Climate Model Intercomparison Project (CMIP) provides a standardized platform for the collection and analysis of GCM model output with the intent of studying natural and anthropogenic climate variability. The most recent release of GCM output for the CMIP Phase 6 provides a new suite of climate model data for use in a wide range of disciplines. While GCMs are often used for the evaluation of long-term variations in climate, there are significant limitations to their use for the study of wildfire impact assessment and fire weather potential. The main limitation is that climate model output is typically of coarse resolution and is subject to significant bias [13]. The mismatch between the spatial scale at which key variability in fire weather variables fluctuate (<10 km), and that of GCM output (>100 km), creates a need to account for these spatial differences.
One method of adjusting the resolution of available GCM output to that which is needed for the evaluation of fire weather variables is statistical downscaling. These methods allow for the reduction in systematic biases inherent to numerical weather prediction systems, as well as an enhancement of grid cell resolution in order to account for the effects of topography and other land surface features [14]. One method of statistical downscaling for use in wildfire risk assessment is bias correction and spatial downscaling (BCSD) [13]. BCSD is a two-step process that involves the bias correction of coarse GCM data using local observations that have been aggregated to the model's resolution [13]. Bias correction techniques include scaling, delta methods, and empirical quantile mapping. In addition, spatial downscaling utilizes techniques such as linear or multiple regression to apply statistical relationships (e.g., slope, y-intercept) derived from high resolution observational data to coarse resolution GCM data in order to improve the resolution [15]. The major assumption made with statistical downscaling methods is that of climate stationarity. It is assumed that the relationships between the coarse and fine scale historical datasets, calculated by the statistical model, will stay the same in the future (with climate change), and that when these models are applied to projected datasets the outcomes will be valid [13][14][15].
In order to track changes in, and probability of, extreme fire weather over time, fire weather indices are often calculated using weather station, reanalysis, or global climate model data. Two commonly used fire weather indices are vapor pressure deficit (VPD) and the Canadian Fire Weather Index System (CFWIS) [16,17]. The CFWIS relies on regular inputs of noon-day temperature, 24 h precipitation accumulation, 10 m wind speed, and relative humidity, and is comprised of six components: three fuel moisture codes and three fire behavior metrics, with the primary purpose of assessing the effect of weather on fuels and potential fires. In general, higher values indicate a more severe risk for ignition or spread of large wildfire [17,18].
On the other hand, VPD, defined as the difference between saturation vapor pressure and actual vapor pressure, is a commonly used metric for the quantification of atmospheric aridity in relation to wildfire risk assessment, and is highly correlated with annual area burned [1,2,19]. VPD accounts for the non-linear (exponential) relationship between saturation vapor pressure and temperature, thus better representing moisture stress and flammability [19]. In fact, maximum VPD in the entire United States is found each year in the desert southwest region, not only affecting agricultural production via crop health but increasing the risk of wildfire ignition and spread [19]. Furthermore, two-thirds of historical (2001-2018) increases in extreme VPD in the western United States is attributable to anthropogenic warming [20].
The objective of this research is to utilize bias correction and statistical downscaling to enhance the spatial representations of CMIP6 simulations for the purposes of assessing the potential impacts of climate change on the occurrence of extreme fire weather. Bias corrected and statistically downscaled daily CMIP6 data (~4 km spatial resolution), covering an historical (1981-2010), mid-century (2041-2070), and late-century (2071-2100) period, will be used as input for two fire weather indices, including the Canadian Fire Weather Index (CFWIS; now referred to as FWI) and vapor pressure deficit (VPD). In order to determine the effects of climate change on the occurrence of extreme fire weather in California ecoregions, we investigate the frequency at which daily index values exceed the average daily historical 95th percentile. We also calculate the extent of each ecoregion that is likely to experience values greater than the 30 year average spatial 95th percentile value for each time period. In addition, time series analysis includes fire weather anomalies relative to an historical period, as well as 30 year Julian-day averages. The diversity of factors that influence fire regime in California, as well as the immense range of biophysical characteristics, require an ecoregion scale analysis. The effects of climate change on extreme fire weather will not be homogeneous, however, our fine scale analysis will allow policy makers and land managers to anticipate changes in climate and determine the areas that are most likely to experience the consequences of increased fire potential.

Study Area
The primary study area is California ( Figure 1). We chose to analyze the effects of climate change on extreme fire weather within ten ecoregions of California [21,22]  Each of these regions is characterized by different dominant ecosystem types that give rise to California's diverse fire regimes, including vegetation, climate, soil types, land use, and hydrology [22]. An ecoregion specific analysis will allow for a nuanced understanding of how climate change may impact extreme fire weather relative to historical conditions within homogeneous regions of particular climate and vegetation types.

Statistical Downscaling Method
Daily near surface maximum air temperature, 24 h precipitation accumulation, 10 m Each of these regions is characterized by different dominant ecosystem types that give rise to California's diverse fire regimes, including vegetation, climate, soil types, land use, and hydrology [22]. An ecoregion specific analysis will allow for a nuanced understanding of how climate change may impact extreme fire weather relative to historical conditions within homogeneous regions of particular climate and vegetation types.

Statistical Downscaling Method
Daily near surface maximum air temperature, 24 h precipitation accumulation, 10 m wind speed, and 2 m relative humidity were acquired for three CMIP6 models (including MIROC6, EC-Earth3, and MRI-ESM2-0), for an historical (1981-2010), mid-century (2041-2070), and late-century (2071-2100) period, from the data repository at the Earth System Grid Federation (https://esgf-node.llnl.gov/search/cmip6/, accessed on 1 January 2021). We derived daily relative humidity from specific humidity, sea level pressure, geopotential height, and near surface air temperature. Future climate projections for a mid-and late-century period were obtained for the Representative Concentration Pathway (RCP8.5), "high emissions" scenario. All bias correction and downscaling methods were completed using a 4 km gridded meteorological climate dataset called gridMET [23]. We used the entirety of the gridMET historical period (i.e., 1979-2014) to train the bias correction algorithm, and in order to avoid training and testing our BCSD output with the same data (i.e., gridMET), we validated the results of the BCSD methodology with two independent observational datasets (see Section 2.3).
A bias correction technique known as empirical quantile mapping was applied separately to the four climate variables for each of the three CMIP6 model outputs. This algorithm adjusted the simulation's (i.e., CMIP6) empirical distribution based on observed patterns (i.e., gridMET) and was used to correct both historic and future climate projections generated by the CMIP6 models [14]. Bias correction methods are often applied to GCM output in order to address inherent biases that arise when modeling earth system processes. The bias corrected GCM data was interpolated to the high-resolution observational grid (i.e., gridMET) using a thin plate spline function and the relationship between the two grids were then estimated to generate regression functions for each grid cell [24]. The linear regression models were then applied to the bias corrected and interpolated GCM data to complete the downscaling.
To elaborate, this spatial downscaling methodology uses regression functions to estimate the linear relationships (in this case the slope and intercept) between the original, high-resolution gridMET data and a "smoothed" estimate, before applying those statistical parameters to the bias-corrected CMIP6 data in a simple linear model to obtain the final downscaled, high-resolution output [25]. For this research, the study area domain was comprised of approximately 70-100 grid cells in the CMIP6 simulations, whereas the gridMET observational data has~81,000 locations over the same area. In the first step, a thin-plate spline interpolation was used to predict CMIP6 values at the gridMET scale, with no addition of topographic or climatological information. This resulted in "smoothed" CMIP6 data that is identical to the original data but sampled at a much greater frequency. In the next step, a linear regression function was derived from the raw and smoothed grid-MET datasets at each grid cell location, for each month of the year (~81,000 grid cells × 12). These month-specific linear regression models were applied to the "smoothed" CMIP6 data to obtain the final downscaled product for each variable.
A square-root transformation was applied to the precipitation data before downscaling to account for its skewed distribution [26,27]. The final bias-corrected and downscaled product of each variable was averaged to obtain an ensemble mean and used in the calculation of the two fire weather indices. Daily maximum temperature, 24 h rainfall, 10 m wind speed, and daily minimum relative humidity were used to calculate FWI, while daily maximum temperature and daily minimum relative humidity are used to calculate VPD. Due to the fact that maximum temperature and minimum relative humidity are used for the calculation of VPD, we make it clear that VPD may be considered "maximum VPD".

Climate Model Simulations of Extreme Fire Weather in California Ecoregions
The BCSD CMIP6 GCM data were used for three primary analyses of extreme fire weather in California, as well as an additional analysis of Julian-day climatology. We used maximum daily surface temperature, minimum relative humidity, and 2 pm wind speed, as daily noon-time data requirements were not available for all models. Prior to the analysis, 1979-2014 BCDS product was validated against two observational analysis datasets, the North American Regional Reanalysis (NARR, psl.noaa.gov/data, accessed on 1 March 2022) and TerraClimate [28,29].
Analysis of extreme fire weather consisted of four main components. First, we calculated the June, July, August (JJA) and the September, October, November (SON) seasonalmean FWI and VPD averages for each ecoregion, as well as for the entire state of California, for an historical (1981-2010), mid-century (2041-2070), and late-century (2071-2100) time period. The historical mean, calculated individually for each ecoregion and for each season, was used as a baseline with which the mid-century and late-century anomalies were calculated. In addition, the relative change (relative to the historical time period) of FWI and VPD in the mid-and late-century periods was calculated. We consider a unique baseline for each ecoregion, in addition to a California state-wide assessment, for each season, in order to properly account for the entire range of climate conditions found across California, and to assess projected changes in each zone relative to their particular local historical climate.
Second, the seasonal, ecoregion specific, daily FWI and VPD averages were used to calculate a 95th percentile value (FW I 95 ; VPD 95 ), which was then used as a threshold for the mapping of the frequency with which fire weather indices exceed that 95th percentile value during the historical, mid-and late-century periods. The maximum occurrence of a given index value above the 95th percentile is 91, or, in other words, all the days in that season in a given year. In order to evaluate the spatial distribution of extreme fire weather throughout the 21st century, the difference in the total number of days per year that exceed the 95th percentile was taken between the mid-and late-century and the historical period. In addition, the ecoregion average number of days per year that exceed the 95th percentile was calculated.
Third, 30 year JJA and SON seasonal average FWI and VPD were calculated for each ecoregion. Considering all grid points within each ecoregion for that 30 year average, the 95th percentiles are calculated and subsequently used as thresholds to determine the extent within each region that will experience FWI and VPD exceeding those values. Frequency distributions of 30 year spatial average FWI and VPD values were plotted for each ecoregion.
Finally, Julian-day averages spanning the entire 30 year periods were calculated for FWI and VPD. Daily climatologies intend to show changing peaks in average maximum values, as well as the lengthening or shortening of each ecoregion's fire season.

Validation of Historical Simulations
Before investigating the impacts of climate change on extreme fire weather, we evaluate the CMIP6 Raw and BCSD ensembles against an observational ensemble that includes NARR and TerraClimate. However, TerraClimate has not released a minimum relative humidity dataset, so in this case, CMIP6 was validated exclusively with NARR. Spatial statistics, including correlation, bias, and root-mean square error (e.g., sCor, sBias, sRMSE) were calculated after resampling all observations and Raw CMIP6 data to the high resolution BCSD (~4 km). In this way, the improved skill of the bias correction and downscaling was demonstrated against a basic resampling. Temporal statistics (tCor, tBias, tRMSE) were calculated using monthly averages spanning the entire historical period . The difference in resolution between the TerraClimate observations (~4 km), raw MIROC6  Spatial correlations between CMIP6 BCSD and the observations were high for all variables, especially for maximum temperature, precipitation accumulation, and minimum relative humidity, with 0.97, 0.88, and 0.87, respectively-all of which were higher relative to the CMIP6 Raw (Table 1). Furthermore, sRMSE and sBias were higher for the CMIP6 BCSD relative to the CMIP6 Raw for all variables except temperature, which showed a higher sBias in the CMIP6 Raw. A substantially higher sCor and sBias in the wind speed resulted between the Raw and the BCSD, with sCor increasing by 0.4 and sBias by 0.38 (Table 2).  Spatial correlations between CMIP6 BCSD and the observations were high for all variables, especially for maximum temperature, precipitation accumulation, and minimum relative humidity, with 0.97, 0.88, and 0.87, respectively-all of which were higher relative to the CMIP6 Raw (Table 1). Furthermore, sRMSE and sBias were higher for the CMIP6 BCSD relative to the CMIP6 Raw for all variables except temperature, which showed a higher sBias in the CMIP6 Raw. A substantially higher sCor and sBias in the wind speed resulted between the Raw and the BCSD, with sCor increasing by 0.4 and sBias by 0.38 (Table 2). tCor stayed the same for maximum temperature, daily precipitation accumulation and minimum relative humidity, but was higher for wind speed in CMIP6 BCSD relative to the CMIP6 Raw. Similar to the spatial statistics, maximum temperature's tBias was the only metric that was not lower in the CMIP6 BCSD.
In addition to the statistical metrics, a monthly climatology of each variable was generated for the CMIP6 Raw, CMIP6 BCSD, and the observational ensemble ( Figure 3). BCSD successfully shifted the summertime temperatures downwards towards the observations, however, temperatures were slightly overestimated in the winter months ( Figure 3a). Both the CMIP6 Raw and BCSD overestimated daily rainfall accumulation in the winter, however, the BCSD more closely approximated the observations throughout the year (Figure 3b). The NARR minimum relative humidity fell between the CMIP6 Raw and BCSD monthly climatology, however, even though BCSD slightly underestimated winter relative humidity, the difference was smaller in all months than it was for CMIP6 Raw ( Figure  3c). CMIP6 BCSD 10 m wind speed was very well simulated in California compared to CMIP6 Raw (Figure 3d).  The NARR minimum relative humidity fell between the CMIP6 Raw and BCSD monthly climatology, however, even though BCSD slightly underestimated winter relative humidity, the difference was smaller in all months than it was for CMIP6 Raw (Figure 3c). CMIP6 BCSD 10 m wind speed was very well simulated in California compared to CMIP6 Raw (Figure 3d).   Overall, JJA and SON anomalies showed increasing average FWI between 2041-2100 (except MWCM), however, relative change during the SON season was greater than JJA in all ecoregions except the Central Valley during the mid-century period. FWI in the MWCM decreased relative to the historical baseline in both seasons, however, decreases during the SON were of much smaller magnitude, and are basically a continuation of historical conditions. In both the mid-and late-century periods, the central and northern California C, EC, SN, and Central Valley ecoregions experienced greater relative change in FWI compared to the southern California ecoregions SCC and SCM.

Fire Weather Anomalies and Relative Change
JJA and SON vapor pressure deficit increased relative to an historical baseline in all ecoregions during both the mid-and late-century time periods ( Figure 5). Relative change of VPD during the fall was greater than in the summer in all ecoregions. JJA VPD during the mid-century period ranged from 14.96 in the SCC to 25.14 in the EC, and 22.08 in the SCC to 39 in the K during the late-century. On the other hand, mid-century SON vapor pressure deficit ranged from 15.68 in the SCC to 27.46 in the K, and from 26.16 in the SCC to 41.91 in the K ecoregion, during the late-century period. Overall, JJA and SON anomalies showed increasing average FWI between 2041-2100 (except MWCM), however, relative change during the SON season was greater than JJA in all ecoregions except the Central Valley during the mid-century period. FWI in the MWCM decreased relative to the historical baseline in both seasons, however, decreases during the SON were of much smaller magnitude, and are basically a continuation of historical conditions. In both the mid-and late-century periods, the central and northern California C, EC, SN, and Central Valley ecoregions experienced greater relative change in FWI compared to the southern California ecoregions SCC and SCM.
JJA and SON vapor pressure deficit increased relative to an historical baseline in all ecoregions during both the mid-and late-century time periods ( Figure 5). Relative change of VPD during the fall was greater than in the summer in all ecoregions. JJA VPD during the mid-century period ranged from 14.96 in the SCC to 25.14 in the EC, and 22.08 in the SCC to 39 in the K during the late-century. On the other hand, mid-century SON vapor pressure deficit ranged from 15.68 in the SCC to 27.46 in the K, and from 26.16 in the SCC to 41.91 in the K ecoregion, during the late-century period. Similar to FWI, JJA and SON relative change in VPD in the mid-and late-century periods increased more in the central and northern California ecoregions of the K, C, EC, SN, and Central Valley compared to the southern California ecoregions SCC and SCM ( Figure 5). While FWI decreased in the MWCM ecoregion relative to its historical baseline (Figure 4j), VPD increased in both the mid-and late-century periods (Figure 5j). Differences between JJA and SON seasonal changes in VPD were greater than for FWI, however lower average historical values during SON led to comparable relative change (to the JJA season) in the mid-and late-century periods. Overall, both JJA and SON vapor pressure deficit is expected to increase linearly throughout 2041-2100, in all ecoregions.

95th Percentile Exceedance Maps
To further investigate the impact of climate change on extreme fire weather, we calculated the number of days that FWI and VPD were greater than or equal to the 95th percentile (calculated using a time series of daily average values for each ecoregion, not per pixel) ( Figure 6). On average, JJA 95 days increased substantially in all ecoregions in both the mid-and late-century periods ( Table 2). The average number of days per year with VPD greater than the 95th percentile doubled in all ecoregions by 2100. In fact, from 2071-2100, over half the total number of summer days in the CCF, C, EC, CV, and NAD ecoregions had a VPD greater than the historical 95th percentile.  Figure 5). While FWI decreased in the MWCM ecoregion relative to its historical baseline (Figure 4j), VPD increased in both the mid-and late-century periods (Figure 5j). Differences between JJA and SON seasonal changes in VPD were greater than for FWI, however lower average historical values during SON led to comparable relative change (to the JJA season) in the mid-and late-century periods. Overall, both JJA and SON vapor pressure deficit is expected to increase linearly throughout 2041-2100, in all ecoregions.

95th Percentile Exceedance Maps
To further investigate the impact of climate change on extreme fire weather, we calculated the number of days that FWI and VPD were greater than or equal to the 95th percentile (calculated using a time series of daily average values for each ecoregion, not per pixel) ( Figure 6). On average, JJA VPD 95 days increased substantially in all ecoregions in both the mid-and late-century periods ( Table 2). The average number of days per year with VPD greater than the 95th percentile doubled in all ecoregions by 2100. In fact, from 2071-2100, over half the total number of summer days in the CCF, C, EC, CV, and NAD ecoregions had a VPD greater than the historical 95th percentile. The largest increases in 95 days occurred in the desert southwest region of California, particularly in the late-century period, where some areas are expected to experience historical 95 on nearly two-thirds of all summer days. Other areas expected to experience large increases in 95 days are the SCC ecoregion east of Los Angeles (34°N, 117°W), the southern portion of the Central Valley and the Central California foothills, the western border of the Sierra Nevada, and scattered areas of the northern California Klamath, Cascades, and Eastern Cascades regions (Figure 6a-c).
JJA 95 days increased in both the mid-and late-century in all ecoregions except the MWCM, which experienced an average decrease in 95 days (Table 2). While the increase in average number of days with 95 is smaller than that of 95 , the SCM, C, EC, and CV double the number of 95 days by 2100. In addition, the greatest magnitude of change of 95 days occurred between the historical and the mid-century periods; the total change between the mid-and late-century was much smaller ( Table 2). The most substantial increases in 95 days occurred in the desert southwest, the CV, the southern tip of the SN, the CCF, and the C and EC in northern California. However, there was a clear pattern of decreased numbers of 95 and 95 days (particularly FWI) in the mid-and late-century periods along the coast of California where the ocean moderates temperature. Furthermore, there was an apparent influence of topography in northern California (K, C, SN ecoregions), where high elevation areas showed little or even negative change, while, for example, the lower elevation, western slopes of the Sierra The largest increases in VPD 95 days occurred in the desert southwest region of California, particularly in the late-century period, where some areas are expected to experience historical VPD 95 on nearly two-thirds of all summer days. Other areas expected to experience large increases in VPD 95 days are the SCC ecoregion east of Los Angeles (34 • N, 117 • W), the southern portion of the Central Valley and the Central California foothills, the western border of the Sierra Nevada, and scattered areas of the northern California Klamath, Cascades, and Eastern Cascades regions (Figure 6a-c).
JJA FW I 95 days increased in both the mid-and late-century in all ecoregions except the MWCM, which experienced an average decrease in FW I 95 days ( Table 2). While the increase in average number of days with FW I 95 is smaller than that of VPD 95 , the SCM, C, EC, and CV double the number of FW I 95 days by 2100. In addition, the greatest magnitude of change of FW I 95 days occurred between the historical and the mid-century periods; the total change between the mid-and late-century was much smaller ( Table 2). The most substantial increases in FW I 95 days occurred in the desert southwest, the CV, the southern tip of the SN, the CCF, and the C and EC in northern California. However, there was a clear pattern of decreased numbers of FW I 95 and VPD 95 days (particularly FWI) in the mid-and late-century periods along the coast of California where the ocean moderates temperature. Furthermore, there was an apparent influence of topography in northern California (K, C, SN ecoregions), where high elevation areas showed little or even negative change, while, for example, the lower elevation, western slopes of the Sierra Nevada showed pronounced positive change. The dramatic increase in VPD 95 and FW I 95 days in the NAD ecoregion is an important climatic signal for human health and well being, as well as for agriculture, but these increases are not likely to translate into increases in annual area burned because of a near total lack of continuous and flammable vegetation.
Historical 95th percentile thresholds were also calculated for each ecoregion for the SON season. The average number of VPD 95 days increased in all ecoregions ( Table 2). In fact, the average number of VPD 95 days tripled by the end of the late-century period in the CCF, K, C, EC, and CV ecoregions, and doubled in the SCC, SN, MWCM, and NAD ecoregions (relative to the historical period). The average number of VPD 95 days in the late-century period was substantially lower than it was in JJA, however, these results still imply that VPD is expected to increase dramatically for the entire duration of the summer and fall in all ecoregions ( Table 2)  Nevada showed pronounced positive change. The dramatic increase in 95 and 95 days in the NAD ecoregion is an important climatic signal for human health and well being, as well as for agriculture, but these increases are not likely to translate into increases in annual area burned because of a near total lack of continuous and flammable vegetation.
Historical 95th percentile thresholds were also calculated for each ecoregion for the SON season. The average number of 95 days increased in all ecoregions ( Table 2). In fact, the average number of 95 days tripled by the end of the late-century period in the CCF, K, C, EC, and CV ecoregions, and doubled in the SCC, SN, MWCM, and NAD ecoregions (relative to the historical period). The average number of 95 days in the late-century period was substantially lower than it was in JJA, however, these results still imply that VPD is expected to increase dramatically for the entire duration of the summer and fall in all ecoregions ( Table 2) Table  2). The C, EC, and CV ecoregions doubled in the number of days experiencing FWI values  Table 2). The C, EC, and CV ecoregions doubled in the number of days experiencing FWI values greater than the historical 95th percentile. Overall changes in the occurrence of FW I 95 days were smaller for the SON season compared to JJA, however, many of the same areas that were impacted in the summer show continued increases in extreme fire weather carrying into the fall, relative to the fall historical 95th percentile (Figure 7f). Some of the most affected areas include the region east of Los Angeles (34 • N, 117 • W), the south western region of the CCF (35-36 • N, 120-121 • W), western border of the Sierra Nevada, and the northern California Cascades and Eastern Cascades (Figure 7). Similar to the summer season, there was a moderating influence on FWI along the coast of California, as well as a topographic influence in the K and SN ecoregions. In fact, nearly the entire state of California, apart from the coastal and high elevation areas, displayed a positive signal in VPD 95 by the end of the 21st century.

Fire Weather Frequency Distributions
The 30 year spatial averages of FWI and VPD were calculated for each of the 10 ecoregions during the historical, mid-and late-century periods. The frequency with which each FWI and VPD value occurred was plotted (as an area), along with the ecoregion specific 95th percentile value and the area that exceeded the threshold value in the three time periods (Figure 8). greater than the historical 95th percentile. Overall changes in the occurrence of 95 days were smaller for the SON season compared to JJA, however, many of the same areas that were impacted in the summer show continued increases in extreme fire weather carrying into the fall, relative to the fall historical 95th percentile (Figure 7f). Some of the most affected areas include the region east of Los Angeles (34°N, 117°W), the south western region of the CCF (35-36°N, 120-121°W), western border of the Sierra Nevada, and the northern California Cascades and Eastern Cascades (Figure 7). Similar to the summer season, there was a moderating influence on FWI along the coast of California, as well as a topographic influence in the K and SN ecoregions. In fact, nearly the entire state of California, apart from the coastal and high elevation areas, displayed a positive signal in 95 by the end of the 21st century.

Fire Weather Frequency Distributions
The 30 year spatial averages of FWI and VPD were calculated for each of the 10 ecoregions during the historical, mid-and late-century periods. The frequency with which each FWI and VPD value occurred was plotted (as an area), along with the ecoregion specific 95th percentile value and the area that exceeded the threshold value in the three time periods ( Figure 8). The total land area that was projected to experience FWI greater than the 95th percentile value increased in all ecoregions apart from Klamath and the Marine West Coast Mountains (Figure 8e,j). The largest increases in area occurred in the Central Valley and the Cascades ecoregions (over seven and six times the historical area, respectively). The total land area that was projected to experience FWI greater than the 95th percentile value increased in all ecoregions apart from Klamath and the Marine West Coast Mountains (Figure 8e,j). The largest increases in area occurred in the Central Valley and the Cascades ecoregions (over seven and six times the historical area, respectively). Furthermore, the increase in area experiencing extreme fire weather in the SCM and the EC, followed by the SCC, CCF, and NAD ecoregions, increased by a factor of four and three, respectively. In many cases, including in the SCM, K, C, EC, SN, CV, and NAD ecoregions, the distribution of average FWI values shifted towards the 95th percentile threshold, indicating that midand late-century average FWI values could be closer to what is now considered extreme. The JJA FW I 95 percentile values ranged from 15.5 in the Cascades to 28.4 in the North American Desert, while the state-wide California JJA FW I 95 was 26.5.
The increase in the amount of land area likely to experience FWI greater than the historical 95th percentile was even greater for SON than JJA, however, the MWCM ecoregion is projected to have a decline in the total land area with extreme fire weather ( Table 3). The ecoregions projected to see the largest changes include the Central Valley, Cascades, and SCM, with ten, eight, and seven times the area under 95th percentile FWI by late-century, relative to the historical period. In addition, the NAD, EC, and SCC ecoregions increased by a factor of six, and the SN by a factor of five. Whereas the area within the Klamath ecoregion with FWI greater than the 95th percentile decreased in the summer, FW I 95 area is projected to double in SON by the end of the late-century. Additionally, the SON FW I 95 percentile values ranged from 8.71 in the Cascades to 16.6 in the North American Desert, while the state-wide California SON FW I 95 was 15.4. Ecoregion increases in JJA VPD 95 area were even more substantial than FW I 95 (Figure 9). The increase in VPD is so significant that the distribution of JJA VPD values shifts in a way that the peak VPD value in the mid-and late-century will exceed that of the historical 95th percentile threshold (Figure 9b-d,f,g,i,k). The largest increases in area exceeding VPD 95 occurred in the Central Valley and East Cascades, with a 17-fold increase, followed by the Cascades (13-fold), and the CCF and SCC (12-fold). The JJA VPD 95 percentile values ranged from 3.12 in the MWCM to 6.98 in the North American Desert, while the state-wide California JJA VPD 95 was 6.6.
The average projected area with extreme SON VPD was also large, implying that conditions in California will be conducive to wildfire for a minimum of half the year, extending the typical fire season into the fall, and contributing to the widespread desiccation of fuels. Unlike with patterns seen in the MWCM FWI values, VPD 95 area was projected to increase in the MWCM in both the JJA and SON by a factor of seven and nine, respectively ( Figure 9j and Table 2). Additionally, the SON VPD 95 percentile values ranged from 2.01 in the East Cascades to 4.28 in the North American Desert, while the state-wide California JJA VPD 95 was 3.97. The average projected area with extreme SON VPD was also large, implying that conditions in California will be conducive to wildfire for a minimum of half the year, extending the typical fire season into the fall, and contributing to the widespread desiccation of fuels. Unlike with patterns seen in the MWCM FWI values, 95 area was projected to increase in the MWCM in both the JJA and SON by a factor of seven and nine, respectively ( Figure 9j and Table 2). Additionally, the SON 95 percentile values ranged from 2.01 in the East Cascades to 4.28 in the North American Desert, while the state-wide California JJA 95 was 3.97.

Fire Weather Index Julian Day Climatologies
Julian day climatologies of VPD indicated an increase in peak VPD in the summer, as well as an earlier start and a delayed end to the fire season in all ecoregions of California ( Figure 10). For example, by 2100, the average onset of high spring VPD (defined here as VPD ≥ 3) is a month or more earlier in the SCC, CCF, K, C, SN ecoregions. Furthermore, the decrease of VPD below the critical value is delayed by a month or more in the C, EC and SN. Our results indicated that climate change will likely increase summer and autumn average VPD and lengthen the amount of time each year that conditions favor large wildfires in California ( Figure 10).

Fire Weather Index Julian Day Climatologies
Julian day climatologies of VPD indicated an increase in peak VPD in the summer, as well as an earlier start and a delayed end to the fire season in all ecoregions of California ( Figure 10). For example, by 2100, the average onset of high spring VPD (defined here as VPD ≥ 3) is a month or more earlier in the SCC, CCF, K, C, SN ecoregions. Furthermore, the decrease of VPD below the critical value is delayed by a month or more in the C, EC and SN. Our results indicated that climate change will likely increase summer and autumn average VPD and lengthen the amount of time each year that conditions favor large wildfires in California ( Figure 10).
Julian day climatologies of FWI indicated that the average summer peak FWI value will increase, and the fire season, defined here by the onset of average FWI greater than 15, will begin earlier and end later in the mid-and late-century periods relative to the historic period. However, the effect of climate change on the growth of the fire season was greater for VPD than it was for FWI (Figure 11). The greatest impact was to the onset of the fire season, with late-century average FWI exceeding a value of 15 anywhere from 6 days (SCC) to one month (NAD) before the historical average. Julian day climatologies of FWI indicated that the average summer peak FWI value will increase, and the fire season, defined here by the onset of average FWI greater than 15, will begin earlier and end later in the mid-and late-century periods relative to the historic period. However, the effect of climate change on the growth of the fire season was greater for VPD than it was for FWI ( Figure 11). The greatest impact was to the onset of the fire season, with late-century average FWI exceeding a value of 15 anywhere from 6 days (SCC) to one month (NAD) before the historical average.  Studies show that the role of temperature and humidity are strongest in FWI, with precipitation and wind speed playing important, but secondary roles [30,31]. However, our results indicated that the increase of VPD in the mid-and late-century periods was much greater than FWI, especially in the MWCM, which experienced a decrease in summer FWI, implying that precipitation and wind patterns can play a strong role in FWI trajectory (Figure 11j). However, the sensitivity of the FWI to incremental changes in input variables is outside the scope of this study.  Studies show that the role of temperature and humidity are strongest in FWI, with precipitation and wind speed playing important, but secondary roles [30,31]. However, our results indicated that the increase of VPD in the mid-and late-century periods was much greater than FWI, especially in the MWCM, which experienced a decrease in summer FWI, implying that precipitation and wind patterns can play a strong role in FWI trajectory (Figure 11j). However, the sensitivity of the FWI to incremental changes in input variables is outside the scope of this study. California's fire weather season length has increased in the past three decades due to increases in surface air temperature, decreases in relative humidity, and longer annual rain-free periods [30,32]. Our results provide further evidence that fire weather conditions conducive to large wildfires will become more extreme and extend later into the fall season in most areas of California by 2100 [3,30]. The peak effects of hot and dry summers will be extended by 30 days into late October and November. This timing will coincide with the Santa Ana wind season (which normally runs from October to March), as well as a projected delay in the onset of the rainy season-producing a high risk period during a time that has normally been outside of California's typical "fire season" for emergency preparedness [30,33]. The catastrophic consequences of severe fire weather towards the end of the typical California fire season were seen in 2020, with the August Complex Fire (August 2020), SCU Lightning Complex (August 2020), Creek (September 2020), LNU Lightning Complex (August 2020), and the North Complex (August 2020)-together accounting for five of the 10 largest wildfires in California history [34,35]. In addition, the increase in land area experiencing extreme VPD and FWI expected by the end of the mid-century (2070) and late-century (2100) period, as well as the extension of the fire season late into the fall, will prime large portions of California for the spread of large wildfires, given an ignition (Figures 10 and 11).

Climate Change Impacts on California Ecosystems
Several conceptual frameworks exist that seek to determine the biophysical controls on wildfire ignition and probability of spread. For example, wildfire has been shown to thrive in the middle of a continuum between aridity (dry) and productivity (wet) [36]. Similar to the way that plants have evolved traits that reflect biotic tradeoffs in favor of various survival strategies, wildfire regime is commonly controlled by the abundance and moisture content of fuels [37][38][39]. For example, the tropics experience large amounts of rainfall, have high biomass production, and contain the necessary fuels for large wildfireshowever, they are rarely dry enough to burn [38,39]. On the other hand, where moisture is limited and biomass production is low, wildfire does not have the necessary fuel to spread. Therefore, wildfire would be most likely to occur and spread in ecosystems that are dry for sufficient periods of time so that vegetation becomes water-stressed, but wet enough to support enough biomass that provide continuous fuels such that large wildfires may burn. These environmental gradients impact all aspects of fire regime, including fire frequency, occurrence, intensity, severity, and total area burned [8,38,39]. This phenomenon implies that the projected impacts of climate change on temperature and atmospheric aridity will have a more significant effect on northern California forested ecoregions sensitive to extended dry periods than on the more consistently hot and dry southern California regions that are adapted to such conditions, and are more likely to lack the necessary abundance and continuity of fuel for large wildfires [9]. These climate-fire relationships are critical for the interpretation of our results, as the most substantial changes to the amount of area experiencing extreme fire weather (Figures 8 and 9), as well as the time added to the onset and end of the fire season (Figures 10 and 11), are projected to occur in forested northern California, where they could have the most dire implications for annual area burned throughout the 21st century.
In southern California, a large variety of topographies, vegetation types, climate regimes, and human populations contribute to a fire regime dominated by frequent wildfires. The precipitation season is confined to a few core months (Dec.-Mar.) with some rain occurring in the shoulder season of November and April, which in combination with the warm and dry summers, leads to an extended annual fire season [32]. In addition, annual Santa Ana winds advected from the eastern deserts towards the coast, move downslope at high speeds and increase the risk of severe wildfire during the fall season. Vegetation types are predominately herbaceous plants (chaparral), grasslands, shrubs, and patches of forest that recover relatively quickly following wildfire. The rapid regeneration of herbaceous fuels after wildfire are a factor in the frequent fire return interval throughout southern California (1-5 years) [40]. Studies show that the impacts of climate change on fire regime in southern California thus far are primarily a result of increasing temperatures and atmospheric aridity, as well as changing patterns of precipitation, that in turn influence the state of the environmental constraints on wildfire activity discussed above [41]. Compared to forested regions of northern California, fire regime in non-forested regions of southern California are less susceptible to shifts in moisture availability and more so to availability of fuel, however, these chaparral and grassland ecosystems are sensitive to rainfall patterns in the year or two prior to fire occurrence because the influx of moisture promotes the growth of vegetation [1,42]. Our results indicate that while the increases in land area experiencing extreme fire weather are smaller than in northern California, and the extension of the fire season is not as substantial, changes to southern California JJA and SON extreme fire weather, in succession, still have the potential to tax vegetation and fire-fighting resources.
Conceptual frameworks help to attribute controlling factors to different fire regimes, given the complex web of climate, human activity, topography, vegetation, hydrology, prior wildfire, fire suppression, and fire management, with their practically infinite combination. Our ecoregion specific results show a future consisting of more persistent extreme fire weather that occurs in a greater land surface area than in the past. Therefore, with the aid of climate-fire frameworks, we may deduce the range of impacts that this increase in aridity and fire potential may have on sensitive ecosystems. While northern and southern California ecoregion's fire regimes are controlled by different factors, the simultaneous increasing trends in atmospheric aridity and fire potential indicate a future in which fire suppression and management resources throughout the state are stretched thin. In order to offset the dangers of extended fire seasons throughout much of California by 2100, alongside the expansion of the wildland urban interface, fire management must implement innovative combinations of new and existing strategies intended to limit the destruction of property and loss of life.

Limitations
Impact studies that utilize global climate model output generally require the use of many models, however, due to the recent release of CMIP6, limited availability of the fire weather variables required to perform the FWI and VPD calculations restricted our ability to increase the number of models used [43]. Despite this, the evaluation of the ensemble-mean of the three models chosen for this study showed good agreement with observational data. In addition, we acknowledge that all statistical downscaling methods are subject to issues of stationarity. To be clear, there are assumptions made regarding the application of historical climatic patterns to projected data, however, the thin plate spline and linear regression method employed here advantageously utilizes data surrounding each grid point in the calculation of its smoothed estimate. This interpolation technique works best with continuous data but it may be adapted to perform well with discontinuous climate variables such as precipitation, and to some extent wind speed. Additionally, it is important to consider the climate-wildfire-vegetation feedback when attempting to predict annual area burned, however, this work focuses on the impacts of climate on extreme fire weather, so the complex relationships among annual area burned, climate patterns, and vegetation distribution are outside the scope of this study [44]. Lastly, it should be noted that projected increases in California's summer and fall VPD and FWI are not synonymous with a corresponding increase in annual area burned, as wildfire is inherently dependent on an ignition source and flammable vegetation [4].

Conclusions
We analyzed the impacts of climate change on extreme fire weather in ten different California ecoregions throughout historical, mid-and late-century periods using two prominent fire weather indices. Climate model output from three CMIP6 models were bias corrected and statistically downscaled using empirical quantile mapping in addition to a simple thin plate spline interpolation method. Historical average baselines for each index were calculated for each ecoregion and used to generate anomalies for 2041-2100. Relative change analysis indicates that FWI and VPD are projected to increase in all ecoregions, except the MWCM, which displayed a decrease in JJA and SON FWI in the mid-and late-century periods. In addition, the average daily 95th percentile value for the historical period (for each ecoregion) was implemented as a threshold in order to calculate the number of summer and fall extreme fire weather days, and to evaluate the spatial distribution of these changes. The frequency of 95th percentile days for JJA and SON VPD and FWI is projected to increase significantly in all ecoregions, except the MWCM, which shows a decline in FWI 95th percentile days. Furthermore, we analyzed the 30 year spatial average of each index and generated frequency distributions to investigate the average change in land area projected to experience extreme fire weather, relative to the historical benchmark. We found that JJA VPD land area exceeding the 95th percentile will increase by a factor of 2 and 17 in different ecoregions by the end of 2100. Lastly, Julian-day climatologies of each index were generated to assess the changing duration of a typical California fire season. We found that the window of extreme fire weather will expand into the spring and late autumn in all ecoregions (apart from the MWCM), extending the amount of time vegetation is exposed to increased atmospheric demand for moisture, and heightening the overall risk for the ignition and spread of large wildfire.
The ecoregion-level spatial scale adopted for this study increases the amount of local information, as well as the resolution with which fire and land managers can implement strategies and counter-measures when addressing issues related to climate change. The spatial downscaling algorithm tested and implemented in this study was computationally inexpensive (compared to dynamical downscaling) and may be readily applied to further climate change scenarios in other Mediterranean and temperate regions [13,15]. High resolution, daily meteorological data is not only useful for wildfire impact studies (and for investigating the post-fire environment), but can be implemented in the research of agricultural productivity, urban heat islands, hydrology, and human health. This intuitive, practical, and simple interpolation method will not only contribute to the improved understanding of the complexities of California's climate-fire dynamics, but will provide the means for more efficiently obtaining climate model data of high spatial and temporal resolution.