Genetic signatures of evolutionary rescue by a selective sweep

One of the most useful models in population genetics is that of a selective sweep and the consequent hitch-hiking of linked neutral alleles. While variations on this model typically assume constant population size, many instances of strong selection and rapid adaptation in nature may co-occur with complex demography. Here we extend the hitch-hiking model to evolutionary rescue, where adaptation and demography not only co-occur but are intimately entwined. Our results show how this feedback between demography and evolution determines – and restricts – the genetic signatures of evolutionary rescue, and how these differ from the signatures of sweeps in populations of constant size. In particular, we find rescue to harden sweeps from standing variance or new mutation (but not from migration), reduce nucleotide diversity both at the selected site and genome-wide, and increase the range of observed Tajima’s D values. For a given rate of population decline, the feedback between demography and evolution makes all of these differences more dramatic under weaker selection, where bottlenecks are prolonged. Nevertheless, it is likely difficult to infer the co-incident timing of the sweep and bottleneck from these simple signatures, never-mind a feedback between them. Temporal samples spanning contemporary rescue events may offer one way forward.


Introduction
The simple models used to predict the genetic signatures of selective sweeps have been incredibly helpful in understanding and identifying population genetic signals of adaptation (e.g., Maynard Smith and Haigh, 1974;Kaplan et al., 1989;reviewed in Stephan, 2019). These models are usually based on constant-sized, Wright-Fisher populations. Meanwhile, many instances of adaptation -and thus selective sweeps -in nature will co-occur with complex demography. In fact, many of the most well-known examples of selective sweeps have arisen following a rather extreme and sudden change in the environment (e.g., after the application of insecticides, Sedghifar et al., 2016, or antimalarial drugs, Nair et al., 2003, which could have simultaneously imposed sharp demographic declines. Attempts to capture such complex demographic scenarios typically impose qualitatively appropriate changes in population size (e.g., Hermisson and Pennings, 2005). Indeed, a number of studies have explored the genetic signatures of selective sweeps during demographic bottlenecks (e.g., Innan and Kim, 2004;Teshima et al., 2006;Wilson et al., 2014). However, these demographies are nearly always chosen in the absence of an explicit population model and independently of evolution.
Here we model selective sweeps in a scenario where demography and adaptive evolution are not independent. In particular we model an instance of evolutionary rescue (Gomulkiewicz and Holt, 1995;reviewed in Bell, 2017), where a sudden environmental change causes population decline that is reverted by a selective sweep. Under this framework, rescue is a simultaneous demographic bottleneck and selective sweep, where each affects the other. First, because the mean absolute fitness of the population changes with the beneficial 1 mmosmond@ucdavis.edu 1 allele's frequency, the depth and duration of the bottleneck depends on the dynamics of the selective sweep, i.e., evolution affects demography. Second, the probability that the beneficial allele establishes depends on the rate of population decline, i.e., demography affects evolution. Together, this feedback between demography and evolution restricts the range of dynamics that are possible, and therefore also restricts the range of genetic signatures we should expect to observe. Our goal here is to describe the range of genetic signatures that (this model of) evolutionary rescue allows, to help elucidate how rescue may obscure inferences of past selection and demography and to identify patterns that could be used to infer rescue in nature.
Most theory on evolutionary rescue to date (reviewed in Alexander et al., 2014) has focused on the probability of rescue. Recently, however, some attention has been given to the dynamics of population size (Orr and Unckless, 2014), the probability of soft sweeps (Wilson et al., 2017), and the genetic basis of adaptation (Osmond et al., 2019) given rescue in haploid or asexual populations. Here we extend this line of thinking to three modes of rescue in diploid, sexual populations, and use coalescent theory and simulations of whole chromosomes to examine the genetic signatures at linked, neutral loci. Our focus is on three common genetic signatures: the number of unique lineages of the beneficial allele that establish (i.e., the softness of the sweep), the pattern of nucleotide diversity around the selected site (i.e., the dip in diversity), and the pattern of Tajima's D around the selected site (i.e., skews in the site-frequency spectrum). We explore three modes of rescue, where the beneficial allele arises from either standing genetic variance, recurrent de novo mutation, or migration. For each mode of rescue we derive the expected forward-time dynamics and resulting genetic signatures and compare these to predictions for populations of constant size as well as to individualbased simulations. Qualitatively, we find that rescue causes faster, harder sweeps that produce wider, deeper dips in diversity and more extreme values of Tajima's D. Due to the feedback between demography and evolution, the effect of rescue on the signatures of selective sweeps, relative to the signatures in populations of constant size, becomes more pronounced as the selection coefficient, and thus the probability of rescue, gets smaller.

Methods and results
Data availability statement Code used to derive and plot all results presented below (Python scripts and Mathematica notebook) is available at https://github.com/mmosmond/rescueCoalescent.git under an MIT licence.

Deterministic trajectories
Consider a population of size N (t), with a beneficial allele, A, at frequency p(t) and an ancestral allele, a, at frequency q(t) = 1 − p(t). Assume non-overlapping generations and let W AA , W Aa , and W aa be the absolute fitness (expected number of offspring) of each genotype. Then with random mating the expected change in allele frequency (equation 5.2.13 in Crow and Kimura, 1970) and population size in one generation is with W (t) the population mean fitness.
Here we are interested in the scenario where a population composed of primarily aa genotypes is declining at some rate d, i.e., W aa = 1 − d. The beneficial allele, A, is then assumed to act multiplicatively with the fitness of the ancestral background, such that W Aa = (1 − d)(1 + hs) and W AA = (1 − d)(1 + s). Throughout the text we assume additivity at the selected locus, h = 1/2, such that with weak selection the allele frequency (c.f., equation 5.3.12 in Crow and Kimura, 1970) and population dynamics can be approximated by . (2) These equations make the feedback between evolution and demography in this model clear. The allele frequency dynamics depend only on relative fitnesses (Equation 1), W i /W (t), and thus can depend on demography only through initial allele frequency (Equation 2). In contrast, demography depends on mean fitness (Equation 1) and thus is strongly influenced by changes in allele frequency (Equation 2).

Conditioning on rescue
Above we have considered only the deterministic dynamics. We are, however, primarily concerned with a stochastic event -the establishment of a rare beneficial allele. In rescue we are particularly interested in only those instances where the beneficial allele establishes (otherwise the population is extinct), which creates a bias away from the deterministic trajectory. As shown in Maynard Smith (1971) (see Orr and Unckless, 2014, for an application to evolutionary rescue), we can approximate this bias by replacing the true initial allele frequency, p(0), in the deterministic predictions (Equations 1-2) by this true value divided by the probability of establishment. In essence, dividing initial allele frequencies by the probability of establishment implies that an allele escaping random loss when rare will initially increase much faster than expected before settling into its deterministic trajectory.
Assuming alleles do not interact (i.e., a branching process), a single copy of an allele that is expected to leave 1 + copies in the next generation with a variance of σ 2 has an establishment probability of (Allen, 2010, p. 172, see also Feller, 1951 For example, in a Wright-Fisher population of size N a rare allele with selection advantage s is expected to leave 1 + s copies, with a variance of 1 + s + O(1/N ). Thus = s and, in a large population with weak selection, σ 2 ≈ 1, such that the probability of establishment is roughly 1 − e 2s (Fisher, 1999, p. 80) or nearly 2s (Haldane, 1927). In our case, = W Aa − 1 = (1 + hs)(1 − d) − 1 ≈ hs − d and σ 2 depends on the particular choice of lifecycle (e.g., see Supplementary text: Probability of establishment and effective population size in the simulated lifecycle). When h = 1 and σ 2 = 1 we recover the probability of establishment in an exponentially growing or declining population, ρ ≈ 2(s − d) (Otto and Whitlock, 1997). As mentioned above, the change in allele frequency (Equation 1), i.e., evolution, depends only on relative fitness and therefore is not influenced by demography. However, now we see that conditioning on rescue will cause the initial allele frequency, and hence evolution, to depend on absolute fitness (through ), and hence on demography. For example, the faster the rate of initial decline, d, the lower the probability of establishment (because the growth rate of the heterozygote is lower). This causes larger effective initial allele frequencies -and hence shorter selective sweeps -in populations that are rescued from more precipitous declines.

Simulations
The simulation details are described in full in Supplementary text: Simulation details. Briefly, the lifecycle described in Supplementary text: Simulated lifecycle, with the addition of a hard carrying capacity at N (0), was simulated in SLiM  with tree-sequence recording . We simulated a 20 Mb segment of a chromosome, with all but one of the center loci neutral, with a per base recombination rate of r bp = 2 × 10 −8 and per base mutation rate at neutral loci of U = 6 × 10 −9 (both inspired by Drosophila estimates). A population was considered rescued when the beneficial mutation was fixed and the population size had recovered to it's initial size, N (0). Pairwise nucleodtide diversity and Tajima's D were calculated directly from the tree sequences using msprime (Kelleher et al., 2016) and tskit (Kelleher et al., 2018).

Forward-time results
Below we derive the effective initial allele frequency given rescue (allowing us to predict the forward-time dynamics; see Conditioning on rescue) and the softness of the sweep given rescue, for three types of rescue: from standing genetic variation, recurrent de novo mutation, or migration. To do so we must first compute the probability of rescue, which we also report. Many similar results have been previously derived: Orr and Unckless (2008) give the probability of rescue from standing genetic variation (haploid) or de novo mutation (haploid and diploid); Orr and Unckless (2014) give the effective initial allele frequencies given rescue from standing genetic variation or de novo mutation in a haploid model; and Wilson et al. (2017) give the probability of a soft sweep given rescue from de novo mutation in a haploid model. Nevertheless, we include the diploid versions of these results below for completeness, and because they are needed to describe the genetic signatures of rescue that follow.

Rescue from standing genetic variation (SGV)
We first consider rescue from genetic variation that is present at the time of the environmental change, ignoring further mutations. To avoid complicating the presentation here we assume the initial number of beneficial alleles is given; treating the initial number as a random variable requires conditioning the distribution on a successful sweep (Hermisson and Pennings, 2017).
Given there are initially κ N (0) copies of the beneficial allele, the number that establish is roughly binomially distributed with κ trials and success probability ρ. The probability of rescue is the probability that at least one establishes, Conditioning on rescue, it is therefore as if the true initial allele frequency, κ/(2N (0)), is inflated by a factor, 1/P SGV rescue , i.e., the beneficial allele quickly jumps to frequency (c.f., equation S1.4 in Orr and Unckless, 2014) after which it follows its deterministic trajectory. This same conditioning applies in a population of constant size (i.e., d = 0). As the decline rate, d, increases, the probability a copy establishes, ρ, and hence the probability of rescue, P SGV rescue , declines, making the conditioning stronger. This causes selective sweeps to get started faster as d increases, implying that, for a given s, rescue sweeps are faster than those in populations of constant size. Figure 1 compares our numerical (Equation 1) and analytical (Equations 2) approximations against individual-based simulations. This shows that the numerical predictions do a reasonably good job in capturing the median simulation dynamics (we use the median as taking the mean of the non-linear trajectories distorts their shape), especially with larger selection coefficients and more initial copies of the beneficial allele. For smaller selection coefficients and fewer initial copies of the beneficial allele we overestimate allele frequency and population size (e.g., panel D). The reasons for this are two-fold. On the one hand we tend to underestimate the probability of establishment with small selection coefficients; beneficial alleles with smaller selection coefficients can exist at low numbers -before securing establishment -for a longer period of time (Maruyama and Kimura, 1974, e.g., with s = 0.13 and d = 0.05 we have W Aa ≈ 1.01, so that such an allele can persist for up to ∼ (W Aa − 1) −1 ≈ 100 generations before establishing or going extinct, Desai and Fisher, 2007). The longer it takes a beneficial allele to establish the fewer copies of the ancestral allele that exist when it does (in the case of rescue, d > 0), increasing the chance the beneficial allele will experience some selection as a homozygote while establishing (where it has roughly twice the probability of survival when d is small given h = 1/2). This is true even when we start with only a single copy of the beneficial allele (as opposed to a model with constant population size, d = 0, where the ancestral allele only decreases in frequency if the beneficial allele increases in number), although this effect is strengthened when there are initially more copies as some copies may reach substantial frequencies in the meantime, causing even greater beneficial allele frequencies during establishment. On the other hand, conditioning on rescue has a larger effect (on effective initial allele frequencies) when rescue is unlikely. This amplifies our underestimates of the probability of establishment discussed above when there are initially fewer copies of the beneficial allele. The analytical approximations, assuming weak selection, provide yet larger overestimates of allele frequency and population size (due to the relatively large selection coefficients used, as expected in populations maladapted enough to require evolutionary rescue).
Note that when the expected number of establishing copies, ρκ, is small, the probability of rescue is roughly ρκ, so that the effective initial allele frequency, p SGV 0|rescue , is independent of the initial number of copies, κ, implying that rescue tends to occur by a hard selective sweep, i.e., only one of the initial copies establishes. For larger values of κ the effective initial frequency is not independent of κ and rescue can occur by a soft selective sweep (Hermisson and Pennings, 2005), where multiple initial copies establish. The probability that multiple copies establish is P SGV rescue − κρ(1 − ρ) κ−1 , where κρ(1 − ρ) κ−1 is the probability of a hard sweep. Given rescue occurs from standing genetic variation, the probability it is due to a soft sweep is therefore In our two examples above this is 0 (when κ = 1) and ≈ 1 (when κ = 100). Between these two extremes we find Equation 6 to provide reasonable estimates for small κ or large s, but to underestimate the probability of a soft sweep otherwise (Figure 2A), when beneficial alleles can persist at low numbers long enough to establish with some non-negligible probability of experiencing some selection as homozygotes (given d > 0). More generally, as mentioned above, the number that establish is binomially distributed, and dividing this by the probability of rescue then provides the distribution given rescue, which has expectation E[number of copies that establish|rescue] = κρ P SGV rescue . (7) Equation 7 also provides reasonable estimates for small κ or large s ( Figure 2B).

Rescue by de novo mutation (DNM)
When there are few copies of the beneficial allele at the time of environmental change rescue may depend on mutations arising de novo at the selected site during population decline. To predict the allele frequency and population size dynamics in this scenario we then need to derive the waiting time until the first successful rescue mutation. The first successful rescue mutation arrives according to a time-inhomogeneous Poisson process with rate, λ(t) = 2N (t)uρ, where 2N (t) = 2N (0)e −dt describes the decline in the number of copies of the ancestral allele. Thus the probability that a rescue mutation has established by time T (i.e., the cumulative distribution of the waiting time) is F (T ) = 1 − e − T 0 λ(t)dt . Taking the limit as time goes to infinity then gives the probability of rescue (c.f., equation 10 in Orr and Unckless, 2008) Following Orr and Unckless (2014), taking the derivative of F (T ) and dividing by the probability of rescue gives the probability distribution function for the arrival time of the first establishing rescue mutation given rescue, f (t). While the first establishing mutation is still rare it will grow exponentially at rate ρ/2, and conditioned on its establishment will very quickly reach 1/ρ copies. Integrating over the possible arrival times then gives the expected number of copies of this successful mutation at time t since the environmental change, ∞ 0 (e (t−τ )ρ/2 /ρ)f (τ )dτ . Dividing by its expected size at time t, exp(ρt/2), and the total number of alleles at the time of environmental change, 2N (0), it is therefore as if the successful mutation was present at the time of environmental change, with frequency  where Γ(z) is the gamma function (equation 6.1.1 in Abramowitz and Stegun, 1972) and Γ(a, x) is the incomplete gamma function (equation 6.5.3 in Abramowitz and Stegun, 1972). The factor 1/ρ > 1 increases the effective initial frequency, because we have conditioned on establishment, while the last factor decreases the effective initial frequency, because we must wait for the mutation to arise. When the expected number of rescue mutations, 2uN (0)ρ/d, is small this expected number cancels out and the last factor becomes approximately d/(d + ρ/2), which is independent of mutational input (Orr and Unckless, 2014). That is, conditioning on unlikely rescue, rescue mutations arise earlier in populations that decline faster. In a population of constant size, N = N (0), mutations arrive at a time-homogeneous rate, λ = 2N (0)uρ and the probability distribution of the waiting time until the first successful mutation is a simple exponential, f (t) = λe −λt . Integrating the expected number of copies of the allele over the waiting times shows that the waiting time factor in a population of roughly constant size is 4N (0)u/(1 + 4N (0)u), which, in constrast to rescue, depends strongly on mutational input but is independent of establishment probability, ρ.
Figure 3 compares our numerical (Equation 1) and analytical (Equations 2) approximations against individual-based simulations. As with rescue from standing genetic variance ( Figure 1), with small selection coefficients we tend to overestimate allele frequencies and population sizes by ignoring beneficial homozygotes and thus underestimating of the probability of establishment. This is an even bigger issue with rescue by de novo mutation, where a weakly selected beneficial allele can not only exist at low numbers for a large number of generations before establishing, but can also arise when the number of ancestral alleles is already considerably lowered (especially so when the mutation rate is small). Our predictions do better with larger selection coefficients, where beneficial alleles quickly establish or go extinct.  We can also calculate the probability of a soft sweep from recurrent mutation. Taking into account the beneficial alleles present at time t, the rate at which additional copies arise and establish is λ(t) = 2N (t)q(t)uρ(t), where 2N (t)q(t) is the number of ancestral alleles at time t and ρ(t) is the probability of establishment, which changes with allele frequency because this influences the genotypes the new allele experiences. Thus the number of mutations that arise and fix is Poisson with rate ∞ 0 λ(t)dt, allowing us to write down an equation for the probability of a soft sweep. To gain intuition we make a very rough approximation, assuming q(t) ≈ 1 while mutations are arriving (i.e., when N (t) is still large), so that ρ(t) ≈ ρ and we get the same Poisson rate we derived above for the first successful mutation, ∞ 0 λ(t)dt ≈ 2N (0)uρ/d, providing us with the probability distribution for the number of mutations that establish. Dividing the expected number of establishing mutations by the probability of rescue, the expected number that establish given rescue is This is analogous to the result in a model with haploid selection (c.f., equation 7 in Wilson et al., 2017). Ignoring density-dependence, our approximation will underestimate the number of establishing mutations when h = 0 since selection in heterozygotes will prolong the persistence of the a allele (creating more opportunity for mutation) and establishment probabilities will rise with the frequency of the beneficial A allele (because more AA genotypes are then created, as discussed in the case of standing genetic variation above). At the same time, if the carrying capacity is reached before fixation, our simple form of densitydependence will introduce additional genetic drift and hence tend to reduce the number of mutations that establish. In the end we find our rough approximation to underestimate the number of establishing mutations ( Figure 4B), suggesting that our underestimate of the probability of establishment has a larger effect than the excess drift brought about by the carrying capacity for these parameter values.
With these same approximations the probability of a soft selective sweep given rescue (i.e., the probability more than one copy establishes) is as in a haploid population (equation 8 in Wilson et al., 2017). As with the expected number of establishing copies (Equation 10), we see this approximation is an underestimate in diploid populations ( Figure 4A). On a related note, Wilson et al. (2017) reasoned that, because soft sweeps are expected when rescue is likely while hard sweeps are expected when rescue is rare, population bottlenecks will tend to be more extreme when rescue occurs by a hard selective sweep (and thus it might be easier to detect soft sweeps from patterns at linked neutral loci, as bottlenecks could hide the signal). Here we show the importance of conditioning on rescue, which roughly equalizes the bottleneck sizes across scenarios with very different probabilities of rescue (e.g., compare Figure 3A and B, or C and D, where the probability of rescue differs by an order of magnitude), potentially making hard sweeps easier to detect due to their greater effect on local gene genealogies (see below).

Rescue by migrant alleles (MIG)
In Supplementary text: Rescue by migrant alleles (MIG) we derive approximations for allele frequency and population size dynamics given rescue by migrant alleles entering the population at a constant per generation rate m. There are close similarities here with rescue by mutation, and so we largely relegate this case to the appendix. However, with migration as opposed to mutation, the rate at which new copies of the beneficial allele arise does not decline with the number of ancestral alleles. This increases the probability that beneficial alleles will establish when the population is closer to extinction, when compared to rescue by mutation, especially under smaller migration rates. This effect is amplified by the fact that the probability of establishment also increases as the number of ancestral alleles declines. Under sufficiently small migration rates, simulations show that the beneficial allele starts to sweep later but increases in frequency faster than the deterministic expectation ( Figure S6). Here we enter a different regime, which we do not explore. We wait until the coalescent has been introduced to explore soft sweeps under migration (see The number of successful migrants).

The structured coalescent
To explore the genetic patterns created by evolutionary rescue we next flip our perspective and think backwards in time, starting from a random sample of chromosomes at the time the beneficial allele fixes. Focusing on a neutral locus that is recombination distance r from the selected site, we are interested in calculating the rate of coalescence, the rates of recombination and mutation off the selected background, and the rate of migration out of the population. If our sample of alleles has k distinct ancestors on the selected background τ generations before fixation, these rates are approximately (table 1 in Hudson and Kaplan, 1988;equation 16 in Pennings and Hermisson, 2006a; see Supplementary text: Deriving the structured coalescent for derivations) where p (τ ) = p(T −τ ), N (τ ) = N (T −τ ), and N e (τ ) = N e (T −τ ) are the allele frequency, census population size, and effective population size, respectively, τ generations before fixation, with fixation occurring in generation T . With slow changes in population size, N (t − 1) ≈ N (t) and the mean number of gametes contributed to the next generation by each gamete in the current generation is 2. The inbreeding (and variance) effective population size, N e (t), is then roughly (4N (t) − 2)/(σ 2 + 2) (equation 7.6.4.3 in Crow and Kimura, 1970), where σ 2 is the variance in the per capita number of gametes contributed to the next generation. Therefore, in a large population, N (t) 1, the ratio of the effective size to the census size is roughly where σ 2 depends on the particular lifecycle (e.g., see Supplementary text: Probability of establishment and effective population size in the simulated lifecycle).

The number of successful migrants
Before moving on, note that Equation 12 (together with Equation 13 when σ 2 is constant) shows that the per generation probability of migration and coalescence depend on population size and beneficial allele frequency in the same way. This similar form implies that the relative rates at which lineages coalesce and migrate at the selected site does not depend on the population size and allele frequency. Pennings and Hermisson (2006a) used this fact to show that, in an ideal population of constant size, the number of unique migrant haplotypes contributing to a present day sample, as well as their proportions, is described by Ewens' sampling formula (pages 334ff in Ewens, 2004) when we replace θ with 2m. Powerfully, this results holds even in non-ideal populations of changing size (as briefly noted by Pennings and Hermisson, 2006a, p. 1081-1082 -including during evolutionary rescue -as long as the relationship between the effective and census population sizes remains the same (i.e., N e (t)/N (t) is constant; we now replace θ with 2mN e /N in Ewens' sampling formula). Thus the softness of a sweep from migration depends only on the migration rate and variance in gamete numbers (σ 2 , which determines N e /N ; Equation 13), and is the same during rescue as it is in a population of constant size. The analogous result for rescue by de novo mutation does not hold (as it does for a population of constant size, Pennings and Hermisson, 2006a), since the rate of mutation is not inversely proportional to population size (Equation 12).
Here we use just two properties of Ewens' sampling formula, the expected number of unique migrants among a sample of size n (page 336 in Ewens, 2004) and the probability this is more than two is (equation 10.9 in Ewens, 2004) Figure 5 shows that these formulas perform very well, even when we sample the entire population (n = 2N (0)).

The timing of coalescence
We now use Equation 12 to calculate the probability that the most recent event is τ generations before fixation and is either coalescence, recombination, mutation, or migration. Letting i, j ∈ {coal, rec, mut, mig}, the probability that i is the most recent event, and occurs τ generations before fixation, is (c.f., equation 6 in Pennings and Hermisson, 2006b) i.e., the waiting time for an inhomogeneous exponential random variable. The approximation assumes the p i (k, τ ) are small enough such that at most one event happens each generation, with small probability, and the changes in the p i (k, τ ) from one generation to the next are small enough that we can approximate a sum across τ with an integral. As a technical aside, to speed computation we analytically solve the integrals and then numerically evaluate Equation 16 under the assumption of weak selection and exponential population growth (Equation 2); in contrast, the simulations impose a hard carrying capacity of N (0), creating a discrepancy when the bottleneck is finished long before the sweep (d s 1). The fixation time, T , is approximated by solving p(T ) = 1 − 1/[2N (0)] using Equation 2, and thus also assumes weak selection. Figures 6-7 and S7 show the probability and timing of the relevant coalescent events for a sample of size 2 at a linked neutral locus, P i (2, τ ), for evolutionary rescue from standing genetic variation, recurrent mutation, and migration. It compares these rescue scenarios to selective sweeps in populations of constant size (d = 0). From this we can make a number of observations: 1) the bottleneck during rescue pushes coalescent times towards the present, so that the distributions of coalescence and recombination times overlap more, 2) the bottleneck during rescue increases the overall probability of coalescence, which reduces the probability of recombination or mutation off the sweep (compare areas under broken curves in Figures 6-7), 3) the bottleneck during rescue pushes migration times towards the present and increase the overall probability of migration (because a larger proportion of the population then descends from a migrant), and 4) the difference between the rescue model and the constant population size model is larger under weaker selection. This latter point nicely illustrates the coupling between demography and evolution in rescue; while weaker selection creates a slower sweep and hence more time for recombination off, it also slows population recovery in the case of rescue, leading to longer and deeper bottlenecks that counteract the additional time provided for recombination.

