Identifying Species Conservation Strategies to Reduce Disease‐Associated Declines

Emerging infectious diseases (EIDs) are a salient threat to many animal taxa, causing local and global extinctions, altering communities and ecosystem function. The EID chytridiomycosis is a prominent driver of amphibian declines, which is caused by the fungal pathogen Batrachochytrium dendrobatidis (Bd). To guide conservation policy, we developed a predictive decision‐analytic model that combines empirical knowledge of host‐pathogen metapopulation dynamics with expert judgment regarding effects of management actions, to select from potential conservation strategies. We apply our approach to a boreal toad (Anaxyrus boreas boreas) and Bd system, identifying optimal strategies that balance tradeoffs in maximizing toad population persistence and landscape‐level distribution, while considering costs. The most robust strategy is expected to reduce the decline of toad breeding sites from 53% to 21% over 50 years. Our findings are incorporated into management policy to guide conservation planning. Our online modeling application provides a template for managers of other systems challenged by EIDs.


Introduction
Emerging infectious diseases (EIDs) threaten a variety of species, and may cause mass mortality, increase local and global extinction rates, and alter communities and ecosystem function (Daszak et al. 2000;Whiles et al. 2013;Langwig et al. 2015). EIDs affect a variety of terrestrial, freshwater, and marine ecosystems and often result in spill-back, influencing domestic animal and human populations (Daszak et al. 2000;Fisher et al. 2012;Thompson 2013). EIDs are particularly threatening to declining or endangered species, leaving small, fragmented populations susceptible to stochastic events (Ginsberg et al. 1995;Cully et al. 1997;Fisher et al. 2012). Despite decades of active research on causal pathogens, uncertainties remain about factors influencing host-pathogen dynamics (e.g., Venesky et al. 2014) and potential effects of in situ management actions (e.g., Garner et al. 2016) that impede species conservation efforts (Scheele et al. 2014;Russell et al. 2017). Disease research alone is inadequate to con-serve species postepidemic (Langwig et al. 2015); a major gap remains in translating existing research to appropriate in situ actions (Scheele et al. 2014). What is needed is a process that explicitly and transparently uses empirical and expert knowledge to predict how the variety of management actions we might consider will influence values (i.e., management objectives) associated with the system, and that further outlines a coherent approach to dealing with tradeoffs when making decisions. Decision analysis (Gregory et al. 2012) is such a process; we demonstrate a decision-analytic process for situations in which little empirical knowledge exists on conservation actions, a common condition when dealing with an EID.
The essential elements of a decision-analytic framework include: (1) identifying management objectives, (2) considering alternative actions, (3) developing predictive models to evaluate actions relative to objectives, and (4) identifying the optimal conservation strategy (e.g., via multiobjective optimization). We use this framework to choose among management actions for an amphibian species (boreal toad, Anaxyrus boreas boreas) threatened by a fungal pathogen (Batrachochytrium dendrobatidis [Bd]). Using the best available science and expert judgment (where research did not exist), we developed a dynamic, spatially explicit, host-pathogen metapopulation model to predict outcomes of potential management actions with respect to conservation objectives. Collectively, these objectives aim to conserve the boreal toad populations in the southern Rocky Mountains (SRM) that have experienced dramatic declines over the last two decades, most likely due to Bd (Muths et al. 2003).
The decision-analytic framework is ideally suited for developing species conservation plans, as it explicitly defines objectives, emphasizes creative generation of alternatives, incorporates multiple sources of uncertainty, and provides a transparent approach to dealing with tradeoffs among objectives. Our work informs a revised SRM Boreal Toad Conservation Plan, prepared by multiple state and federal wildlife and land management agencies. The conservation plan defines the problem and identifies specific objectives of maximizing the persistence and number of boreal toad breeding populations across the historic range over a 50-year time horizon, while considering financial costs (Crockett 2017). We use our metapopulation model to forecast landscape-level breeding occurrence for a baseline conservation strategy (no additional actions) and sets of potential actions (strategies) to identify an optimal strategy, given the management objectives. Our practical and synthetic approach demonstrates how conservation agencies can identify the most effective management strategy while incorporating scientific knowledge and uncertainty, and should serve as a guide for other species challenged by EIDs. Conservation practitioners could benefit from adopting our elicitation, modeling, and optimization process to tackle species affected by EIDs.

