The role of spatial structure in at-risk metapopulation recoveries

Metapopulations are often managed as a single contiguous population despite the spatial structure underlying their local and regional dynamics. Disturbances from human activities can also be spatially structured with mortality impacts concentrated to just a few local populations among the aggregate. Scale transitions between local and regional processes can generate emergent properties whereby the whole system can fail to recover as quickly as expected for an equivalent single population. Here, we draw on theory and empirical case studies to ask: what is the consequence of spatially


INTRODUCTION
The interplay between demographic processes across spatial scales can play an important role in how spatially structured populations, like metapopulations, recover from disturbances (Hanski, 1998;Melbourne & Chesson, 2006;Thomas & Kunin, 1999).For example, dispersing silver-spotted skipper butterflies, Hesperia comma, re-colonize extirpated patches, which buffers their metapopulation from widespread extinction (Davies et al., 2005;Hanski & Thomas, 1994).However, recoveries also depend upon local population dynamics, not just dispersal from neighboring patches.For example, local overfishing in Pacific herring Clupea pallasii arises from feedback between local density-dependent population regulation, dispersal, and harvest location, which can result in localized collapses even as other populations remain abundant (Benson et al., 2015;Okamoto et al., 2020).The interplay between local and metapopulation processes can challenge effective conservation and management because they shape context-dependent pathways toward recovery (Chesson, 2012;Zelnik et al., 2019).
Metapopulation recovery dynamics arise from feedback across spatial scales between dispersal and productivity, often regulated by density dependence (Howe & Davis, 1991;Melbourne & Chesson, 2006;Pulliam & Danielson, 1991).Compensatory density dependence in local patches results in per capita increases in productivity as abundances drop below local carrying capacity (Rose et al., 2001).Dispersal can subsidize local abundances and rescue extirpated patches, while lowering abundance and potentially triggering compensatory dynamics in source patches (Thomas & Kunin, 1999;Watkinson & Sutherland, 1995).Metapopulation recoveries thus span a continuum that ranges from high productivity and low dispersal, where patches depend on their own internal dynamics for recovery ("isolated recovery regime"), to low local productivity and high dispersal, where patches rely on external rescue effects via dispersal ("rescue recovery regime") (Zelnik et al., 2019).The balance between local regulation and regional dispersal can control source-sink dynamics, which ultimately structures spatial recovery regimes and metapopulation resilience, defined as the ability of metapopulations to compensate and to recover from disturbances (Hanski, 1998;Levins, 1969;Thomas & Kunin, 1999;Zelnik et al., 2019).
The spatial distribution of disturbances from harvest, habitat losses, or perturbations among local populations can further alter metapopulation resilience (Fahrig & Merriam, 1985;Zelnik et al., 2019).In fisheries and wildlife management, harvesting across several populations at once (e.g., mixed-stock fishing) is a common form of disturbance, with sustainable targets determined through broad-scale assessments of the entire metapopulation or aggregate (McCullough, 1996;Moore et al., 2021;Ricker, 1958).If total impacts are evenly distributed among local populations relative to their recovery potential, then spatial recovery patterns may be uniform as the aggregate-scale targets reflect the sum of local patches, creating a tight correspondence between monitoring signals and ecological dynamics (Cumming et al., 2006;Zelnik et al., 2019).However, uneven or asymmetric disturbances among patches can weaken the feedback between management and metapopulation dynamics, allowing disturbances or harvest to appear sustainable when they are not (Cumming et al., 2006;McCullough, 1996;Ricker, 1958).We illustrate this in a simple two-patch scenario, where a 50% loss in aggregate abundances is either shared evenly among both patches or completely extirpates a single patch (Figure 1).Although the total abundance lost to the whole system is equal, the scenario with local extirpation erodes the resilience of the entire aggregate by maintaining local density-dependent regulation, which can slow recoveries, aggregate recovery thus depends upon external rescue effects via dispersal.Thus, simple differences in the spatial pattern of disturbance can produce emergent outcomes such as masked local collapses, reduced patch occupancy, and lost metapopulation resilience (McCullough, 1996;Okamoto et al., 2020;Zelnik et al., 2019).
Scale mismatches between management and ecological processes can lead to the emergence of unexpected spatial dynamics that challenges management (Cumming et al., 2006;Okamoto et al., 2020).For example, the IUCN Red List often sets recovery targets at aggregate scales that span entire metapopulations or even species (Hutchings et al., 2012;Pfab et al., 2011).Monitoring an entire resource portfolio may signal high productivity (Appendix S1: Figure S1) that appears sustainable on the aggregate, but can mask local collapses (Okamoto et al., 2020;Ying et al., 2011).Alternatively, risk assessments may narrowly focus on impacts on local dynamics while overlooking regional processes, including the potential spatial accumulation of risks from organisms dispersing through impacted patches (Moore et al., 2015).Such mismatches in scales of management and relevant ecological processes can decouple local population dynamics from aggregate performance metrics, leading to unexpected ecological or management outcomes that can be overlooked if the spatial structure is not accounted for in monitoring programs (Fullerton et al., 2016;Melbourne & Chesson, 2006;Post et al., 2021).
Advances in theoretical and statistical ecology have allowed practitioners to address many challenges that arise from scale dependence in metapopulation dynamics.For example, metapopulation theory has been foundational to spatial management designs, like marine protected area networks, as optimal designs can depend on dispersal patterns and local dynamics (Aalto et al., 2019;Halpern, 2003;Sanchirico & Wilen, 2005).Furthermore, statistical advances like spatial-temporal hierarchical models improve the estimation of the scale transitions between local and aggregate dynamics (Anderson et al., 2022;Thorson, 2019;Thorson et al., 2019), and have been increasingly used in quantitative assessments to inform F I G U R E 1 Metapopulation recovery regimes (black) depend on how disturbances impact two local patches (green and orange) at three time periods (t 1 = pre-disturbance phase, t 2 = disturbance phase, and t 3 = recovery phase) and whether impacts are equally (a) or unequally (b) distributed.Lower panels demonstrate density-dependent population dynamics at each time period.(a) Two abundant local patches are equally impacted by the same disturbance reducing both patches to 50% capacity, and each recover from local density dependence alone.(b) A disturbance of the same magnitude but unequally applied to the two patches during t 2 depletes one patch to 0% (orange) while maintaining the other patch at 100% capacity (green).Recovery of the depleted patch (orange) during t 3 depends on the potential for rescue effects from the green patch.Without rescues (e.g., if disturbances also affect connectivity, like dams), then metapopulation resilience may be permanently lowered (* indicates this context dependence) leading to a 50% loss to the metapopulation.spatial management, particularly in data-rich contexts (Cadrin et al., 2018;Kerr et al., 2016;Punt et al., 2015).Despite these advances, issues continue to arise that challenge management including sparse monitoring of spatial structure to inform assessments (Post et al., 2021) and uncertainties of the context dependence underlying spatial recovery regimes that may emerge from the interplay of local and regional processes (Zelnik et al., 2019).
Here, we combined applications from a theoretical metapopulation model and reviews of empirical case studies to examine how feedback between local and regional processes mediate spatial recovery regimes and how these emergent patterns can challenge management and conservation.Overall, we propose that managing successful recoveries of disturbed metapopulations requires an understanding of the context and scale dependence resulting from the interplay of key demographic processes like local dynamics driven by density dependence, regional dispersal, habitat networks and topography, and spatial disturbances (Figure 2).To demonstrate this, we first used a simulation model to explore how the interplay between these demographic processes affects metapopulation recovery criteria, such as recovery rate.We then illustrated many real-world challenges faced by management through a literature review of three species listed under the US Endangered Species Act (ESA): Florida Everglade snail kites (Rostrhamus sociabilis), California and Alaska sea otters (Enhydra lutris), and Snake River Chinook salmon (Oncorhynchus tshawytscha).Together, these case studies spanned a range of spatial processes, recovery outcomes, and management challenges germane to at-risk metapopulations.With this understanding, we present some advice to practitioners tasked with managing metapopulation recoveries and identify opportunities for applied research to continue tackling these challenges.