Genetic signatures at linked neutral loci
We now use the P i (k, τ ) to describe patterns of genetic variation at linked neutral loci in a random sample of chromosomes at (or not long after) the time of fixation. Under rescue from standing genetic variance, we assume each of the κ initial copies has independently arose via mutation in the recent past, which assumes the allele was sufficiently deleterious before the environmental change (c.f., Prezeworski et al., 2005). We neglect the migration case as this requires a number of assumptions about the history of the metapopulation, e.g., how and when the sweeps occurred in the neighbouring patches, historical migration rates, etc.

Genetic diversity
One classic pattern of genetic variation produced by a selective sweep is a dip in genetic diversity around the selected site (Maynard Smith and Haigh, 1974;Kaplan et al., 1989). Here we consider the average number of nucleotide differences between two randomly sampled sequences, π (Tajima, 1983), focusing on sequences of length 1 (i.e., heterozygosity).
We first consider our expectation for π at a site that is far enough away from the selected site to be unaffected by the sweep; this provides us with an expectation for the genome-wide average, which is determined by the mutation rate at neutral loci and the population bottleneck. In particular, our expectation for π at such a site in a population of constant effective size, N e , is simply θ = 4N e U (Watterson, 1975), with U the per base per generation mutation rate at neutral loci. Ignoring neutral mutation input during the bottleneck, the π at a sufficiently loosely linked site is this neutral expectation times the probability a sample of size two does not coalesce during the period of interest (in our case, from the time we take the sample at the time of fixation, t = T , until the time of environmental change, t = 0). This is (c.f., equation 4 in Slatkin andHudson, 1991 andequation 7 in Griffiths andTavare, 1994) We next consider sites that are more closely linked to the selected locus, and are thus directly affected by the selective sweep. To keep the analysis simple, we assume that if one of the sampled alleles recombines or mutates off the beneficial background before the two samples coalesce then it is as if both samples were on the ancestral background from the start and therefore coalesce with each other as if they were at an unlinked locus (Equation 17). This assumption will be appropriate with large recombination and mutation rates, where both samples will quickly recombine or mutate off the sweep, and with small recombination and mutation rates, where coalescence will typically occur first. At moderate recombination and mutation rates we will tend to underestimate diversity as the second recombination or mutation event will then take place further back in the past, restricting the time over which the samples can coalesce. Finally, there is also the possibility of no events occurring in the history of the sample during the sweep; if the sweep arose from mutation (or migration) or from a single copy of the beneficial allele (κ = 1) then the sample must coalesce, otherwise (i.e., if we start with more than one beneficial copy, κ > 1), under our assumption that each copy of the beneficial allele among the standing variance has a unique mutational origin, the samples are independent draws from a neutral population. Ignoring neutral mutations during the sweep, a simple approximation for π at any location in the genome is then where P off (k, T ) = T 0 [P rec (k, τ ) + P mut (k, τ )]dτ is the probability of recombination or mutation before coalescence during the sweep, P ∅ (k, T ) = exp − j T 0 p j (k, τ )dτ is the probability that no events have occurred in the history of the sample during the sweep (j ∈ {rec, mut, coal}), and δ κ>1 is 1 if the sweep arose from more than one copy of the beneficial allele (κ > 1) and 0 otherwise. Figures S2-S3 compare our predictions of π after evolutionary rescue against both simulations and the constant population size scenario (d = 0). While our predictions qualitatively match simulations, the tendency of our deterministic approximations to overestimate population size when the probability of rescue is small (Figures 1 and 3) causes us to overestimate diversity in these cases. To correct for this, in Figures  8-9 we replace our prediction for E[π|unlinked] (Equation 17) in the rescue scenario with the observed mean diversity level (as a technical aside, because we only simulate a 40cM chromosome, much of which is affected by the sweep, we do not use sites within 5cM of the selected locus in our average; including these sites would have very little effect in larger genomes). A very similar result was achieved by converting the observed population sizes to effective population sizes and computing the sum in Equation 17 (results not shown). Using the observed mean diversity level is justified by the fact that genome-wide diversity can be measured directly from data and is highly variable across populations (Tajima, 1983). In fact, unless the population was sampled both before and after the selective sweep (or we have good estimates of its mutation rate and long-term effective population size), the amount of background pairwise diversity tells us very little about recent population size changes and all the information is contained in relative diversity (the diversity in a window divided by genome-wide diversity). Figures S4-S5 show that our predictions of relative diversity (E[π]/E[π|unlinked) closely match that observed in simulations. Figures 8-9 and S4-S5 show that evolutionary rescue has three main effects relative to the constant population case: 1) rescue can greatly reduce genome-wide diversity under sufficiently weak selection, where bottlenecks are long and deep (by increasing T and decreasing N e (τ ); Equation 17), 2) rescue tends to deepen dips in diversity when soft sweeps are possible, i.e., it hardens soft sweeps (by decreasing P off (2, T ) at the selected site, r = 0; Equation 18), and 3) rescue generally produces wider dips in diversity due to excess coalescence during the sweep (by decreasing P off (2, T ) at a given r > 0; Equation 18).

