Modeling the effects of spatially explicit patterns of climate and fire on future populations of a fire-dependent plant

Many plant species are likely to face population decline or even extinction in the coming century, especially those with a limited distribution and inadequate dispersal relative to the projected rates of climate change. The obligate seeding California endemic, Ceanothus perplexans is especially at risk, and depending on how climate change interacts with altered fire regimes in Southern California, certain populations are likely to be more at risk than others. To identify which areas within the species’ range might need conservation intervention, we modeled population dynamics of C. perplexans under various climate and fire regime change scenarios, focusing on spatially explicit patterns in fire frequency. We used a species distribution model to predict the initial range and potential future habitat, while adapting a density-dependent, stage-structured population model to simulate population dynamics. As a fire-adapted obligate seeder, simulated fire events caused C. perplexans seeds to germinate, but also killed all adults in the population. Our simulations showed that the total population would likely decline under any combination of climate change and fire scenario, with the species faring best at an intermediate fire return interval of around 30–50 years. Nevertheless, while the total population declines least with a 30–50 year fire return interval, the effect of individual subpopulations varies depending on spatially explicit patterns in fire simulations. Though climate change is a greater threat to most subpopulations, increased fire frequencies particularly threatened populations in the northwest of the species’ range closest to human development. Subpopulations in the mountainous southern end of the range are likely to face the sharpest declines regardless of fire. Through a combination of species distribution modeling, fire modeling, and spatially explicit demographic simulations, we can better prepare for targeted conservation management of vulnerable species affected by global change.

Many plant species are likely to face population decline or even extinction in the coming century, especially those with a limited distribution and inadequate dispersal relative to the projected rates of climate change. The obligate seeding California endemic, Ceanothus perplexans is especially at risk, and depending on how climate change interacts with altered fire regimes in Southern California, certain populations are likely to be more at risk than others. To identify which areas within the species' range might need conservation intervention, we modeled population dynamics of C. perplexans under various climate and fire regime change scenarios, focusing on spatially explicit patterns in fire frequency. We used a species distribution model to predict the initial range and potential future habitat, while adapting a density-dependent, stage-structured population model to simulate population dynamics. As a fire-adapted obligate seeder, simulated fire events caused C. perplexans seeds to germinate, but also killed all adults in the population. Our simulations showed that the total population would likely decline under any combination of climate change and fire scenario, with the species faring best at an intermediate fire return interval of around 30-50 years. Nevertheless, while the total population declines least with a 30-50 year fire return interval, the effect of individual subpopulations varies depending on spatially explicit patterns in fire simulations. Though climate change is a greater threat to most subpopulations, increased fire frequencies particularly threatened populations in the northwest of the species' range closest to human development. Subpopulations in the mountainous southern end of the range are likely to face the sharpest declines regardless of fire. Through a combination of species distribution modeling, fire modeling, and spatially explicit demographic simulations, we can better prepare for targeted conservation management of vulnerable species affected by global change. KEYWORDS fire, Ceanothus, chaparral, climate change, species range, population model, population dynamics