Rationale of key processes
We developed a numerical model to explore the context dependence of metapopulation recoveries along gradients of multiple local and regional demographic processes.We hypothesize that the spatial recovery regimes of disturbed metapopulations emerge from how four key processes interact and transition across scale (Melbourne & Chesson, 2006; Figure 2): 1. Local density-dependent population dynamics determines the amount of surplus produced as local populations drop below carrying capacity.Higher compensation (i.e., per capita productivity) may allow patches to self-recover at faster rates while at low densities (Bowlby & Gibson, 2020;Rose et al., 2001).Similarly, higher carrying capacities can reduce the strength of density dependence over time, all else equal, helping to avoid local patches from slowing metapopulation recoveries (i.e., local bottlenecks).2. Dispersal rates control the proportion of organisms that move between neighboring patches and create sourcesink dynamics.Dispersal from neighbor patches can rescue disturbed patches but at the loss of individuals in source patches that may also be slow to recover (Bradford & Braun, 2021).However, spatial recoveries can become limited if organisms disperse into patches already at their carrying capacity and are lost due to density dependence.Increased dispersal can improve metapopulation recovery rates but can either increase extinction risks through population synchrony and spatial-temporal correlations (Earn et al., 2000) or buffer against risks via stabilizing portfolio effects (Moore et al., 2021).3. The spatial arrangement of patches within a habitat network dictates which new habitats dispersing individuals move into, thereby shaping the extent of rescue potential (Earn et al., 2000;Fahrig & Merriam, 1985).For example, networks may be arranged linearly (e.g., river or coastline) or they may be highly interconnected (e.g., forests or meadows).Thus, "connectivity" is shaped by both network complexity and dispersal rate.Well connected habitat networks tend to enhance resilience compared with networks with weak or sparse connections (Campbell Grant et al., 2007;Sanchirico & Wilen, 1999), but excessive connectivity can increase synchrony and erode resilience (Earn et al., 2000).
4. The spatial allocation of disturbances will also impact recovery patterns by controlling the balance between local and regional dynamics and influencing whether patches can be rescued (Atlas et al., 2015;Okamoto et al., 2020;Zelnik et al., 2019).Disturbances allocated evenly across patches may achieve quicker recovery than metapopulations with uneven disturbances.

Model description
Following the above concepts, we defined our metapopulations as a set of P local populations for a species with a 1-year generation time with population dynamics that follows basic birth (i.e., recruitment), immigration, mortality, and emigration processes typical to metapopulation and portfolio theory (Anderson et al., 2015;Bowlby & Gibson, 2020;Fullerton et al., 2016): where N i,t was the number of adults in patch i at time t (i.e., the remaining individuals after recruitment, dispersal, and death have occurred), R i,t the number of within-patch recruits, P P j¼1 j ≠ i ωp i,j R j,t was the number of recruits immigrating into patch i from any other patch, ω set the proportion of local recruits to disperse, and d i,t the proportion of post-dispersal recruits that died from the disturbance regime.Density-independent dispersal was shaped by a distance-decay function p i,j such that migrating recruits most often moved into their closest neighbor patch (Appendix S1: Section S1.1.3).Local recruitment at time t depended on adult densities at t − 1 following a Beverton-Holt function (Walters & Martell, 2004): where α i was the recruitment compensation ratio, β i was the carrying capacity for patch i, and ϵ i,t represented lognormal recruitment deviates used to generate stochasticity.We focus on the Beverton-Holt function rather than the Ricker function to model density dependence for two main reasons: first, asymptotic rather than over-compensatory density dependence was more relevant to our case studies and, second, the Beverton-Holt was more straightforward to standardize demographic parameters among scales and scenarios.Both the magnitude of variance and the spatial-temporal correlations in local recruitment deviates depended on the scenario (Appendix S1: Section S1.1.5).
Natural resource managers often monitor and manage resource portfolios at the scale of the metapopulation (Anderson et al., 2015;Moore et al., 2021;Okamoto et al., 2020).Hence, metapopulation patterns arise from the aggregation of local patch dynamics.We therefore defined metapopulation adults A t as the adult abundances summed across patches: with metapopulation recruits M t as: Expected adult recruitment E M t ð Þ among the whole aggregate set of populations, presuming no spatial structure, follows: where α was the compensation ratio averaged across the metapopulation, β was the carrying capacity summed across all populations, and ϵ t were lognormal recruitment deviates generated using identical scale and variance parameters as in Equation (2).We then compared the value in Equation ( 5) to the realized, actual value in Equation (3) to calculate potential shortfalls of ignoring spatial structure (i.e., monitoring as a single contiguous population).
We developed a numerical model of metapopulation dynamics integrating the above processes and explored emergent spatial and temporal recovery patterns across 9504 ecological scenarios, with each scenario replicated 100 times to evaluate the consequence of stochastic processes (Appendix S1: Section S1.3).Briefly, our ecological scenarios tested gradients in (1) heterogeneity in local patch productivity (i.e., compensation ratio), (2) heterogeneity in local patch carrying capacity, (3) dispersal rates, (4) the spatial allocation of disturbances among patches, (5) topography of the habitat network, (6) stochastic lognormal recruitment dynamics, (7) temporal correlations in recruitment dynamics, and (8) spatial correlations in recruitment dynamics (e.g., Moran effects).
The total disturbance was the same in all scenarios: we reduced the metapopulation by 90% in a single year following the IUCN criterion for critically endangered status (COSEWIC, 2019;IUCN, 2015).Disturbance scenarios were developed to test gradients in how spatial variation in the disturbance regime among local patches shapes recovery dynamics (see details in Appendix S1: Sections S1.3 and S1.3.1).In each scenario, the proportional loss to each patch was a random draw and each scenario varied in the randomization algorithm used to determine how those losses were distributed.Scenarios included (i) all patches were equally vulnerable to a 90% loss in local abundance such that the total loss summed to a 90% aggregate loss (hereafter termed "uniform"; Appendix S1: Figure S7), (ii) disturbance impacts were concentrated to a random subset of local patches and each patch among that subset had an equal probability of receiving the same proportional loss of individuals such that the total loss summed to a 90% aggregate loss (hereafter termed "local even"; Appendix S1: Figure S8) or (iii) disturbance impacts were concentrated to a random subset of local patches such that total extirpation of each patch among that subset would result in a 90% loss on the aggregate (hereafter termed "local uneven"; Appendix S1: Figure S9).These scenarios provide a useful range of how spatial disturbances can impact metapopulations in different ways, for example, spatial asymmetries in how harvest is prosecuted among mixedstock fisheries can manifest in relatively uniform or uneven impacts on local populations (Moore et al., 2021).For (ii) and (iii), we developed a numerical algorithm to find a set of conditions whereby disturbances (either a proportion or extirpation) to a random subset of local patches summed to a ~90% loss in abundance to the whole metapopulation and left at least one local patch undisturbed for potential rescue effects.Correspondingly, the local, even disturbance scenarios tended to be slightly less severe versions of local, uneven scenarios (all else being equal).After this pulsed disturbance, metapopulation recovery dynamics were allowed to arise naturally from their disturbed state.