Species background
The boreal toad is a subspecies of the widely occurring western toad (A. boreas) of North America. Once common in the mountainous regions of southeastern Wyoming, Colorado, and northern New Mexico, boreal toad populations declined noticeably beginning in the late-1970s (Carey et al. 2005), with multiple local extirpations linked to Bd (Muths et al. 2003). The SRM metapopulation may be an evolutionarily significant segment of the boreal toad (Goebel et al. 2009) and belongs to a clade being evaluated by the U.S. Fish and Wildlife Service for listing under the Endangered Species Act (USFWS 2012).
The fungus Bd causes chytridiomycosis through infection of keratinized epithelial cells. Amphibians die when high infection loads disrupt osmotic and electrolyte regulation, causing cardiac arrest (Voyles et al. 2009;Vredenburg et al. 2010). Mortality rates vary by species, life stage, infection intensity, the virulence of the fungal strain, and environmental factors (Van Rooij et al. 2015). Observational studies have documented the influence of environmental features (Beyer et al. 2015) and natural disturbances (Hossack et al. 2013;Roznik et al. 2015) on amphibian-Bd dynamics, but few mitigation strategies have been field tested and no single strategy is appropriate for all species threatened by Bd (Converse et al. 2016;Garner et al. 2016). Collectively, these limited studies form the foundation of knowledge regarding potential in situ strategies available to managers facing Bd-related declines. However, this knowledge alone is not adequate for identifying the most appropriate action to take in any particular situation.

Metapopulation model
We built a metapopulation projection model based on the model proposed by Converse et al. (2016), which uses a multistate occupancy estimation framework (Nichols et al. 2007). Boreal toad ("toad," hereafter) sites include 125 wetland complexes distributed across 11 mountain ranges that have historically supported breeding (Figure 1). In a given year, each site may be in one of the four exclusive states, defined by the status of toads (A) and Bd (B): AB = toad breeding occurs and Bd is present. A0 = toad breeding occurs and Bd is absent. B0 = Bd is present and toad breeding does not occur. 00 = neither toad breeding nor Bd occurs.
Sites may change state from year t to t+1 based on probabilities of colonization (γ ) and extirpation ( ). These dynamic parameters depend on the state of the site at time t. For example, given that toad breeding occurs at a site without Bd (state A0) in year t, A0 is the probability that no breeding will occur the following year; alternatively, γ B A is the probability that Bd will colonize a site at time t+1, given that toad breeding occurred at time t. Thus, the probability that a site changes from state A0 to state AB is γ B A (1 − A0 ), describing the joint probability that Bd colonizes a site where toad breeding occurs at t, and that toads continue to breed at the site in t+1. State transitions can be described via a matrix ( Figure 2).
We estimated colonization and extirpation parameters from monitoring data collected at 82 of the 125 sites within the SRM over a 10-year period (2001-2010) using a dynamic two-species occupancy model that accounts for false-negatives (Table 1; Richmond et al. 2010;Dugger et al. 2016). Parameter estimation details are provided in Appendix S1. Then, we integrated these empirical param- The logit function projects probabilities (θ) onto the real number line, where logit(θ) = log( θ 1−θ ). Notation definition includes: = extirpation probability; σ = spatial colonization scalar; and ξ = baseline (nonspatial) colonization probability of Bd. The mean probability (inverse-logit) for AB , A0 , and B0 = BA is 0.341, 0.040, and 0.036, respectively. See Appendix S1 for details on parameter estimation and Figure S1 for a visualization of each spatial scalar probability relating to colonization of toads and Bd. b This parameter only applies to the colonization probability of Bd.
eters with expert-elicited management effects (described below) to predict toad and Bd dynamics under various conservation strategies.