Tajima's D
Finally, we consider Tajima's D statistic (Tajima, 1989), which measures the relative excess (positive D) or deficiency (negative D) of intermediate frequency polymorphisms, relative to the standard neutral model (i.e., constant population size, neutral evolution). Quantitative predictions of Tajima's D require one to consider samples of size greater than 2, which quickly becomes complicated with selection and complex demography. Instead, here we discuss the expected qualitative patterns, based on intuition from the analysis presented above, and compare these to simulation results.
First, hard selective sweeps tend to produce star-like gene genealogies, with most samples coalescing near the beginning of the sweep and recombination allowing a few samples to coalesce much further back in time (Kaplan et al., 1989). Hard sweeps therefore produce an excess of low frequency polymorphisms (Wakeley, 2009, p. 120), leading to negative D (Braverman et al., 1995). The larger the selection coefficient the more star-like the genealogy (less time for coalescence or recombination during the sweep), and thus the more negative D when conditioned on a hard sweep. However, with sufficient standing genetic variation or rates of recurrent mutation or migration, larger selection coefficients will tend to cause softer selective sweeps  (Pennings and Hermisson, 2006b). As linkage to the selected site decreases so too does this skew in genealogies. In the case of a constant population size, D should asymptote to the neutral expectation of zero. In the case of rescue, however, the bottleneck will cause an excess of intermediate frequency polymorphisms (Wakeley, 2009, p. 120), and therefore D should asymptote at some positive value (more positive with more severe bottlenecks).
These patterns are borne out in simulations (Figures 10-11), where we see that rescue consistently causes positive background D. When sweeps are guaranteed to be hard (κ = 1), D around the selected site is similarly negative in both rescue and under a constant population size. When there is some possibility for a soft sweep (all cases but κ = 1), rescue tends to harden the sweep (by reducing probabilities of establishment, increasing coalescence, and reducing mutational input) and thus produce lower values of D at the selected site. Together these patterns cause rescue to stretch out or even invert the pattern of D observed in populations of constant size: under rescue, D tends to be greater away from the selected site and lower at the selected site.