Metapopulation recovery metrics
Here, we focused on characterizing the spatial recovery dynamics of disturbed metapopulations using four quantitative metrics: (1) recovery rate (i.e., the time for the metapopulation to reach ≥pre-disturbance abundance for five consecutive years), (2) patch occupancy (i.e., the proportion of patches with post-disturbance abundances ≥10% pre-disturbance), (3) relative production (i.e., the ratio between the summed adult abundances across all patches in Equation 3 and the expected adult abundances of an equivalent-sized contiguous population in Equation 5), and (4) rate of nonrecovery (i.e., the percent of 100 stochastic simulations where the metapopulation failed to recover (as defined in 1) within 50 years after the disturbance; Appendix S1: Section S1.4).We used these recovery metrics as they provide multiple indices of metapopulation recoveries beyond typical measures of population viability or extinction risks (Wolf et al., 2015).For example, the recovery rate measures how quickly a metapopulation recovers to a target baseline.Patch occupancy measures how many local patches recovered to a minimum value (10% pre-disturbance) or, conversely, how many remain collapsed 25 years after disturbance.Relative production provides an indicator of the impacts of spatial structure on adult recruitment compared with nonspatial expectations (i.e., the consequence of monitoring at the aggregate scale despite population regulation occurring at local scales), measured 25 years after disturbance.The rate of nonrecovery quantifies the frequency of metapopulations failing to reach their target baseline due to stochastic population dynamics.We then used hierarchical clustering analyses of the above quantitative metrics across all ecological scenarios to evaluate whether there were common recovery dynamics (i.e., recovery regimes) with all analyses implemented in R version 4.3.0(Hennig, 2014;Murtagh & Legendre, 2014;R Core Team, 2023).Hierarchical clustering analyses were useful for our purposes in providing statistical support for the unique recovery regimes that emerged from the simulations, based on a dissimilarity matrix of the above four recovery metrics, without imposing strict grouping criteria a priori (Appendix S1: Section S1.6).The clusters described qualitatively similar modes among an otherwise continuous range of recovery outcomes.We then assessed how ecological and disturbance conditions were associated with these modes using ordered logistic regression in the "MASS" package (Venables & Ripley, 2002) where clusters (ordered from "best" to "worst" based on recovery metrics) were the dependent variable and the ecological and disturbance scenarios as predictors.Hence, combining clustering analyses with ordered logistic regression allowed us to infer both how many and what recovery regimes commonly emerged from the interactions of ecological and disturbance scenarios, helpful for crafting qualitative heuristics or "rules of thumb" in decision-making.

Overview of case studies
We reviewed available published literature from three species listed under the US ESA (Everglade snail kites, California and Alaska sea otters, and Snake River Chinook salmon) to highlight real-world challenges in managing spatial recoveries in at-risk metapopulations.These case studies spanned a gradient of common spatial processes and recovery outcomes to at-risk metapopulations that were directly relevant to our theoretical model.Overall, we reviewed peer-reviewed articles, government technical reports, and species recovery plans (NMFS, 2017a(NMFS, , 2017b;;USFWS, 2003USFWS, , 2013USFWS, , 2019)).In brief, each case study represents a highly disturbed metapopulation with a patchy habitat network of local populations that vary in productivity and are intermittently linked via dispersal.While disturbance typologies varied among species (e.g., overharvesting, habitat fragmentation, or degradation), each of their respective spatial disturbance regimes led to relatively uneven impacts distributed among local populations that emerged as substantial declines to the metapopulation sufficient enough to trigger federal listing.The relevant federal management authorities for each case study (US Fish and Wildlife Service [USFWS] for sea otters and snail kites, and National Marine Fisheries Service [NMFS] for Chinook salmon) initially established quantitative recovery targets at aggregate spatial scales and were subsequently challenged by unexpectedly slow recoveries after listings, with spatial structure and scale transitions playing an important role in each case.Each case study was also supported by relatively data-rich monitoring programs that allowed for scientific insight into some of the spatial processes critical to their respective metapopulation dynamics and recovery plans.In the Results below, we will first overview each case study's history, spatial recovery dynamics, and management challenges and link how these outcomes relate to results from our theoretical model.