Introduction
Many plant species are likely to face population decline or extinction in the coming century due to global change, including climate change, ongoing land-use change, invasive species, and altered disturbance regimes such as fires and floods (Komatsu et al., 2019;Driscoll et al., 2021). Species with limited distributions, low rates of dispersal relative to projected rates of habitat shifts due to climate change (Corlett and Westcott, 2013), and those lacking physiological adaptations to disturbance (Keith, 2012;Miller et al., 2019) will be the most vulnerable to these threats. Changes in precipitation and temperature due to climate change are expected to impact the climate water deficit (Flint et al., 2013;Mitchell et al., 2016) -the quantified evaporative demand that exceeds available soil moisture for plants (Stephenson, 1998). In response, suitable habitat for many plant species will shift, fragment, reduce, and in some cases expand (Lutz et al., 2010;Restaino et al., 2016). In addition to habitat changes, climate and land use change may also exacerbate changes in some episodic environmental processes under which many plants have evolved, e.g., fire regimes, hydrological regimes, and storm events (Johnstone et al., 2016). In turn, changes in disturbance regimes such as fire could rival climate and land-use change as primary biodiversity threats (Turner, 2010;Bowman et al., 2011), with significant consequences for species in fire-prone ecosystems (Syphard et al., 2009;Keeley et al., 2012).
Of particular concern are Mediterranean-type ecosystems (MTEs) which world-wide are among the highest in plant species richness and endemism, including many species with restricted ranges (Cowling et al., 1996). MTEs are characterized by cool, wet winters and warm dry summers and are found adjacent to cold ocean currents on the west coasts of continents at mid-latitudes. All five Mediterraneanclimate regions-Eurasian Mediterranean, southeastern Australia, South African Cape, Chile, and California-are biodiversity hotspots of global significance (Myers et al., 2000). Changes in climate and fire in MTEs are driven by spatially explicit processes, show spatial dependence, and occur in specific geographic locations and across characteristic landscape configurations (Franklin et al., 2021). Biodiversity in MTEs is particularly vulnerable to changes in climate and fire regime Doblas-Miranda et al., 2015); climate change may also increase fire risk in MTEs (Pausas et al., 2004;Keeley and Syphard, 2016).
The California Floristic Province (CFP) is a global MTE and biodiversity hotspot (Baldwin, 2014), with 20% of all vascular plant species found in the United States, 30% of which are endemic (Master et al., 2000). Species distributions among key plant functional types in the CFP are closely tied to climate, especially winter minimum temperatures (Ewers et al., 2003) and drought severity (Davis et al., 1999) which is closely related to temperature. These climate variables have already changed in the 20th century and are projected to change more in the future (Hayhoe et al., 2004;Cayan et al., 2008;Sun et al., 2015). It is possible that 20th century climate change has already had a measurable impact on plant species distributions in this region (Rapacciuolo et al., 2014).
California is also one of the most fire-prone regions in the world due to its unique combinations of climate, vegetation types, land use patterns, and fire management practices (Kelly et al., 2020;Syphard and Keeley, 2020). The state hosts a diversity of natural fire regimes (i.e., the characteristics of wildfires that occur in an area over long time periods; Syphard and Keeley, 2020), most of which have been altered over the past century due to a complex combination of these factors, which are changing independently and synergistically in different ways across the state (Keeley and Syphard, 2018;Williams et al., 2019;Syphard and Keeley, 2020;Li and Banerjee, 2021). The primary cause of fire ignitions in California is anthropogenic, although natural lightning-caused fires do occur in higher elevations, largely in the eastern part of the state (Keeley and Syphard, 2019). The highest fire frequency tends to occur at intermediate levels of development and population density, near the wildland-urban interface (Syphard et al., 2007;Syphard and Keeley, 2020;Jenerette et al., 2022). However, area burned, and the incidence of large fires tends to occur in more remote locations, with fire behavior governed by biophysical factors, such as topography, climate, and fuel conditions (Keeley and Syphard, 2021;Li and Banerjee, 2021). Changes in the timing and location of ignitions, along with increases in the number of large fires, are occurring to various extents across the state, often in different directions (Keeley and Syphard, 2021). The complex and important consequence of these changes to biodiversity is that plants and animals have evolved with and adapted to the long-term historical fire regimes in their range (Keith, 2012;Miller et al., 2019). Changes to these fire regimes, in both magnitude and direction, are likely to result in disturbances that exceed species' various adaptations.
Several studies have shown important cumulative and spatially heterogeneous effects of multiple threats on plant species persistence in MTEs (Keith et al., 2008;Conlisk et al., 2012;Regan et al., 2012;Bonebrake et al., 2014). These studies indicate that increased fire frequency can be more detrimental to some plant functional typesthose "obligate seeders" dependent on fire to break seed dormancy and a certain minimum time between fires for seed bank developmentthan habitat contractions due to climate change or urban growth . Threats and species traits could present complex interactions; it was shown that for an obligate seeding shrub, habitat fragmentation increased population persistence under frequent fires despite habitat loss because it decoupled fires across the landscape, reducing the size of individual fires (Regan et al., 2010). But this result held only up to a point; once habitat loss surpassed a threshold, the negative effects of reduced habitat through fragmentation outweighed the beneficial effects of risk spreading of frequent fire (Regan et al., 2010). Spatial location and context play a large role in the risk of decline of fire-dependent plants (Syphard et al., 2013;Bonebrake et al., 2014). For instance, for two obligate seeding species, under a historical fire regime, the ranking of threats, in terms of greatest to least reductions in population abundance, depended both on the driver of habitat change (urban growth or climate change) and where the species were located (range characteristics; Syphard et al., 2013). Nevertheless, for both species studied, high fire frequency caused populations to decline to a greater extent than projected habitat change under climate change scenarios, highlighting the complex interaction of plant functional type, spatial context, and fire frequency.
Modeling and mapping species' range shifts or reductions in suitable habitat under land use or climate change projections have been the major focus of investigations into the impacts of global change on biodiversity (e.g., Warren et al., 2018). They are necessary components of identifying climate and fire refugia for vulnerable species and for managing landscapes to ensure an appropriate amount, pattern, and timing of fire to promote species persistence (Kelly et al., 2020). However, in models where the effects of fire on plant population Frontiers in Ecology and Evolution 03 frontiersin.org persistence are explicitly considered, most studies assume a consistent fire regime across the study extent and do not spatially vary the average fire frequency in ways that account for the effects of fireregime drivers, such as proximity to human populations, spatially explicit climate, or topography. Studies that do alter the timing of fire to simulate temporal changes in fire frequency do not vary this aspect of the fire regime across a landscape (e.g., Conlisk et al., 2012). Hence, to provide more meaningful projections of plant population trends to support management decision making, it is necessary to give greater consideration to how fire regimes vary spatially.
In this paper, we use an integrated modeling framework to link climate induced changes in species distributions with population dynamics and spatially explicit fire frequency for Ceanothus perplexans, a native plant species from an iconic functional type across California and most MTEs: obligate seeding woody plants. In obligate seeding plants, soil-or canopy-stored seeds require the smoke and heat from fire to germinate, an evolutionary adaptation to fire regimes which requires a minimum number of fire free years to accrue a seed bank to maintain population viability . Ceanothus perplexans is an ideal case study on which gain insights on the effects of spatially explicit fire regimes on plant biodiversity in MTEs because it is a characteristic species in chaparral shrublands, a significant plant community in California. While all MTEs support different iconic plant communities, they share the defining characteristics of sclerophyllous, evergreen, deep-rooted shrubs, an outcome of convergent evolution (Rundel, 2018). Ceanothus perplexans also occurs in a landscape of intense urban growth and other land-use changes, its life history is dependent on fire, it is a poor disperser, and population reductions due to climate change are likely in the next century requiring conservation intervention (Schwilk and Keeley, 2012). We link downscaled climate projections for California and species distribution models with an established population model and historical fire records to investigate how including spatially explicit fire frequency compares and contrasts with spatially uniform fire frequency across the range of the species. We ask three main questions: i. Are there areas of the species range that are more (or less) affected by fire and climate change? ii. How do the population-level effects of spatially explicit fire frequency compare with uniform fire frequency? iii. Are some areas of the species range more likely to be appropriate as climate and/or fire refugia for C. perplexans than others?