Discussion
Here we have explored genetic signatures of evolutionary rescue by a selective sweep. By allowing demography to depend on the absolute fitness of the genotypes that comprise the population we explicitly invoke a feedback between demography and evolution. This feedback restricts the range of dynamics, and thus the signatures, that one should expect to observe. We find that, because the probability of establishment for an allele with a given selective advantage is reduced in declining populations (Equation 3; see also Otto and Whitlock, 1997), selective sweeps causing rescue are expected to be harder than those in populations of constant size when sweeps arise from standing genetic variance or recurrent mutation (Figures 2 and 4; consistent with Wilson et al., 2014 andWilson et al., 2017). Further from the selected locus, the demographic bottleneck experienced during rescue increases rates of coalescence relative to mutation and recombination (Figures 6-7), creating wider dips in diversity and lower diversity genome-wide (Figures 8-9; consistent with Innan and Kim, 2004). Tajima's D captures both the hardening of the sweep and the bottleneck, causing D to generally reach both higher and lower values under rescue (Figures 10-11). These differences between evolutionary rescue and standard sweeps all become larger under weaker selection (i.e., when the heterozygote has a smaller growth rate, s/2 − d 1) as the slower sweeps that result imply deeper, longer bottlenecks during rescue. In contrast to standing variance or mutation, when sweeps arise from a constant rate of migration demography has no affect on the number of beneficial alleles that establish (as briefly noted by Pennings and Hermisson, 2006a) and thus rescue has no affect on the hardness of the sweep ( Figure 5). Further, because the rates of coalescence and migration are both inversely proportional to the number of beneficial alleles N (t)p(t) at the selected site (Equation 12, Figure S7), the distribution of the number and frequency of migrant haplotypes spanning the selected site is given by Ewens' sampling formula (Ewens,  1972), with mutation replaced by migration (Pennings and Hermisson, 2006a). As we move away from the selected site recombination breaks apart these migrant haplotypes, leading to patterns of nucleotide diversity that depend on the migration rate and history of the two populations. If the migration rate is low we should expect an excess of migrant haplotypes at the site of the selective introgression and, if the divergence between the migrant and focal population is high, the so-called "volcano" pattern of diversity (Setter et al., 2019), where diversity is maximized at an intermediate distance from the selected site due to a more balanced presence of both migrant and non-migrant alleles. Evolutionary rescue has been explored theoretically (e.g., Gomulkiewicz and Holt, 1995;Uecker and Hermisson, 2016;Anciaux et al., 2018) and observed repeatedly in both experiments (e.g., Bell and Gonzalez, 2009;Lindsey et al., 2013;Ramsayer et al., 2013) and in host-pathogen systems in nature (e.,g., Wei et al., 1995;Feder et al., 2016). More recently, a number of studies have used genetic data to suggest that evolutionary rescue has occurred in the wild, including crickets becoming song-less to avoid parasitoid flies (Pascoal et al., 2018, reviewed in McDermott, 2019, killifish deleting receptors to tolerate pollution (Oziolor et al., 2019), hares moulting brown instead of white to avoid predation in snowless winters (Jones et al., 2018), bats altering hibernation to survive white-nose syndrome (Gignoux-Wolfsohn et al., 2018), and tall waterhemp evolving herbicide resistance (Kreiner et al., 2019). In nearly all of these cases there is strong evidence of a recent selective sweep by a very beneficial allele. Genetic evidence for a demographic bottleneck, on the other hand, is generally lacking, although genome-wide reductions in nucleotide diversity and increases in Tajima's D, relative to non-stressed populations, are sometimes detected (Oziolor et al., 2019). This begs the question of whether one can infer evolutionary rescue from genetic data alone, which would greatly help in assessing the relevance of rescue in nature. Strong support for rescue would come from the coincident timing of a sweep and bottleneck, which is difficult given the imprecise time estimates from a genetic sample collected from a single time point. Sampling before and after the potential rescue event would therefore be highly advantageous in determining co-occurrence. However, it should be noted that even if the sweep and bottleneck appear to have co-occurred, this correlation in timing does not imply it was caused by a feedback between demography and evolution. It is, of course, difficult to say in any case -without observing replicate populations go extinct or performing experiments -whether extinction would have occurred (or will occur) without adaptive evolution, as required by the strict definition of evolutionary rescue. We therefore need more experiments (such as Rêgo et al., 2019) that explore the genetic consequences of verified rescue to confirm the theoretical results presented here and help develop a robust signal of rescue to compare patterns from natural populations to. A strength of the above analysis is that we have explicitly modelled a feedback between demography and evolution, restricting the range of genetic signatures we consequently expect to observe. To take a recent example, Harris et al. (2018) have claimed that the lower reductions in genetic diversity within HIV populations adapting to less efficient drugs (as observed by Feder et al., 2016) could be due to weaker bottlenecks or slower sweeps rather than sweeps being softer, i.e. arising from multiple mutations. Fortunately, in this case genetic time-series data were available to show that the ability of HIV to reliably adapt on a short time-scale necessitates mutation rates and selection coefficients that imply adaptation by soft sweeps is likely (Feder et al., 2018). Formally modeling a feedback between demography and evolution also helps narrow the relevant parameter range. For example, under a haploid version of the model explored here (as is applicable to HIV) the minimum population size during rescue by new mutations is N (0)s(2N (0)s) −d/s /(s − d) (equation 22 in Orr and Unckless, 2014). Thus, for a given N (0) and s the minimum population size is e ln S/(2s), where S = 2N (0)s, implying that the minimum population size consistent with the model is roughly proportional to 1/s. The imposed feedback between demography and evolution therefore precludes simultaneously slow sweeps and large bottlenecks (for example negating the two smallest bottleneck sizes in figure 3A of Harris et al., 2018 and constraining one to the upper right portion of figure 3A -B in Feder et al., 2018). While it is very likely in this case that soft sweeps are indeed the cause of the pattern (Feder et al., 2018), incorporating an explicit model of how demography and evolution interact could help focus future debates.
The model presented here is but one model of evolutionary rescue, which involved a number of important assumptions. One of these is that the beneficial allele acts multiplicatively with the ancestral background, so that its marginal fitness is affected by the decline rate of the ancestral genotype. This in turn caused the dynamics of the sweep, once started, to depend only on relative fitness (Equation 1) while also making the probability of establishment (Equation 3) depend on the initial rate of population decline. If, instead, the absolute fitness of the heterozygote and mutant homozygote were independent of the initial decline rate, say 1 + hs and 1 + s , then the reverse would be true; the dynamics of the sweep would depend on the initial rate of population decline while the probability of establishment would not. We expect these effects would, however, largely cancel out. In any case, because our results depend primarily on the absolute and relative fitness of the heterozygote, the alternative model just described may closely match the model analyzed in detail here when hs is replaced by (hs + d)/(1 − d). A second key assumption we have made is that the beneficial allele acts additively with the ancestral allele at that locus (h = 1/2). Alternative forms of dominance will impact our results. At one extreme, a completely recessive beneficial allele (h = 0) is unlikely to establish, making rescue nearly impossible in outcrossing populations (Uecker, 2017). At the other extreme, complete dominance will greatly increase the probability of establishment and rescue (Uecker, 2017), as well as population mean fitness and thus population size. All else equal, we therefore expect rescue to have less effect on the signatures of selective sweeps relative to those in populations of constant size when the rescuing allele is more dominant. Given that the marginal fitness of the beneficial allele will not depend on allele frequency under complete dominance, the model will behave much more like a haploid model, where simple predictions of allele frequency and population size are more accurate (Orr and Unckless, 2014). Finally, it is of course possible to model rescue under much more complex lifecycles and population structure (e.g., as expected for the evolution of malarial drug resistance; Kim et al., 2014), at least using simulations. More complex lifecycles, such as those of parasites like Plasmodium and HIV, could cause the bottleneck to have additional impacts on the resulting genetic signature. For example, our populations are obligate sexual out-crossers, such that the probability of recombination does not depend on the population size (c.f., Equation 12), as everyone must mate with some one. However, with selfing and/or facultative sex (genetic exchange), rates of recombination could be lower at lower population densities, which would increase the impact of bottlenecks on resulting genetic signatures.
Evolutionary rescue is only one example of a myriad of processes where demography and evolution feedback on one another. This approach -combining forward-time eco-evolutionary models with coalescent theory to predict genetic signatures -could be used in many other scenarios. For instance, adaptive colonization of new habitat (a.k.a., adaptive niche expansion) is a closely related process for which a similar approach has already been taken (Kim and Gulisija, 2010). As in the case of rescue, explicitly modelling the feedback between demography and evolution in adaptive niche expansion changes the expected signatures left behind by selective sweeps as compared to Wright-Fisher populations. Such an approach is interesting from a conceptual point-of-view, improving our understanding of how eco-evolutionary dynamics affect genetic signatures. But further, given the computational power and simulation platforms available today, it is no longer necessary to restrict oneself to Wright-Fisher populations; researchers may now simulate under much more ecologically-realistic models. establishment and the rate of coalescence; see Supplementary text: Simulation details) as well as the same census population size (affecting initial allele frequency, p(0)) as in the case of rescue.