Theoretical model
Metapopulation recoveries emerged from contextdependent feedback between local and regional processes.For illustration, we show two example recovery regimes for metapopulations that varied only in network complexity when disturbances were locally uneven (Figure 3).We first examined a metapopulation with linear network complexity composed of 16 identical patches, a dispersal rate of 1%, and low stochasticity (0.1%).This lower dispersal rate was chosen to allow metapopulation recovery dynamics to emerge from a balance between local dynamics and regional dispersal.As shown in Figure 3a,b, the linear network has a slowed and modular recovery pattern where local recoveries cascaded outward from the last remaining patch (e.g., similar to California sea otters; see Empirical Case Studies below).In contrast, subregions within the dendritic network (Figure 3c,d) were able to recover quickly as dispersal quickly rescued neighboring patches (e.g., similar to Snake River Chinook salmon, see Empirical Case Studies below), similar to recovery in a river network.We further demonstrated the role of these ecological and disturbance processes in a walkthrough that iteratively layers complexity to an example metapopulation in Appendix S1: Section S1.3.2.
The spatial allocation of disturbance was a strong determinant of metapopulation recovery patterns (Figure 4; Appendix S1: Figures S16-S21).Locally even disturbances tended to reduce recovery rates by 16% compared with uniform disturbances and increased the risk of local patch depletions (i.e., spatial contraction; Appendix S1: Figure S16a).The poorest recovery outcomes emerged when the disturbance was locally uneven.In these scenarios, recovery rates were slowed by 42% (from 0.82 to 0.47 on average; Appendix S1: Figure S16a), relative production was reduced by 28% (Appendix S1: Figure S16b), long-term patch occupancy was reduced by 59% (Appendix S1: Figure S16a), and the risk of nonrecovery increased almost seven-fold from 2.3% to 18% (Appendix S1: Figure S16c).By altering aggregate production from reduced patch occupancy, localized disturbances eroded metapopulation resilience.In other words, through changes in source-sink dynamics, metapopulations under local disturbances acted less than the sum of their parts: the more uneven the local disturbances, the worse these effects.Conversely, uniform disturbances typically allowed for resilient metapopulation recovery dynamics.
Interactions between dispersal, network structure, and variability in local dynamics also affected metapopulation recoveries in several ways (Figure 4 and Appendix S1: Figures S17-S21).First, recovery rates increased with increased dispersal (Figure 4).However, this relationship was sensitive to disturbance, with stronger effects of dispersal under local, uneven disturbances and weaker effects for uniform disturbances.Second, more linearized networks had slower recovery times than more connected networks, as rescue effects were slow to cascade through the whole system (Figure 3 and Appendix S1: Figures S16 and S17).The slowed recovery of linearized networks was most common when disturbances were localized (Figure 4a).Furthermore, variability in local patch demographic rates (both productivity and carrying capacity) slowed metapopulation recoveries but high dispersal rates tended to dampen this effect (Appendix S1: Figure S17).Last, interactions between the magnitude of recruitment variability and the strength of spatial-temporal correlations affected recovery dynamics, particularly for linear and dendritic habitat networks under locally uneven disturbances (Appendix S1: Figures S18-S21).For example, scenarios with high recruitment stochasticity showed less variation in recovery rates than scenarios with low stochasticity because of stabilizing effects generated by asynchronous dynamics (Appendix S1: Figure S18).However, when the magnitude of recruitment variability and the strength of spatial-temporal correlations were both high, there was less variation in most recovery metrics compared with scenarios with low spatial-temporal correlations and both patch occupancy and relative produced were reduced (Appendix S1: Figures S18-S21).Conversely, when recruitment stochasticity was low, high spatialtemporal correlations in recruitment increased the variation in recovery rates than for equivalent scenarios with low spatial-temporal correlations.
A few common types of spatial recovery patterns emerged from the interplay between ecological and disturbance conditions (Table 1; Appendix S1: Figure S23).Overall, clustering analyses of our model results found evidence for five common outcomes (Appendix S1: Figure S23).We defined these modes as (Table 1 and Appendix S1: Table S2): (1) resilient recovery: metapopulation recovered to pre-disturbance abundances relatively quickly with all local patch occupancy was widespread, (2) slow recovery: metapopulation recoveries were either slowed (compared with resilient recoveries), had reduced patch occupancy, or reduced relative production, (3) hidden collapses: the metapopulation recovered to pre-disturbance abundances but many local patches were still collapsed, (4) lost capacity: metapopulation recovery rates were slow, many local patches were still collapsed, long-term relative production was ≤80% pre-disturbance, and the risk of nonrecovery from stochastic processes was relatively high, and (5) critical risk: the metapopulation failed to return to recover, abundances remained ≤10% of pre-disturbance, and the risk of nonrecovery due to stochastic processes was highest.
In general, the five outcomes we identified spanned a continuum of recovery patterns that ranged from better (resilient) to worse (critical risks).Furthermore, the probability of any particular recovery outcome was modulated by the interplay between ecological and disturbance conditions (Figures 5 and 6; Appendix S1: Figures S24--S27).For example, resilient recoveries were highly probable for any metapopulations subjected to uniform disturbances (Figure 5; Appendix S1: Table S2).When metapopulations were subjected to locally even disturbances, the probability for resilient recoveries was increased under conditions of high dispersal, homogenous local demographic rates (α and β in the Beverton-Holt model), low spatial-temporal recruitment variability, and gridded habitat networks (Figure 5).However, resilient recoveries were highly improbable for metapopulations subjected to locally uneven disturbances, where outcomes instead ranged from slow recovery (at best) to critical risks (at worst).While a relatively infrequent outcome (Appendix S1: Figure S24), the probability of critical risks increased under conditions of locally uneven disturbances, low dispersal, variable local demographic rates, high spatial-temporal recruitment variance, and linear habitat networks (Figure 6).These dynamics illustrate some of the context dependence underlying how local and regional processes shaped metapopulation recoveries.
Collectively, these simulations revealed the overarching role of spatial structure in metapopulation recoveries and illuminated several factors that affected how at-risk metapopulations are likely to recover.While based on theory, these results provide insights into the real world by helping practitioners identify the underlying ecological context and knowledge gaps likely to drive or limit recoveries in focal ecosystems.Below, we explore case studies of three US ESA-listed species to illustrate how these insights can be applied to species metapopulation recoveries.