Spatially explicit colonization
We derived spatially explicit colonization probabilities for toads and Bd as a function of a site's proximity to other occupied sites using our empirical colonization estimates (Converse et al. 2016). We used the Euclidian distance between unoccupied and occupied sites, specified as a Gaussian dispersal function (Equation (1)), as a predictor of colonization. This approach is empirically supported via estimates of toad movements between breeding Figure 2 State transition probabilities for a dynamic host-pathogen metapopulation model, specified using colonization (γ ) and extirpation ( ) parameters for four states, AB, A0, B0, and 00; "A" indicates the presence of boreal toad breeding, "B" indicates the presence of Bd, and "0" indicates an unoccupied site.
sites (E. Muths, unpublished data). For Bd, we considered colonization probability to be a joint function of distance to neighboring Bd-occupied sites and a distanceindependent baseline probability (ξ ; Equation (1)). If ξ is greater than zero, all sites have some probability of being colonized, regardless of their proximity to Bd-occupied sites. Including a positive ξ captures uncertainty about how environmental or anthropogenic factors influence landscape-scale spread of Bd (Venesky et al. 2014;Kolby et al. 2015). The probability of Bd colonization (γ B A s,t = γ B0 s,t ) for site s in year t is the joint probability of colonization as a function of distance (d) to occupied sites in the set J in year t, weighted by the Bd-specific spatial scalar, σ B , and the baseline ξ ( Figure S1) Toad colonization (γ AB s,t and γ A0 s,t ) is modeled in a similar fashion, with a scalar value of σ A0 or σ AB (Table 1), but with no baseline colonization (ξ = 0; colonization of toads is strictly dependent on neighboring occupied sites; Figure S1). We derived spatial scalars using a simulated (calibration) inverse prediction method (Converse et al. 2016). The probability distribution of ξ was estimated via expert elicitation (Table 1; Figure S1).

Management actions and expert elicitation
We considered two types of management actions: (1) those that alter colonization and extirpation probabilities of Bd or toads and (2) toad reintroductions. We used expert elicitation (Martin et al. 2012) to identify potential management actions and quantify the associated effects on toad and Bd dynamics (i.e., extirpation and colonization probabilities). Experts were either members of the Boreal Toad Conservation Team or researchers.
We elicited the effects of actions and their uncertainty using a four-point process to help guard against expert overconfidence (Speirs-Bridge et al. 2010). For each action, experts considered how the action would proportionally change one or more model parameters and independently judged the most likely effect value, lowest and highest possible values, and their confidence that the true value was between their lowest and highest value ( Figure S2). Following a modified Delphi process (Martin et al. 2012), experts initially recorded their four-point values independently, then shared thoughts and experiences with the group, after which they could privately revise their final values. We used the four-point values to specify a probability distribution for each expert and action (Burgman 2005; Appendix S1). We derived an integrated probability distribution for each action, equally weighting all experts (Casella & Berger 2002); these integrated distributions incorporated within and among expert uncertainty. Last, we elicited the value for the distance-independent colonization probability (ξ ) for Bd ( Figure S1) and the relative financial costs of each action.
Toad reintroductions were incorporated in our model assuming typical management practices; harvested wild eggs are reared in hatchery facilities and tadpoles are released in selected unoccupied sites (Muths et al. 2014). Tadpoles are released annually for 10 years, then a new set of reintroduction sites is selected and the process is repeated. Reintroduction sites were initially absent of toads and Bd, and became occupied (e.g., 00 → A0) if reproduction occurred, usually 6-9 years after a cohort is released, corresponding to female age of maturity (Carey et al. 2005). The probability that a single cohort successfully reproduces was elicited from experts and differed depending on whether Bd was present. The probability of reproduction in a given year depends on the number of sequential cohort releases and the timing of potential reproduction overlap among cohorts ( Figure S3). Reintroduction sites were selected using one of the four selection strategies (see Table S1 description).