Supplementary text: Rescue by migrant alleles (MIG)
It is also possible for rescue to occur through the immigration of beneficial alleles. Assuming that the number of migrant alleles that replace a resident allele each generation is Poisson with mean m, the waiting time until the first successful migrant is exponential with rate λ = mρ. Given the population is expected to persist for log(N (0))/d generations, the probability of rescue is therefore roughly the probability the first successful migrant arrives by then, Dividing the waiting time distribution by the probability of rescue then gives the waiting time distribution conditioned on rescue, f (t). Using the approach we have taken in Rescue by de novo mutation (DNM), the effective initial frequency of the beneficial allele given rescue is (S4) When the migration rate is small this last factor is nearly independent of m (analogous to the mutation case). In a population of constant size the waiting time to the first successful migrant allele is simply exponential with rate λ = mρ, giving a waiting time factor 2m/(1 + 2m), which is strongly dependent on m but independent of ρ. Figure S6 compares our numerical (Equation 1) and analytical (Equations 2) approximations against individual-based simulations. We see the predictions do fairly well for larger values of m, but very poorly for small m. In the latter case the first successful migrant allele tends to arrive once the population is so small that beneficial homozygotes are regularly produced during establishment, causing us to greatly underestimate the probability of establishment (and thus overestimate allele frequencies and population size), as well as the rate of allele frequency increase. Figure S7 shows the timing of migration relative to recombination and coalescence (Equation 16). As with rescue from standing genetic variance or mutation (Figures 6-7), the bottleneck increases the overall coalescence rate and shifts its timing closer to fixation, overlapping more with recombination. Migration scales with coalescence (Equation 12) and is thus similarly increased and shifted.