Spatial modularity in the recovery of Everglade snail kites
The Florida Everglade snail kite was listed as a critically endangered species under the initial US ESA of 1967 (USFWS, 2016).Snail kites inhabit much of Southern and Central Florida's large swathes of wetlands with their largest historical population in the Greater Everglades (Figure 7).Beginning in the late 1800s, habitat losses from development and drought in the Greater Everglades led to a decline in snail kites, which ultimately prompted their listing under the original ESA.In 1997, the snail kite monitoring program was refined to account for spatial structure, leading to improved estimates of local and aggregate population sizes.In 1999, the USFWS set quantitative recovery criteria based on abundances for the entire metapopulation.In 2019, recovery criteria were amended to establish spatially explicit targets due to new insights from the improved monitoring framework (USFWS, 2019).The snail kite metapopulation began an unexpected increase in abundance in the mid-1990s, but this increase The probability of resilient metapopulation recoveries along gradients of habitat network structure, dispersal rates (low ≤0.001; high >0.001),spatial distribution of disturbance, heterogeneity in local demographic rates (α was local patch productivity and β was local patch carrying capacity), and spatial-temporal recruitment variation (high variation was for scenarios where ρ = 0.6 and σ = 0.1; low variation was for scenarios where ρ ≈ 0 and σ = 0.001, ρ denoted the scale of spatial-temporal correlations and σ was the scale of lognormal recruitment variance) based on ordered logistic regression.Note that the probability of resilient metapopulation recoveries for uniform and local, uneven treatments were unaffected by spatial-temporal recruitment variance and have overlapping circles and squares.
was not equally distributed among regions (Reichert et al., 2016).In the lake-dominated northern region, a sudden increase in prey resources following the introduction of an invasive freshwater snail spurred a localized population resurgence.The apparent rebound of snail kite abundances exceeded prescribed recovery baselines, sparking debates about whether snail kites had recovered and should be downlisted (Martin et al., 2007).However, the kites in the largely contiguous wetlands of the Everglades remained in a state of collapse.This regionalized snail kite recovery revealed that limited connectivity and nest successes would lead to critical bottlenecks and reduced long-term viability of the metapopulation (Bowling et al., 2012;Fletcher et al., 2013Fletcher et al., , 2015; Martin The probability of long-term critical risks in metapopulation recoveries along gradients of habitat network structure, dispersal rates (low ≤0.001; high >0.001),spatial distribution of disturbance, heterogeneity in local demographic rates (α was local patch productivity and β was local patch carrying capacity), and spatial-temporal recruitment variation (high variation was for scenarios where ρ = 0.6 and σ = 0.1; low variation was for scenarios where ρ ≈ 0 and σ = 0.001, ρ denoted the scale of spatial-temporal correlations and σ was the scale of lognormal recruitment variance) based on ordered logistic regression.Note that the probability of long-term critical risks in metapopulation recoveries for uniform and local, even treatments were unaffected by spatial-temporal recruitment variance and have overlapping circles and squares.
F I G U R E 7 Metapopulation dynamics, recovery outcomes, and management challenges emerge from interplay between local dynamics, regional dispersal, and spatial structure for three species listed under the US Endangered Species Act.Map and time series for Florida Everglade snail kite recovery dynamics recreated from Reichert et al. ( 2021) with text supporting key historical points from case study reviews and recovery plans (Fletcher et al., 2013(Fletcher et al., , 2015;;Reichert et al., 2021;USFWS, 2019).Map and time series for California sea otter recovery dynamics generated from USGS monitoring data provided by M. Tim Tinker, with historical footnotes from case study reviews and recovery plans (Lafferty & Tinker, 2014;Tinker et al., 2016;Tinker, Tomoleoni, et al., 2019;USFWS, 2003).Map and time series for Snake River Chinook recovery dynamics recreated from Northwest Fisheries Science Center (2015) and Thorson et al. (2014) respectively.Historical footnotes on Snake River Chinook time series supported by case study reviews and recovery plans (Achord et al., 2003;Fullerton et al., 2016;NMFS, 2017a;Walters et al., 2013).NMFS, National Marine Fisheries Service.et al., 2006).This was in part because the Everglades were a historical stronghold for snail kites and represented their largest potential habitat source, but limited connectivity between the regions would keep the northern recovery from spreading southward toward the Everglades, which remained impacted by drought and development (Martin et al., 2006(Martin et al., , 2008)).
The mismatch between the management and ecological scales of Florida's snail kites revealed that their apparent aggregate recovery was masking poor status in one key region (Reichert et al., 2021).Compared with the northern region, snail kites from southern regions show stronger site fidelity, less dispersal, and reduced nest successes made worse by the impacts of local drought and disturbance (Fletcher et al., 2013;Martin et al., 2006).Yet, the resurgence of the northern regions masked the continued decline in the Everglades, misinforming assessments of their recovery in the early 2000s, similar to our simulated metapopulations when conditions of localized disturbances and low dispersal rates allowed for hidden local collapses to occur despite high aggregate abundance (Appendix S1: Table S2).The improved monitoring program better aligned the spatial scales between population dynamics and recovery potential, but initial USFWS recovery criteria were set for the metapopulation as a whole.Hence, despite increasing abundances, the regionalized dynamics among Everglade snail kite populations would continue to bottleneck the long-term viability of the metapopulation and slow recovery.This case study reveals two important insights: (1) spatial complexity in population structure and monitoring designs can mask hidden collapses within metapopulations, and (2) using shifting historical baselines to define recovery criteria can lead to overly optimistic estimates of recovery (suggesting that downlisting may have been an inappropriate recommendation).Despite these debates, snail kites retained a protective listing and appear headed toward a more widespread but slow recovery (Reichert et al., 2021).(1967,1973) (Ralls et al., 1996).In North America, sea otters have naturally re-colonized or been re-introduced to much of their former range, and now inhabit coastal nearshore waters in California, Washington, British Columbia and Alaska.Population structure is regionalized into subspecies (Northern sea otters E. lutris kenyoni from Washington north and Southern sea otters E. lutris nereis in California) and local stocks within each subspecies (Davis et al., 2019;USFWS, 2003USFWS, , 2013)).For southern sea otters, the fur trade harvest reduced their abundance and distribution to a mere handful of kelp forests along southcentral California (Figure 7-California sea otter).In 1982, the USFWS began setting recovery targets for the whole subspecies, which encompassed all of the California coastline (Ralls et al., 1996).

Spatially disparate recoveries in sea otters
Beginning in the 1990s, a puzzle was posed for sea otter management: why were California sea otters recovering slower than the Alaska subspecies, despite concurrent federal protections of both stocks (Figure 8; Estes, 1990)?Specifically, California sea otter population growth rates (~3%-5% per year) were four times lower than the growth rates (~20%) of recovering populations of the northern subspecies in Alaska and British Columbia (Esslinger & Bodkin, 2009;Estes, 1990;Gregr et al., 2008).This slow recovery pressured the managing authority, USFWS, to consider shifting sea otter recovery baselines and criteria for eventual delisting (Ralls et al., 1996).Fortunately, new insights into sea otter spatial dynamics revealed that this apparent mystery may stem from a mismatch between the regional scale of management (subspecies or "stock") and the local scale of population regulation.For example, tagged otters were found to exhibit low dispersal, with lifetime home ranges of just 5-25 km of coast for most individuals (Tarjan & Tinker, 2016) with occasional longer dispersals when local patches reached high densities (Tinker et al., 2008).As a result, density-dependent population regulation occurred at the scale of tens of kilometers, a much finer scale than the 500 km coastline of the entire California stock (Tinker et al., 2021).Local population regulation in California led to the slowed and sequential recovery, as food-limited populations near carrying capacity slowly re-colonized neighboring patches of abundant prey resources, which subsequently led to rapid local population growth rates more similar to Alaska trends (Lafferty & Tinker, 2014).
The slow recovery of California sea otters was initially unexpected by management, in part, because northern sea otters in regions such as Southeast Alaska had recovered so quickly (Estes, 1990;Laidre et al., 2001).Despite almost identical demographic processes in local patches, the metapopulation dynamics of California and Southeast Alaska were entirely different (Figure 8).Spatial network structure emerged as one of the primary differences underlying the two recovery regimes (Tinker, 2015;Tinker et al., 2021) as seen in the difference between Figure 4a,d from our simulation models.Specifically, California sea otters were contracted by the fur trade harvest to a single source patch and, similar to our modeled linear networks (Figure 3a,b), their dispersal opportunities were limited to immediate neighbor patches either north-or southward (Lafferty & Tinker, 2014;Tinker, Gill, et al., 2019).Intriguingly, the recovery of California sea otters began from a single source patch and continued along a linear coastal network, but is now bottlenecked at the northern and southern extremes of the range that appears driven, in part, by increased mortality from the recent resurgence of great white sharks (Carcharodon carcharias; Moxley et al., 2019;Tinker et al., 2016).Meanwhile, the Southeast Alaska metapopulation was re-established from seven source patches, each dispersing and expanding among a complex network of islands, inlets, and bays (Tinker, Gill, et al., 2019), similar to our modeled grid networks.This case study reveals three important insights: (1) that accurate expectations of metapopulation recovery depend upon monitoring that accounts for spatial structure, (2) that habitat network complexity plays a key role in shaping recovery regimes, and (3) that spatially contracted metapopulations can become vulnerable to emergent properties that further limit their recovery.Despite these challenges, California sea otters maintained their protected status under the US ESA and recovered in much of their range in Central and Southern California.

The persistent collapse of Snake River Chinook salmon
After nearly a century-long precipitous decline, wild Chinook salmon (O.tshawytscha) in the Snake River (a tributary of the Columbia River) were listed as a threatened population under the US ESA in 1992 (NMFS, 2017a).In the late 1800s, an estimated 1.5 million spring and summer run Chinook salmon returned to the Snake River each year.By 1992, abundance had declined to <1% of historical baselines, driven by cumulative disturbances from dams, forestry, land use, overfishing, and ecosystem changes that affected the entire Columbia River basin (Matthews & Waples, 1991).Notably, the construction of 12 hydroelectric dams along the river blocked migration to spawning habitats and led to widespread spatial contractions (Matthews & Waples, 1991;Raymond, 1979Raymond, , 1988)).The management authority, NMFS, set multiple pathways for recovery criteria that recognized the spatial structure of Chinook salmon within the watershed (NMFS, 2017a): (1) >3000 total spawners in at least two patches (i.e., a Distinct Population Segment or Evolutionary Significant Unit), (2) >4500 total spawners in the metapopulation (i.e., the Major Population Group), and/or (3) 80% patch occupancy in the lower portion of the watershed.
Despite considerable effort toward their recovery, the persistent low abundance of Snake River Chinook salmon continues to challenge management with only a modest rebound in recent decades (Isaak et al., 2007).For example, commercial and recreational harvest of Chinook salmon has been reduced by more than 90% (NMFS, 2014) and more than 500 million USD per year has been invested in research, enhancement, and restoration of the Snake and Columbia Rivers (Northwest Power and Conservation Council, 2017).Since the collapse of the metapopulation, the number of recorded spawning sites has decreased, with many historically occupied areas remaining unused (Fullerton et al., 2016).In years with high adult returns, Chinook salmon expand their spatial distribution but still tend to cluster around core spawning areas.Moreover, river spawning fidelity is high (within 10 km for wild fish), which limits dispersal (Isaak et al., 2007;Isaak & Thurow, 2006).In theory, low spawner densities should trigger compensatory improvements to freshwater growth and survival.Yet, while metapopulation densities remain low, the spatial clustering of adult Chinook salmon in select streams triggers density-dependent regulation that limits the subsequent growth and survival of juveniles (Thorson et al., 2014;Walters et al., 2013).Thus, the reduced resilience and decreased connectivity maintained by spatial contraction and site fidelity prevent rescue effects and spatial recovery across the Snake River (sensu Atlas et al., 2015).Despite ongoing challenges to their widespread recovery, the recovery of some local populations within the Snake River has put pressure on NMFS to potentially downgrade recovery targets of Chinook salmon in the face of uncertain historical baselines (Thurow et al., 2020).
The persistent collapse and nonrecovery of Snake River Chinook salmon continue to challenge management despite effective monitoring of the spatial structure within the watershed (Northwest Fisheries Science Center, 2015).As shown in our simulated results, such localized disturbances to metapopulations with low dispersal rates among dendritic habitat networks, like the Snake River watershed, can lead to prolonged risks of nonrecovery and local collapses (Figure 7; Appendix S1: Table S2).This case study reveals two important insights: (1) broad disturbances to metapopulations with low dispersal can manifest as spatially uneven local impacts, potentially triggering severe spatial contraction and local bottlenecks that dampen widespread recovery, and (2) low dispersal and local density dependence within disturbed metapopulations can undermine management even with data-rich monitoring and informed recovery targets.Despite local recoveries in a few streams, Snake River Chinook salmon remain at risk and have not seen more widespread recovery across the watershed.