Study species
Ceanothus perplexans Trel. (previously Ceanothus greggii var. perplexans Jeps.), Cupped leaf ceanothus, is an obligate seeding evergreen shrub primarily distributed in the southwestern region of the California Floristic Province (CFP; Hickman, 1993). Fire events generally kill all individual C. perplexans plants while stimulating the germination of seeds in the seed bank. Because germination rarely occurs in the absence of fire for this obligate seeder, C. perplexans are typically found in even-age cohorts (Keeley, 1992;Zammit and Zedler, 1993). Moreover, because this species has very limited dispersal ability, it relies largely on an extant seed bank to recover from fire events (Holl et al., 2000). Fire typically kills all seeds that do not germinate, and seeds should reduce to about 5% of the original value after 20 years with no additional input (Quinn, 1994;Regan et al., 2010). While the maximum lifespan of this species is uncertain, their survival and fecundity vary with age, typically declining above 80 years of age (Zammit and Zedler, 1993).

Current and future habitat distribution
We used species distribution models (SDMs) to estimate habitat suitability maps of C. perplexans at 270 × 270 m spatial resolution, restricted to the CFP ( Figure 1) for current (1981-2010 climatic averages) and future climate conditions. These SDMs were developed as part of a larger project assessing the distributions of 106 plant species in the California Floristic Province under climate change (Rose et al., 2023). Details of the SDM workflow are documented in that paper and in Supplementary Table S1 (following the ODMAP protocol; Zurell et al., 2020), and are outlined here. We created these SDMs using presence-absence data from vegetation surveys between 1980 and 2019 (provided by California Department of Fish and Wildlife Biogeographic Information and Observation System and Dr. James Thorne, University of California, Davis) and a set of predictor variables (Supplementary Table S1) related to climatic, hydrologic, terrain, and soil factors associated with plant distributions in waterlimited ecosystems like California (Stephenson, 1998;Kueppers et al., 2005;Dubuis et al., 2013). We used these initial habitat suitability maps (where the suitability of each cell is a value between 0 and 1) to represent the full extent of the species' distribution in 1995.
For future species distributions, we used climate projections from the Basin Charaterization Model (Flint et al., 2013) driven by each combination of two representative concentration pathways of greenhouse gas emissions scenarios (RCPs) under two different climate models. We chose these four projections to represent a range of emissions outcomes and climate models, with RCP 4.5 being moderate emissions and RCP 8.5 being high emissions, and a range of climatic conditions in California, with the CNRM-CM5 (Centre National de Recherches Meteorologiques climate model version 5) climate model projecting relatively cooler and wetter and HadGEM2-ES (Hadley Centre Global Environmental Model version 2-Earth system configuration) relatively hotter and drier conditions for California (Pierce et al., 2018). With these projections, we used the same climatic predictor variables from our initial habitat suitability maps, but projected to the years 2055 (midcentury, based on 2040-2069 climate normals) and 2085 (end-of-century, averaged between 2070 and 2099). We interpolated habitat suitability for years between 1995 and 2055 and years between 2055 and 2085. Because the species has negligible dispersal ability, we restricted future suitability maps to the extent of the initial distribution.
The distribution of C. perplexans notably extends slightly into Nevada, Arizona, and Baja California, Mexico, beyond the full extent of our model (Hickman, 1993). We did not include these areas at the edge of its range for lack of comparable high quality environmental and species occurrence data. By excluding this portion of the range of C. perplexans, we potentially introduce bias to our analysis. In particular, if temperatures in in those areas are Frontiers in Ecology and Evolution 04 frontiersin.org warmer than the rest of the range, then the species could potentially tolerate higher temperatures than our SDMs estimate and our future projections might overestimate habitat loss on the trailing edge of the distribution. However, the southernmost areas where C. perplexans is found, in the mountains of northern Baja California, experience a Mediterranean climate very similar to southern California. Including them is unlikely to introduce substantially new climate space to the SDM as was demonstrated in Regan et al. (2012). Overall, these limitations in our SDMs predicting future habitat suitability should not affect our focus on the effects of fire on population dynamics within the core of the species range.