Supplementary text: Deriving the structured coalescent
Let the allele frequency and population size τ generations before the present be p (τ ) and N (τ ). Following Pennings and Hermisson (2006a), we artificially subdivide the time within a generation, to be able to identify any period between two successive events ( Figure S1). We now go about deriving Equation 12.

Migration
The number of migrant alleles that arrive each generation is Poisson with mean m. Given that there are 2N (τ − 1)p (τ − 1) beneficial alleles in the next generation, the probability that any one is a new migrant where the approximation assumes the number of beneficial alleles changes little from one generation to the next. The probability that at least one of k beneficial alleles is a migrant is 1 − (1 − P ) k , which, with rare migration, is approximately For a given probability of being replaced by a migrant allele, the rate of migration in a diploid model is half that of the haploid model (equation 15 in Pennings and Hermisson, 2006a, replacing M with m) as there are twice as many resident alleles.

Mutation
The number of beneficial alleles after mutation, 2N (τ − 4/5)p (τ − 4/5), is the number before mutation plus the number of new mutants Because the population size does not change during mutation, N (τ − 4/5) = N (τ − 3/5), the frequency of beneficial alleles after mutation is simply The probability a beneficial allele is a new mutant is therefore u[1 − p (τ − 3/5)]/p (τ − 4/5), which, using the previous equation, is equivalent to The probability that at least one of k beneficial alleles is a new mutant is 1 − (1 − P ) k , which, when mutation is rare, is approximately ku[1 − p (τ − 4/5)]/p (τ − 4/5). With little change in allele frequency from one generation to the next this is This is equivalent to the haploid result (e.g., equation 5 in Pennings and Hermisson, 2006a) as both the mutation rate and number of alleles are multiplied by the ploidy level, which cancels.