DISCUSSION
Interactions between local population dynamics, connectivity, and spatial disturbances can generate difficultto-predict metapopulation recovery patterns.Here, we revealed how these emergent outcomes challenge the long-term management of at-risk metapopulations using modeled simulations supported by case study reviews of three species listed under the US ESA.Our findings provide several insights germane to managing successful recoveries of disturbed metapopulations.We found that spatial complexity in ecological and management scales can result in misleading expectations in spatial recovery regimes.For example, spatial mismatches between management and ecology can mask hidden local collapses within metapopulations (Okamoto et al., 2020;Ying et al., 2011).We also found that the spatial structure of both the habitat network and, most importantly, the disturbance regime itself, can undermine monitoring and management efforts.For example, localized uneven disturbances in habitat networks with low connectivity (e.g., linear coastal networks) can trigger a spatial contraction and perpetuate local bottlenecks that prohibit multiple recovery criteria from being achieved.In other words, metapopulations can act less than the sum of their parts in certain contexts, particularly when localized disturbance regimes trigger feedbacks between local density dependence (which limits compensatory productivity) and low connectivity (which limits external rescues).Furthermore, the case studies we present demonstrate that spatially contracted metapopulations can become vulnerable to new pressures that further compromise their recovery (e.g., sea otter recovery bottlenecks created by great white shark recovery; Moxley et al., 2019;Tinker et al., 2016).Ultimately, these emergent recovery patterns can increase political pressure on managers to adjust recovery criteria or shift target baselines for managing atrisk metapopulations (Bergstrom et al., 2009).Despite political pressures resulting from slow recoveries, all three species reviewed in our case studies maintained more protective listings and both snail kites and sea otters appear to be on the road to recovery (Lafferty & Tinker, 2014;Reichert et al., 2021).
Recent work by Zelnik et al. (2019) demonstrated that disturbed metapopulations exhibited three general recovery regimes-isolated recovery, rescue recovery, or mixed-although a fourth pattern can also emerge: persistent collapse.Our results further suggest that metapopulations under rescue recoveries (i.e., collapsed patches relying upon rescue effects from neighbors) can be sensitive to the configuration of their habitat networks.In particular, metapopulations with large maximum distances between patches (e.g., the two endpoints of a coastal network) were sensitive to unevenly distributed spatial disturbances.In such cases, recovery relies on external rescue effects, but local patches are unlikely to receive that rescue within any short-term time period.Zelnik et al. (2019) found similar results where the effective system size controlled whether metapopulation recovery depended upon internal or external processes.Smaller system sizes increase dispersal and connectivity and allow the metapopulation to be well mixed, which improves recovery rates but can generate new forms of risk by synchronizing local dynamics.Recovery timelines also depend, in part, on whether dispersal is density independent, because density-dependent dispersal can slow the rate of recovery as dispersal only becomes significant as local patches fill up (e.g., California sea otters).Dispersal patterns are shaped by complex behaviors, like collective decision-making and learned migration, that can lead to density dependence, with the directionality either weakening or strengthening metapopulation resilience by altering population synchrony and source-sink dynamics (MacCall et al., 2019;Rogers et al., 2018;Yeakel et al., 2018).Density dependence can maintain population bottlenecks for longer than the assumptions in our model (Einum et al., 1998), suggesting that spatial recovery regimes may be sensitive to a wider array of density-dependent processes.
In practice, many of our model assumptions do not capture the vast diversity of spatially structured populations and their management systems.Overall, our model most directly pertains to understanding metapopulations in their recovery phase, as we modeled spatial disturbances that were curtailed after a single pulsed event.However, many metapopulations are impacted by a mix of pulsed and persistent disturbances (e.g., harvest, land-use development) that can be uniformly or asymmetrically distributed among local populations or even individuals (e.g., size-selective fisheries).While the asymmetric disturbance regimes we explored generally led to worse recoveries, it may be that asymmetries could lead to resilient recoveries if disturbances were aligned with per capita recovery potential, that is, harvesting more productive populations while conserving the less productive (Moore et al., 2021).Furthermore, we designed our model for a semelparous species with a 1-year generation time that neglected important size and age structure dynamics.Studies on the efficacy of spatial protections revealed that, for example, size and age structure can create storage effects (sensu Warner & Chesson, 1985) within protected areas that increase metapopulation resilience to a mix of pulsed, persistent, or stochastic disturbances (Micheli et al., 2012;Nowlis & Roberts, 1999).Our model also lacked behavioral and evolutionary processes that allow for metapopulations to adapt to disturbances (Dunlop et al., 2009;Rogers et al., 2018).For example, adaptive responses to selective harvest can improve spatial conservation, but these responses depend on dispersal distances (Mee et al., 2017).Conversely, some behavioral or local adaptations can induce threshold responses and depensatory dynamics (e.g., Allee effects) that weaken recoveries (Aalto et al., 2019;Lande, 1998;Rogers et al., 2018).Modern risk assessments can account for many of these above dynamics by tailoring the model design to specific life histories, management levers, spatial structure, disturbance regimes, and data availability (Benson et al., 2015;Cadrin et al., 2018;Thorson et al., 2019).Despite these simplifications, our model builds off previous work to further illuminate context-dependent pathways in spatial recoveries (Okamoto et al., 2020;Zelnik et al., 2019).
Multiple scales must be considered by managers designing conservation plans because recovery regimes are governed, in part, by both local and regional processes.For example, planning optimal spatial conservation areas, like marine protected areas, depends partly on how organisms move through habitat networks (Fullerton et al., 2016;Sanchirico & Wilen, 2005).We also found that recovery regimes can generate different ecological outcomes that affect the spatial and temporal provisioning of resources along a path to recovery.Local productivity and dispersal can define the appropriate scale for monitoring to reliably track both local and aggregate risks among metapopulations (Cadrin et al., 2018;Fletcher et al., 2013;Okamoto et al., 2020).If spatial structure is unaccounted for, management improvements may be necessary because spatial mismatches between management and ecology can impose inequitable spatial distribution of sustainable natural resources (Okamoto et al., 2020;Stier et al., 2020).Hence, understanding and managing these disparate recovery regimes is not just of theoretical interest, but has real-world consequences for resource management.For example, hidden local collapses within Pacific herring fisheries that may otherwise be sustainably managed on the aggregate can disproportionately affect human communities that rely on the harvest of local stocks (Okamoto et al., 2020).As well, understanding recovery regimes can help managers predict differences in expected recovery timelines to trigger delisting criteria for at-risk species (Lafferty & Tinker, 2014;Tinker, Gill, et al., 2019).
The spatial and temporal scales governing metapopulation dynamics can exceed or outpace the ability of management to confront ecological change.The policy focus of laws, like the ESA and MMPA in the United States or the Species At Risk Act in Canada, is typically on extinction risks and promoting recovery at the level of species, subspecies, or distinct population segments.Regulatory agencies thus rely upon the bestavailable science to inform listing criteria and recovery plans.Understanding extinction risks and achieving recovery goals for metapopulations may further depend on spatial structure, as illustrated here.When management and policies are nonspatial but ecosystem dynamics are, there can be unrealistic expectations for recovery that challenges management (Stier et al., 2020).For example, the recovery plan for California sea otters listed recovery at 3090 total otters across coastal California, but until 2010 there was limited understanding of the spatial processes underlying their slow recovery (Lafferty & Tinker, 2014).Unexpected recovery dynamics can put political pressure on management to shift recovery baselines, lower criteria for delisting, or authorize targeted disturbances, like harvest, even when species or ecosystems are still at risk (Bergstrom et al., 2009;Pauly, 1995).Furthermore, there can be downstream consequences when conservation fails to consider spatial structure; for example, the MMPA prohibits the take of any sea otters in Alaska (except for limited harvests by coastal Alaska Native peoples), with equivalent protections in British Columbia (Canada), but these policies may inadvertently jeopardize commercial and subsistence shellfish fisheries because of the cascading impacts of sea otter recoveries in these regions (Burt et al., 2020;Salomon et al., 2015).The role of spatial structure in recoveries can be addressed by scientists and policymakers working together to fill data limitations through improved monitoring (see Everglade snail kites) and updating recovery plans accordingly (see Snake River Chinook salmon).
Recovery plans for at-risk metapopulations are challenged by uncertainties that arise from mismatches between ecology and monitoring (Connors et al., 2019).Here, we explored a rich interplay between local and regional processes to illuminate how a few key factors can shape qualitatively disparate recovery patterns.We thus recommend a few heuristics for managing at-risk metapopulations: 1. We suggest that spatial structure be more frequently recognized in listing criteria and recovery plans.Our findings, and that of others (Okamoto et al., 2020;Zelnik et al., 2019) highlight that spatial processes, including dispersal, habitat networks, and disturbances, can lead to many emergent outcomes that slow recovery and challenge management.2. We suggest that scientific monitoring explores and understands the potential role of spatial structure to informing recovery plans (Connors & Atnarko Sockeye Recovery Planning Committee, 2016;Fullerton et al., 2016).Practical tools and technologies such as genomics and genetic stock identification, spatial microchemistry, and telemetry can help to detect and refine the spatial resolution relevant for monitoring that better matches the underlying ecological dynamics (Stier et al., 2020).This would allow assessments to quantify the relative contributions between local and regional processes to metapopulation dynamics to target recovery bottlenecks.Even a qualitative understanding of such processes could help to identify pathways of metapopulation recoveries (Zelnik et al., 2019).3. We suggest that recovery targets be designed to balance local-level and aggregate-level performance metrics, particularly if people rely upon local populations as resources (Okamoto et al., 2020).Recovery plans often aim to improve population abundances on the aggregate, but this strategy may inadvertently expose human communities to risks associated with a lack of access to local resources (Löfqvist et al., 2022).Given that recovery targets are a management construct, both social values and ecological processes can be used to inform these targets.Alternative plans could target an equal distribution of recovery targets across local populations (i.e., equitable utility), prioritize core habitats or source populations (Bradford & Braun, 2021), or maximize the minimum among local populations to diversify the social benefits provided by restoration (i.e., raising the floor in Rawls maximin principle; Rawls, 1971).4. We recommend that management and recovery plans consider the historical and ongoing spatial disturbances as a key determinant of recovery (e.g., harvest allocation in mixed-stock fisheries; Moore et al., 2021).For example, we found that metapopulations subjected to localized and uneven disturbances can experience spatial contraction.In such a case, protecting productive patches and promoting connectivity between local patches might improve recovery.
Conversely, the same metapopulation under a localized and evenly distributed disturbance may experience recruitment declines, suggesting that reducing the overall disturbance (e.g., total harvest) or spreading the disturbance across more subpopulations might improve recovery.
While we are unlikely to fully understand the full suite of social or ecological dynamics underlying many metapopulations, more flexible management frameworks that incorporate the above precautionary provisions may be more robust to this incomplete knowledge.Proactive management using such frameworks could lower the risk of ecological dynamics induced directly or indirectly by a proposed disturbance regime (e.g., cumulative effects from overharvest and habitat disturbance).Given that most anthropogenic impacts on metapopulations are likely to vary across space, such a framework can employ flexible management to better hit moving targets by redefining management goals or reassessing the distribution of risks during the recovery windows for at-risk metapopulations.