Population model
We simulated population dynamics of C. perplexans using stagebased matrix models adapted from published demographic models (Regan et al., 2010, Syphard et al., 2013 Figure 1). These studies provided survival, transition, and fecundity parameters. Generally, Simulations of Ceanothus perplexans were a combination of a species distribution model, spatially explicit fire patterns, and a matrix-based population model. (A) To build the spatial patch structure of the model, the initial and future suitable habitats were built from a species distribution model. Future climate models came from two general circulation models and two emissions scenarios. We split suitable habitat into smaller patches where we simulated subpopulation dynamics. (B) Using historical records of fires in California, we estimated a spatially explicit map of mean fire return intervals. We divided this into patches the size of simulated subpopulations and used this as a baseline for spatially explicit fire heterogeneity and modified it to reflect changing fire regimes. (C) We simulated the population dynamics in each patch stage-structured demographic model. Seeds only germinated after fire events. (D) Combining these together, we compared final populations among several combinations of climate and fire models. each population follows the stage-based matrix model during each time step, with demographic stochasticity around the mean survival and fecundity.
Originally an age-structured matrix model with 2-year time steps, we converted the models into a stage-structured model with 2-year time steps. Before age 6, C. perplexans had lower survival and no fecundity compared with older age classes. To retain this dynamic, we preserved the age structure of all ages up to age 6. After this, we pooled adult C. perplexans survival and fecundity into stages. In the model, these stages were ages 6-12, 13-31, 32-56, 57-81, and 82+. Survival for all ages in each stage class was the same within the original model, so we used these survival probabilities for this model. Fecundity changed with each age class in the original model based on a polynomial fit of seed production (Syphard et al., 2013). We used fecundity at the lowest age of each stage. We determined transition probabilities between stages such that the mean time between stage classes was equal to the number of years between the age classes that these stage represented. For example, there were 7 years between the first adult stage (ages 6-12) and the second adult stage (ages 13-31), so the probability of a first stage adult surviving but remaining in this stage was 5 7 5 s and the probability of surviving but transitioning to second stage adult was 2 7 5 s (recalling that each time step is 2 years). Overall, the 2-year time step stage structure model for these species featured 9 stage classes: seeds, germinants, seedlings aged 2-3, seedlings aged 4-5, adults in stage 5 (equivalent to ages 6-12), adults in stage 6 (equivalent to ages 13-31), adults in stage 7 (equivalent to ages 32-56), adults in stage 8 (equivalent to ages 57-81), and adults in stage 9 (equivalent to ages 82+). The base transition matrix for C. perplexans was  Table 1. The standard deviation matrix that determined variation around the mean transition probabilities was set to 0 2 . g ij (where g ij is the transition from stage i to stage j) for 0.5 ij g £ and 0 2 1 . -( ) g ij for g ij > 0 5 . . For fecundity rates, the standard deviation was set to 2 f i .
Each time step in which a fire occurs, all C. perplexans plants in the population die. However, as obligate seeding species, a fraction of the seeds, r = 0 017167 . , germinate after a fire and survive to the next 2-year time step.
Though we analyzed populations by relative proportional change, carrying capacities and size classes were based on the maximum of 150 largest size adult-equivalent plants per hectare for C. perplexans (Syphard et al., 2013). The relative weight of each size class to the population density were based on age-specific data on canopy area from Zammit and Zedler (1993). We used a best fit function from Syphard et al. (2013), such that the larger individuals contributed more to population density and the largest individuals contributed 100 times more density than germinants, and seeds did not contribute to population density. Overall, the relative contribution to density was 0, 0.01, 0.1567, 0.2945, 0.3928, 0.5141, 0.7710, 0.9779, and 1 for size classes of increasing size, respectively.
To determine the carrying capacity of each 270 × 270 m grid cell over time, we estimated the total habitat suitability of the area from 1995 to 2085. As habitat suitability changed over time, so did the carrying capacity of the simulated patches. We grouped cells into larger patches, defined as distinct clusters of 270 × 270 m grid cells were within a two-cell radius of each other. However, using these criteria, the entire metapopulation consisted of several smaller patches and a small number of very large patches. Because the smallest patches required substantial computational resources,