Coalescence
Considering k beneficial alleles at the time of census, and ignoring any migration or mutation, the probability of at least one coalescence event is then

28
where N e (τ −2/5) is the effective population size at the time of syngamy. When allele frequency and effective population size changes little from one generation to the next this is roughly .
(S11) This is half the rate observed in a haploid model with the same population size (equation 5 in Pennings and Hermisson, 2006a) as there are twice as many alleles in a diploid population.

Recombination
Consider a neutral locus at recombination distance r from the selected site. Assuming weak selection such that the survivors of viability selection remain in Hardy-Weinberg proportions, the number of alleles linked to the beneficial allele after recombination is The first term on the right hand side is the number currently linked to the beneficial allele multiplied by the probability of being in a beneficial homozygote plus the probability of being in a heterozygote but not recombining. The second term on the right hand side is the number not currently linked with the beneficial allele times the probability of being in a heterozygote and recombining onto the beneficial background. The probability an allele on the beneficial background after recombination was not there before is then which, because recombination does not change allele frequency or population size, is P = 2N (τ − 1/5)[1 − p (τ − 1/5)]p (τ − 1/5)r 2N (τ − 1/5)p (τ − 1/5) = [1 − p (τ − 1/5)]r. (S14) The probability at least one of k alleles on the beneficial background recombines off is 1 − (1 − P ) k , which, when recombination is rare, is approximately kr[1 − p (τ − 1/5)]. Assuming allele frequency changes little through one bout of selection this is kr[1 − p (τ )]. Finally, assuming migration, mutation, and coalescence are rare, the probability that none of k beneficial alleles migrates or mutates times the probability none coalesce times the probability at least one of the k linked alleles recombines off is roughly (table 1 in Hudson and Kaplan, 1988) p rec (k, τ ) = kr[1 − p (τ )].