Management strategies and outcomes
No single action is likely to eliminate Bd or its effects on toad populations (Garner et al. 2016). Using the elicited probability distributions associated with each action, the conservation team developed a set of 35 management strategies which combined disease, habitat, and reintroduction actions (Table S1); some strategies focused on minimizing reproductive failure at toad-occupied sites, while others aimed to increase solar radiation to increase toad persistence and local Bd extinction. All management strategies called for agencies to implement disinfection practices for employees, researchers, and fire-fighting personnel to reduce transmission of Bd. Strategies also differed in the number of reintroductions per year (0, 2, or 4) and how reintroduction sites were chosen.
We evaluated the i = 35 management strategies and a baseline strategy (no additional actions) according to the management objectives of maximizing toad persistence and distribution, while considering financial costs. We measured toad persistence in several ways: using the expected number of active toad breeding sites (ToadSites), the probability of metapopulation extinction (0 breeding sites), and quasi-extinction (<5, <10, and <20 breeding sites) in 50 years. Toad distribution was measured as the expected number of occupied mountain ranges in 50 years (OccMtns). We balanced tradeoffs among objectives to identify optimal strategies across all permutations of objective weights (w o ) for toad persistence, OccMtns, and the relative financial costs among strategies (Costs). We scaled each objective attribute to a mean of zero and unit variance, summed the weighted value for all possible permutations of objective weights, where w o ∈ [0, 1] and {O = 3} {O = 1} w o = 1, and identified the optimal strategy as the one that maximized the weighted value: Optimal strategy w1,w2,w3

Model initialization and implementation
We implemented our metapopulation model in R (R Development Core Team 2016; any use of trade, firm, or product names is for descriptive purposes only and does not imply endorsement by the U.S. Government). We projected all strategies with actions implemented each for 50 years with 4,000 iterations. The initial states of sites (AB, A0, B0, 00) were based on empirical estimates from 2014, which included 55 active toad breeding sites across eight mountain ranges. We provide an online graphical interactive Web-app for researchers and managers to investigate additional strategies using this model, https://borealtoad.shinyapps.io/ SRM-BT-MODEL-Public/.

Expert elicitation actions
We elicited effects of 43 actions from boreal toad experts (Tables 2 and S2). The actions identified as most influential in altering toad breeding or Bd differed for colonization and extirpation processes. For example, experts believed that re-establishing beavers (Castor canadensis) in unoccupied drainages was the best single action to increase toad colonization, while their perceived best action to decrease toad extirpation was to restrict all ground disturbing activities within a specified distance of active breeding sites. Experts were relatively pessimistic about reducing Bd colonization, but their perceived best strategy was to limit visitor use at breeding sites, and their perceived best action to increase Bd extirpation was to remove trees and sedges to increase solar radiation at breeding sites (Table 2).

Predicted outcomes and evaluating strategies
The probability of SRM metapopulation extinction in 50 years was very low, regardless of the management strategy employed (baseline: 0.03 and ࣘ0.02 for all other strategies; Table S3; Figure S4). Quasi-extinction (<5 sites) was moderately high under the baseline strategy (0.09), but was expected to be low under any other strategy (ࣘ0.04; Figure S5). In contrast, quasi-extinction (<10, <20 sites) probabilities were relatively high under the baseline strategy (0.23 and 0.50, respectively) and alternative management strategies had a range of outcomes (<10 sites: 0.03-0.13, <20 sites: 0.15-0.38; Table S3; Figures S6 and S7). After 50 years, the number of toad sites and occupied mountain ranges under the baseline strategy declined by 53% (ToadSites = 25.7 sites, OccMtns = 3.8 mountain ranges; Table S3; Figure 3). Any of the alternative strategies yielded fewer declines, with best case values of 18% decline in OccMtns and 19% decline in ToadSites.
Predicted metapopulation extinction and quasiextinction probabilities (<5, <10 sites) were relatively small and varied little among strategies, making them ineffective at discriminating among strategies. Quasiextinction (<20 sites) probability was more useful at differentiating among strategies and was highly correlated (R 2 = 0.92) with ToadSites, which we used to investigate objective tradeoffs. Several strategies performed well, yielding the highest expected number of occupied toad sites (40-45) distributed among six to seven mountain ranges ( Figure 3); other strategies performed worse, but were less costly.
Model results identified six optimal strategies, across all permutations of objective weights (Figure 4). The baseline strategy was only optimal if the objective weight on Costs was ࣙ0.80. Strategy C.5 (Table S1) was the optimal strategy across the most combinations of objective weights (66%), provided Costs weight was ࣘ0.55. Strategy C.5 included the maximum number of reintroductions and incorporated a number of actions to reduce Bd colonization and protect toad sites, increasing toad colonization and persistence (Table S1). Strategy C.5 was predicted to reduce the expected decline of toad-occupied breeding sites and mountain ranges to 21% (26 to 44 sites) and 17% (3.8 to 6.6 mountain ranges), respectively, relative to the estimated 53% decline in each under the baseline strategy ( Figure S8).