Symbol
Value Description Frontiers in Ecology and Evolution 06 frontiersin.org and they were likely to have little effect on overall metapopulation results, we excluded patches where the total habitat suitability of all cells in the patch (calculated in the SDM) was <5 (meaning that patches were a minimum of 5 cells, or 0.135 hectares). These large populations masked a great deal of information on the spatial heterogeneity of climate change. Additionally, wildfires necessarily occur on the scale of an entire patch in this model, so ecoregionspanning fires would be too large to be realistic. Because of this, we placed a lattice across California, separating the entire metapopulation further. This lattice artificially created a series of vertical and horizontal lines of unsuitable habitat, 4 cells wide with 30 cells of potentially suitable in between. This created a maximum of 30 × 30 grid cell patches (6,561 ha). While patches at this scale are still larger than most of the historic wildfires in California, this only represents the maximum patch size, so most fires would be substantially smaller. Whenever a population exceeded the carrying capacity, the model temporarily decreased the survival and transition probabilities, causing a gradual decline in abundance back to carrying capacity (Conlisk et al., 2012). To reduce the high variability in population size caused by fire-induced germination, population sizes of each patch only counted the number of reproductive adults.
We did not simulate dispersal because Ceanothus species have very limited dispersal and are unlikely to disperse between patches on the scale we modeled (Holl et al., 2000;Regan et al., 2010).

Fire
We modeled fire occurrence in each patch using a discrete time Weibull hazard function. This distribution uses two parameters b (scale) and c (shape), and as long as c > 1 fires are more likely to occur the longer it has been since the last fire event, emulating fuel accumulation. The hazard function is the probability that a fire occurs in the patch at time step t given that the last fire occurred T t ( ) time steps ago. Thus, for each patch k (given b k , c k , and T k (t)), the probability of a fire occurring at time step t is During each time step t and each patch k, we generated a uniform random number between 0 and 1, and if the number was less than l k k T t ( ) é ë ù û , a fire occured in that patch. When a fire occurs, individuals in the population in patch k have modified survival, transition, and fecundity values for the time step in which the fire occurs.
The b parameter is a unitless value that determines the scale of the average time between fires; with higher values of b , there is a longer average time between fires. To create spatially heterogeneous fire probabilities, we derived the initial value of the b parameter for each patch based on historical fire data in California. Using polygons of fire perimeters that occurred between 1878 and 2020, 1 we calculated 1 https://frap.fire.ca.gov/frap-projects/fire-perimeters/ the mean number of years between each fire event for each 270 × 270 m grid cell, using this as the mean historical fire return interval for that cell (Figure 1). For the grid cells where fires did not occur over this time period, we assigned a fire return interval as 142 years (the entire timespan of this dataset). In our model, fire occurs at the scale of an entire patch. To estimate a baseline value of b for a patch, we averaged the mean time between fires across all cells in each patch across the entire range of C. perplexans.
The c parameter value determines the fuel-dependence of the fire, such that fire likelihood changes exponentially with time since the previous fire relative to c . The fuel-dependence of an area is likely to depend on the density and types of plant species. We set c = 1 42 . , consistent with previous models of Ceanothus species in mixed chapparal in southern California where fire is thought to result from a combination of weather and fuel accumulation (Polakow et al., 1999;Regan et al., 2010;Syphard et al., 2013).
To understand the importance of spatial patterns in fire, we compared two general patterns of fire frequency distribution around the species' range. We compared population simulations under the spatially explicit heterogeneous fire frequency distribution (described above) with a control baseline with the same mean fire return interval, but all patches having identical fire return interval. To understand the importance of overall fire frequency, we also compared simulations over a range of overall mean fire return intervals by modifying the b parameter for each patch proportionally the same for eadch patch. Overall, we compared mean fire return intervals varying at 5-year intervals between 5 and 120 years. Under these uniform fire frequency distributions, we set all patches to have equal values of b . For heterogeneous distributions, we multiplicatively adjusted the baseline b values of each patch such that the average fire return interval was equal to the desired value.