F
I G U R E 2 Interplay between local productivity, dispersal between patches, and network complexity shape context dependence of metapopulation recoveries under alternate spatial disturbance regimes.Large gray horizontal arrows in panels (a.1-a.3)highlight how a suite of modulating factors-like productivity, dispersal, and network complexity (high$low; black vertical arrows)-affect recovery dynamics (panels b, c).Panels (b, c) show recovery dynamics for a metapopulation (black line) compromised of two local patches (green and orange lines) under either equal (panels b.1-b.3;where both patches reduced to 50% capacity equally) or unequal disturbance regimes (panels c.1-c.3; with one patch reduced to 0% and the other remaining at 100%).The dashed lines and gray shaded regions within panels (b.2) and (c.2) indicate how the interplay between modulating factors can promote faster or slower recoveries, depending on feedback between productivity, dispersal, and network complexity (gray horizontal arrows).

F
I G U R E 3 Spatial recovery regimes for two example metapopulations (solid black line) composed of 16 local patches (dashed lines) with linear (a, b) and dendritic (c, d) networks impacted by local uneven disturbances in year 50.Metapopulation processes were otherwise identical with low stochasticity, 1% dispersal rates, and the same number of local patches each with identical productivity and capacity.Relative abundance (scaled blue to red by relative abundance at year 100 in panels a, c) measured as the ratio of total abundance N at time t to carrying capacity K.