Discussion
Identifying effective strategies is challenging in most species conservation programs because of a lack of clarity about management objectives, an insufficient set of candidate actions, or uncertainty about the species' ecology or effects of actions. EIDs have driven wildlife species and populations to extinction and extirpation, and fungal pathogens are particularly challenging to combat following emergence (Fisher et al. 2012;Langwig et al. 2015). Given the dire situation of global amphibian   Table S2. a Indicates effect differences in the same action when Bd is present and absent, respectively.
declines, effective conservation actions that ameliorate the effects of Bd are needed.
We encourage conservation practitioners in other systems to take a decision-analytic approach to management in the face of EIDs. This process requires that decision makers and stakeholders define specific quantifiable objectives, identify a robust set of actions for further consideration based on experience, expert knowledge or assessing the situation from unique viewpoints; and then develop a predictive modeling framework. The framework integrates available empirical and expert knowledge, including uncertainties about system dynamics and the effectiveness of actions, and then quantitatively evaluates management strategies in terms of identified objectives (Gregory et al. 2012;Converse et al. 2016). Models may be developed initially with expert Figure 3 Predicted outcomes of a baseline (no additional action) strategy and 35 potential management strategies to reduce the impacts of Bd on landscape-level breeding occurrence of boreal toads within the SRM. Outcomes for each action include: the expected number of active toad breeding sites, distributed among mountain ranges, and relative cost. For instance, the baseline strategy (no action, in black) has the lowest cost, but also the lowest expected number of breeding sites and occupied mountain ranges (25.6 and 3.8, respectively). In contrast, Strategy D.4 has relatively high costs, buts is predicted to have the most number of breeding sites and occupied mountain ranges (44.4 and 6.6, respectively).

Figure 4
The optimal strategy across all permutations of objective weights for three fundamental objectives: maximize the predicted number of active toad breeding sites, maximize toad breeding distribution among mountain ranges, and minimize relative costs. The location of each cell corresponds to the weight for the number of toad-occupied mountain ranges (x-axis) and the number of active toad breeding sites in the SRM (y-axis) in 50 years. The number within each cell is the weight for relative costs. Collectively, the three weights sum to 1. The color of the cell represents the optimal action for a given combination of weights.
knowledge, and then updated as empirical information becomes available. Expert knowledge is commonly used in science and conservation when data are lacking, immediate action is warranted, or the complexity of the problem exceeds the feasibility of collecting all necessary data, which all apply to management problems involving EIDs (Martin et al. 2012). This approach directly addresses calls for conservation responses to emerging wildlife diseases (e.g., Scheele et al. 2014;Langwig et al. 2015;Garner et al. 2016), describing a clear and defensible process to determine optimal actions, given existing uncertainties.
For the SRM toad, it is reassuring that our best empirical data suggest that the metapopulation will not go extinct in the near future, but the expected rate of decline in toad breeding sites, if no additional actions are implemented, is substantial. According to our projections, a loss of occupied breeding sites and reduced distribution of toads is inevitable even under optimal management strategies. Increasing the number of reintroduction sites could bolster the population, but this action would be costly and would not affect Bd dynamics (a primary driver of toad population dynamics). To substantially reduce declines, the dynamics of this host-pathogen system would have to shift to a greater degree than we predict will occur under the strategies we considered (e.g., a change of annual toad colonization probability of >0.4 or annual toad extirpation probability in the presence of Bd of <0.05). Our results assume that metapopulation dynamics will remain unchanged, except when influenced by the proposed management strategies. However, other factors may ameliorate or exacerbate the challenges faced by toads, including demographic stochasticity in small populations (Lande 1993); temporal variation in hydrological conditions (Amburgey et al. 2012); or rapid selection of adaptive tolerance to Bd (Converse et al. 2016;Savage & Zamudio 2016).
We suggest that researchers need to focus on in situ experiments in natural settings to test how management actions impact host-pathogen dynamics. Such studies are currently lacking and are sorely needed (Garner et al. 2016). In the short term, effective actions identified here will be prioritized in conservation efforts. Implementation of one of the optimal strategies in conjunction with continued toad-Bd monitoring will provide a means to evaluate the effect of the strategy on metapopulation dynamics. These experiments can be guided by a decisionanalytic process, which highlights the uncertainties that most impede decision making (Runge et al. 2011). As sequential monitoring and decision making occur, adaptive management will help decision makers learn about those components of the system that are directly pertinent to making better decisions (Williams 2011).
The future remains uncertain for many species threatened by EIDs. Addressing these serious conservation challenges requires an approach that synthetically and transparently applies our best understanding of system function to the problem of determining how to act to conserve species. We present a framework for integration of empirical and expert knowledge in conservation decision making, which could be adapted to other conservation efforts focused on assessing in situ management strategies for host-pathogen systems. Importantly, our framework encourages the incorporation of empirical information collected after strategies are implemented, thus progressively improving decisions and, hopefully, conservation outcomes.