Simulation and analysis
We used the software RAMAS GIS ® (Akçakaya and Root, 2005) to connect habitat suitability maps and fire occurrence probabilities with population models and to run our population simulations. We ran 10,000 replicate stochastic simulations of each combination of the two climate models (CNRM-CM5, HadGEM-ES2, plus a set of simulations with no climate change), two emissions scenarios (RCP 4.5, RCP 8.5), 24 mean fire return intervals (5-120 years), and two fire frequency distributions (spatially explicit, uniform), with 240 variations in total.
We determined initial population sizes by simulating population dynamics for 90 time steps without changes in habitat suitability and mean fire return interval equal to the baseline. Because obligate seeding dynamics can dramatically change local population sizes from year-to-year, there can be high temporal and spatial variability in population densities (especially with seeds and seedlings), even when the population is near carrying capacity. Because of this, it was difficult to assign initial population densities for our simulations. Therefore, we initially set each population with 109,300 germinant seedlings per 270 × 270 m cell. For each replication, we used the population in each patch at the end of the burn-in period of 90 time steps as the initial population. We then ran the population model over the 90-year climate change simulation. At this point, we adjusted the fire return interval of Frontiers in Ecology and Evolution 07 frontiersin.org each population to match the fire frequency we were testing, and the habitat suitability changed with climate change in each time step.
To understand how the independent variables of climate models, emissions scenarios, and mean fire return intervals affect the final population sizes of individual subpopulations in the model, we used a linear model to calculate the percentage of deviance explained by each independent variable on the final abundance for each subpopulation. To simplify the analysis, we used a binary comparison for mean fire return intervals (in particular short: 20 years vs. long: 80 years).

Results
The total final abundance of Ceanothus perplexans was greatest in simulations with intermediate fire return intervals with a mean between 30 and 50 years ( Figure 2). This peak around 30-50 year fire return intervals occurred with all variations of climate model, emissions scenario, and whether or not the spatial distribution of fire frequencies was uniform or heterogeneous. Shorter fire return intervals were more likely to cause extinction of a subpopulation than longer fire return intervals ( Figure 3). Expectedly, higher emissions scenarios resulted in much lower relative abundances. Populations under the HadGEM-ES2 climate model (which generally projected hotter and drier conditions; Pierce et al., 2018) resulted in somewhat higher relative abundances than populations under the CNRM-CM5 climate model.
Even with optimal fire return intervals, a majority of individual subpopulations decreased in abundance relative to initial subpopulation abundances (Figure 3); the decrease in final abundance for longer fire return intervals was the result of decreased abundance rather than local extinction. When fire return intervals were uniform, subpopulations only went extinct from fire when the fire return interval was 5 years (Supplementary Figure S1).
Though the final relative abundance of total metapopulation did not strongly depend on whether fire probabilities were uniform or spatially explicit (Figure 3), the final relative abundance of subpopulations did depend on geographical patterns in fire return intervals (Figure 4). In particular, when fire return intervals were uniform across the species' distribution, all subpopulations had a similar final abundance. However, when fire return intervals were spatially explicit, there was greater variability in the final relative abundance of populations. Moreover, this variability was greater with shorter mean fire return intervals ( Figure 4A) than with longer fire return intervals ( Figure 4B).
With spatially explicit fire return intervals, emissions scenario (RCP 4.5 vs. RCP 8.5) was generally the most important factor for explaining the deviance of the final relative subpopulation abundance ( Figure 5A). However, the climate model (CNRM-CM5 vs. HadGEM-ES2) and the fire return interval (20 years vs. 80 years) explained the greatest amount of deviance in determining the relative abundance of several subpopulations. The fire return interval was important for subpopulations on the northwestern and southwestern edge of the species' initial range ( Figure 5B). When fire return intervals were uniform across the species' range, they were unlikely to have a large effect on the relative abundance of individual subpopulations (Supplementary Figure S2).