F
I G U R E 4 Metapopulation recovery rates along gradients of dispersal rates (shown on a logarithmic scale), spatial distribution of disturbances, and habitat configurations including (a) linear, (b) dendritic, (c) star, and (d) grid networks.The lines represent the median recovery rate across all simulated scenarios for a given disturbance regime (9504 total simulated scenarios that tested interplay between heterogeneity in local demography, spatial disturbances, dispersal rates, habitat network, spatial-temporal correlations in stochasticity, and scale of recruitment stochasticity).
After their near-extinction from the North Pacific fur trade from 1750 to 1911, sea otters were protected by international treaty and are classified on the IUCN Red List as endangered, and protected in the US under the Marine Mammal Protection Act (MMPA) of 1972 and, in certain regions (SW Alaska and California), by the US ESA

F
I G U R E 8 Comparison of observed (points) and predicted (lines) recovery patterns between California and Alaska sea otter metapopulations.Years since initiation of recovery corresponds to the years since management began recovery efforts: 1947 for California when a remnant colony in Big Sur became legally protected, and 1969 in Alaska when sea otters were translocated to seven patches.Estimated mean abundances (95% CIs in shaded regions) in California and Alaska from logistic growth models in Tinker et al. (2021) and Tinker, Tomoleoni, et al. (2019), respectively.The horizontal dashed line corresponds to the California recovery target of 3090 sea otters.Panel (a) shows California metapopulation recovery dynamics during 1983-2019 when spatial structure became increasingly recognized.