Acknowledgments
We are incredibly grateful to all the members of Boreal Toad Conservation Team. This is contribution number 593 of the U.S. Geological Survey Amphibian Research and Monitoring Initiative (ARMI).

Supporting Information
Additional Supporting Information may be found in the online version of this article at the publisher's web site: Appendix S1. Additional information on the empirical metapopulation parameter estimates and expert elicitation. Table A1. Model selection results for the multistate (conditional two species) dynamic occupancy model fit to boreal toad breeding (A) and Bd (B) detection data collected in the SRM from 2001 to 2010. Table A2. Model-averaged parameter estimates from a dynamic two-species occupancy model for breeding boreal toad breeding (A) and Bd (B) in the SRM. Figure A1. Two triangular distributions using an expert's lowest, highest, and most likely values, along with their certainty that the true value lies between the lowest and highest values. Figure A2. An example of an expert elicitation using the Beta distribution. Figure A3. An example of an expert elicitation using the truncated normal distribution. Figure A4. An example of an expert elicitation using the triangular dsistribution. Figure S1. Boreal toad and Bd distance-dependent colonization expressed via Equation 1 (A). Figure S2. An example how an expert's four elicitation points are used to define a probability distribution. Figure S3. The cumulative probability of successful reintroduction (reproduction occurred) if sequential cohorts of tadpoles are released at a site starting in years 1, 11, 21, 31, and 41 and the probability of a single cohort reproducing is 0.4. Figure S4. Predicted outcomes of 36 strategies to improve breeding boreal toad occurrence over 50 years, including metapopulation extinction probability (no sites with toad reproduction), the expected number of mountain ranges with toad reproduction, and the relative financial costs of each strategy. Figure S5. Predicted outcomes of 36 strategies to improve breeding boreal toad occurrence over 50 years, including quasi-metapopulation extinction probability (<5 active toad breeding sites), the expected number of mountain ranges with toad reproduction, and the relative financial costs of each strategy. Figure S6. Predicted outcomes of 36 strategies to improve breeding boreal toad occurrence over 50 years, including quasi-metapopulation extinction probability (<10 active toad breeding sites), the expected number of mountain ranges with toad reproduction, and the relative financial costs of each strategy. Figure S7. Predicted outcomes of 36 strategies to improve breeding boreal toad occurrence over 50 years, including quasi-metapopulation extinction probability (<20 active toad breeding sites), the expected number of mountain ranges with toad reproduction, and the relative financial costs of each strategy. Figure S8. Animated time-series of the predicted distributions for the number of breeding boreal toad sites and expected number of occupied mountain ranges under the baseline and C.5 strategies. Table S1. Description of the management strategies (combined actions) evaluated using the metapopulation projection model. Table S2. All considered management actions and elicited proportional mean changes in model parameters. Table S3. Predicted outcomes of boreal toad breeding occurrence from a meta-population model that integrates empirical estimates and expert-elicited management action effects.