Discussion
In the coming century, the fate of Ceanothus perplexans and other fire-adapted plant species will ultimately depend on changes in both the climate and fire regime, and the rate of these changes will not be geographically uniform (Keeley and Syphard, 2016). Our simulations showed that the total abundance of C. perplexans near the end of the century depended strongly on changes in the fire regime with an optimal fire return interval of around 30-50 years, supporting previous results (Regan et al., 2010;Syphard et al., 2013). While the effect on the total population was similar regardless of the fire frequency distribution across the landscape, we found that individual subpopulations depended much more on the spatial distribution of fires. In other words, some populations of C. perplexans are likely to be more threatened by climate change while other populations are more likely to be threatened by too-frequent fire.
Ceanothus perplexans benefits from fire as long it does not occur too frequently. Ceanothus perplexans require fire for seed germination (i.e., they rarely germinate without fire), but fire generally kills all stems in the population (Keeley and Zedler, 1978). Though it is not optimal, populations of C. perplexans can persist in environments with more than 80 years between fires (Regan et al., 2010;Syphard et al., 2013). However, fire return intervals that are more frequent than the optimum are more detrimental to populations, as subsequent fires can kill plants before they reach maturity and produce seeds. Especially with obligate seeding plants that depend on fire but cannot tolerate fires occurring too often, the effect of shifting fire regimes is likely to be complex and nonlinear. This complexity can be even further complicated because obligateseeding species recovery can be hindered by other factors associated with climate and land use change, such as drought (Jacobsen and Pratt, 2018), and competition with non-native annuals that are replacing chaparral in many human-dominated parts of the state . Because of this, we found that spatially explicit patterns in fire frequency can have a strong effect on individual subpopulations without having a large effect on the total metapopulation. An overall average around the optimum of 30-50 year fire return intervals could mean that many subpopulations will individually have both shorter or longer fire return intervals. This means that an obligate seeder would require different management strategies throughout different areas of its range depending on whether they have high or low fire return intervals.
Especially due to increased ignitions associated with population and urban growth in southern California, fires have exceeded their historical frequency in recent decades (Safford and Van De Water, 2014), and may continue to increase with climate, land use, and vegetation type change (Keeley and Syphard, 2021;Fusco et al., 2022). Past studies suggest that climate change and fire will have a larger effect than human development on future populations of C. perplexans (Syphard et al., 2007(Syphard et al., , 2013. If fire probabilities increase throughout the species' range, many areas will be more at risk of local extinction. Assuming that fires become more frequent throughout the entire species' range, areas with historically lower fire return intervals might see smaller reductions in population size, if not an increase in population, as fire frequencies move closer to the species' optimal fire return interval. Our model assumes that the Frontiers in Ecology and Evolution 08 frontiersin.org parameters controlling the population dynamics of all populations are identical throughout the species' range. However, populations of plants could be locally adapted to local fire regimes (Schwilk, 2003;Kane et al., 2022). Given C. perplexans' relatively long generation time and limited dispersal ability, low gene flow could create conditions for local adaptation. Populations in areas that historically have higher fire frequencies could have higher germination rates, and populations in areas with lower frequencies could have higher survival rates for adults and seed banks. In this situation, any shift in fire regime could reduce population sizes regardless of direction. However, anthropogenic influence of fire regimes means that populations might not be adapted to local fire conditions in evolutionary time (Keeley and Pausas, 2022). Our model considered the novel inclusion of heterogeneous fire probabilities over space, but the simulated fires always had the same size and intensity. In particular, all simulated fires occurred at the scale of entire subpopulations, ranging from a few to 6,561 hectares. Though our maximum-sized fires might be considered large in a historical context, they are well below the size typically considered a megafire (Linley et al., 2022), and large wildfires in California are becoming more frequent (Keeley and Syphard, 2021). Model limitations did not allow the size of fires to change, but we can speculate how this would affect populations of C. perplexans. Even with increased fire frequencies, smaller and uncorrelated fires might spread out the risk of driving an entire population of an obligate seeder extinct. In other words, the "immaturity risk" of multiple fires occurring in rapid succession and killing the population before it reaches maturity and produces seeds could be spread out among multiple subpopulations when fires are smaller (Regan et al., 2010). If large wildfires in become more common in California, or if extreme fire-weather conditions drive fire spread any benefits of smaller population sizes might be negated.
Though our model suggests that C. perplexans is at risk of population decline from both climate change and increased fire frequencies, and the causes of each are difficult to counteract, spatially explicit models can help to identify priority areas for conservation within a species' range. We found that populations across California are at risk from climate change and many areas in the southern end of their distribution within the United States are likely to decrease under most variations on climate change that we modeled, often driving them to local extinction regardless of fire frequencies. Conversely, changes in fire frequencies were more important for populations near the northwestern edge of the distribution. While most of the distribution of C. perplexans is in higher elevation chaparral in remote protected national forests, the northwestern edge is in the San Jacinto Mountains, nearer to high human population density. Bordering urban areas, fire has been historically more frequent and could continue to become more frequent from human ignition events. Thus, conserving C. perplexans could involve fire suppression in the northwestern part of its distribution, but fire management might not benefit the majority of the southern and central areas of the distribution in California where climate change has a larger effect. One potential solution is the assisted migration of the population farther north and west to the Transverse Ranges, though those regions are likely to undergo their own climate and Final abundance of the total metapopulation relative to the initial abundance of the total metapopulation across different mean fire return intervals. Color represents different emissions scenarios in the climate model, left/righ represents different climate models, and symbol represents how fire frequency was distributed around the landscape.
Frontiers in Ecology and Evolution 09 frontiersin.org fire regime changes. One limitation of these results is that we did not model the distribution of the species outside of the CFP in the United States, despite this species existing in Baja California, Arizona, and Nevada (Hickman, 1993). This limitation could bias our SDM and our estimates of habitat change (overestimating the effect of climate change on the southern edge of our modeled area), but this bias is likely to be small when considering similar climates between high elevations in California and Baja California Final abundance of each of the 149 different subpopulations relative to the initial abundance of those subpopulations across different mean fire return intervals under spatially explicit fire frequency distributions. Individual subpopulations are stacked such that subpopulations that increase in abundance are at the bottom, while subpopulations that decrease in abundance or go locally extinct are at the top. Climate models and emissions scenarios are arranged in columns and rows, respectively.

FIGURE 4
Final abundance of each subpopulation relative to the initial abundance of the subpopulation across different mean fire return intervals (A/B), spatial patterns in fire frequency (color), climate model (columns), and emissions scenario (rows).
Frontiers in Ecology and Evolution 10 frontiersin.org (Regan et al., 2012). When considering populations of C. perplexans outside of California, it will also important to understand differences in historical and future fire regimes. In particular, historical differences in fire regimes between Baja California and southern California (Minnich and Chou, 1997) could result in different effects on the fate of populations growing in either region. It can be difficult to understand the complexities of fire and its effects on ecosystems, especially when it comes to predicting future fire dynamics under climate change, but we can hope to make smart conservation decisions under various outcomes with what we know (Keeley and Syphard, 2016). For obligate seeding species, we know that they are best adapted to certain optimal fire and climate conditions; those conditions can be relatively broad, considering that wildfires do not occur on a regular timetable. However, as we enter unprecedented conditions, building ecological models under a range of climatic and fire regime conditions help us prepare for conservation actions we might need to take.

Data availability statement
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/Supplementary material.

Author contributions
GB and HR concieved and designed the study and wrote the manuscript. GB constructed the population models and performed the analysis. MR, SV, and JF constructed and ran species distribution models. AS provided expertise on fire history, data, and modeling. MR, SV, JF, and AS edited the manuscript. All authors contributed to the article and approved the submitted version.

Funding
This research was supported by National Science Foundation grant #1853697 and California Strategic Growth Council, Climate Change Research award #CCR30009.

Acknowledgments
The University of California, Riverside occupies the ancestral lands of the Cahuilla, Tongva, Luiseño, and Serrano peoples and we are grateful for the opportunity to conduct this research on their homelands. We thank James Thorne and the California Department of Fish and Wildlife for providing vegetation survey data as well as Alan and Lorrie Flint for providing the Basin Characterization Model hydroclimate data. Percent deviance of each independent variable in explaining the relative final abundance of each subpopulation relative to the initial abundance of the subpopulation under spatially explicit fire return intervals. (A) Distribution of explained deviance among all subpopulations and independent variables. (B) Spatial representation of explained deviance from independent variables across the estimated range of C. perplexans (excluding 15 small subpopulations outside of this map in the northwest).