key: cord-0480770-dq9oha6z
authors: Zhang, Xiyun; Ruan, Zhongyuan; Zheng, Muhua; Zhou, Jie; Boccaletti, Stefano; Barzel, Baruch
title: Epidemic spreading under pathogen evolution
date: 2021-02-19
journal: nan
DOI: nan
sha: 0ffc010e0ecb6df511206595544f9fc7b220a768
doc_id: 480770
cord_uid: dq9oha6z

Battling a widespread pandemic is an arms race between our mitigation efforts, e.g., social distancing or vaccination, and the pathogen's evolving persistence. This is being observed firsthand during the current COVID-19 crisis, as novel mutations challenge our global vaccination race. To address this, we introduce here a general framework to model epidemic spreading under pathogen evolution, finding that mutations can fundamentally alter the projection of the spread. Specifically, we detect a new pandemic phase - the mutated phase - in which, despite the fact that the pathogen is initially non-pandemic (R0<1), it may still spread due to the emergence of a critical mutation. The boundaries of this phase portray a balance between the epidemic and the evolutionary time-scales. If the mutation rate is too low, the pathogen prevalence decays prior to the appearance of a critical mutation. On the other hand, if mutations are too rapid, the pathogen evolution becomes volatile and, once again, it fails to spread. Between these two extremes, however, we observe a broad range of conditions in which an initially sub-pandemic pathogen will eventually gain prevalence. This is especially relevant during vaccination, which creates, as it progresses, increasing selection pressure towards vaccine-resistance. To overcome this, we show that vaccination campaigns must be accompanied by fierce mitigation efforts, to suppress the potential rise of a resistant mutant strain.

Another crucial aspect, absent in the present treatment based on pre-mutated strains, is the fact that pathogen evolution is responsive. As we tighten our mitigation, either through prophylactic measures [24] [25] [26] or via pharmaceutical interventions, we induce a selective pressure for mitigation resistance. For example, if one enforces social distancing to push the reproduction rate R 0 below the pandemic threshold, the pathogen becomes naturally pressured towards higher transmissibility. Similarly, if one employs therapeutic treatment to expedite recovery, natural selection will push the pathogen to higher drug-persistence.

To address this, we introduce here an evolving pathogen model, designed to capture the delicate interplay between the pathogen's spread and its developing fitness. The evolution, a random walk in fitness space, is driven by the pathogen's mutation rate. At the same time the natural selection, in which the fitter strains proliferate, is pushed by the epidemiological parameters, characterizing how fast a mutated strain propagates. Together, we identify a rather broad set of conditions -the mutated phase -in which a non-pandemic pathogen will eventually reach an evolved pandemic state.

We find that besides the classic epidemiological parameters, i.e. infection/recovery rates, two additional components factor in -the mutation rate governing the evolutionary time-scales, and the number of infected individuals, which determines the likelihood of a critical mutation to occur within the relevant time-frame. Therefore, the current prevalence ρ(t) of the pathogen has direct impact on its projected spread. This has significant implications pertaining to our two main mitigation strategies • Social distancing suppresses the reproduction number R 0 to below the pandemic threshold. However, if many have already been infected, i.e. ρ(t) is large, then a stricter suppression may be required to avoid the emergence of a critical mutation. This indicates that the projection of the spread, and hence its mitigation, depends on the current prevalence ρ(t) -a hysteresis phenomenon, unobserved in the classic modeling frameworks • Vaccination campaigns create strong evolutionary pressure towards a vaccine resistant mutation, whose risk, once again, is directly related to the current pathogen prevalence. Hence, to succeed, we show that vaccine roll-out, no matter how efficient, must be coupled with fierce suppression via social distancing.

Consider a social network of N individuals, linked through the adjacency matrix A ij , a random network with average degree k. At t = 0 the network experiences an outbreak, which then spreads via the susceptible-infected-susceptible 18 (SIS) dynamics. In the classic SIS formulation, the projected spread is driven by two time-independent parameters: the recovery rate µ and the infection rate β, whose ratio R 0 = kβ/µ, the reproduction number, determines the state of the system -pandemic (R 0 ≥ 1) or healthy (R 0 < 1). Here, however, the pathogen is allowed to evolve, therefore these parameter may change over the course of the spread. This is captured by the individual recovery rate

where the fitness F i (t) stands for the level of mutation of the pathogen carried by individual i at time t, hence the unmutated pathogen has F i (0) = 1. The above equation models the fact that (i) each individual i may carry a distinct version of the virus; (ii) this version my gradually change in time t due to mutations. The smaller is µ i (t), the higher is the transmissibility of the pathogen, as described by the evolving reproduction number

Indeed, a low rate of recovery µ i (t) extends the duration of the infectious state, providing individual i with more opportunities to infect their peers. Hence, as the r.h.s. of (2) indicates, pathogens with increased F i (t) exhibit higher reproduction, and therefore they spread more efficiently than their lower fitness competitors.

While mutation may also impact the transmissibility of the pathogen directly (by altering the value of β), in (1) we only evolve µ i (t), as indeed only the ratio R i (t) affects the spread in the SIS framework, and not the specific values of µ or β. Therefore, for simplicity, here we only track the pathogen evolution through µ i (t), and its subsequent R i (t), and set β stationary.

Complementary to this, in Supplementary Section XXX we examine the case of β-mutations.

The spread is driven by the following set of transitions: the process of infection between a pair of individuals i and j is modeled by

in which a susceptible (S) individual i interacts with their infected (I) neighbor j (A ij = 1) at rate β. This leads to both individuals becoming infected. The newly infected individual i inherits j's pathogen, and hence in (4) we set i's fitness at the time of infection equal to that of j. Both fitness parameters, F i (t) and F j (t) may later change via mutation. Next, we consider the process of recovery

in which an infected individual I i transitions to S i at the evolved recovery rate µ i (t) of (1). Finally, the process of mutation follows

capturing a random walk with variable step size δ i (t), i.e. a sequence of random shifts in fitness, caused by small discrepancies in the pathogen's reproduction. Note that F i (t) is prohibited from becoming negative, as, indeed, a below zero fitness in (2) is meaningless. The case where F i (t) does approach zero corresponds to µ i (t) → ∞ in (1), a limit in which recovery is instantaneous, and hence the pathogen is unfit for reproduction. Such strains will be rapidly eliminated from the pathogen pool.

The magnitude of each mutation step is extracted from a zero-mean normal distribution, namely δ i (t) ∼ N (0, σ 2 ). Consequently, in the limit where σ = 0, we have δ i (t) = 0 at all times, mutations are suppressed, and Eqs. (3) -(6) converge to the classic SIS model, with R i (t) = R 0 as the constant reproduction number. In contrast, as σ is increased, significant mutations become more frequent and the pathogens rapidly evolve. We therefore vary σ to control the mutation rate of the pathogens.

Taken together, our modeling framework accounts for the dynamics of infection and recovery (SIS) under the effect of pathogen mutation. As the spread progresses, pathogens evolve via Eq. (6), blindly altering their epidemiological parameters at random. Natural selection, however, will favor the positive mutations, in which δ i (t) > 0. Indeed, such mutations lead to higher fitness, reducing the recovery rate µ i (t), and consequently increasing R i (t). Such pathogens, with increased R i (t), will proliferate more rapidly, and will eventually dominate the population.

Critical mutation. Consider an outbreak of a pathogen with R 0 < 1, i.e. below the epidemic threshold. This can be either due to the pathogen's initial sub-pandemic parameters, or a result of mitigation, e.g., social distancing to reduce β. In the classic SIS formulation, such pathogen with fail to penetrate the network. However, in the presence of mutations (σ > 0) the pathogen may potentially undergo selection, reach R i (t) > 1, and from that time onward begin to proliferate. This represents a critical mutation, which, using (2), translates to

the critical fitness that, once crossed, may lead an initially non-pandemic pathogen to become pandemic. The smaller is R 0 the higher is F c , as, indeed, weakly transmissible pathogens require a long evolutionary path to reach pandemic spread.

Next, we analyze the spreading patterns of our evolving pathogens, seeking the conditions for the appearance of the critical mutation.

To examine the behavior of (1) -(6) we constructed an Erdős-Rényi (ER) network with N = 5, 000 nodes and k = 15, as a testing ground for several epidemic scenarios (Fig. 1) . Each scenario is characterized by a different selection of our model's three epidemiological parameters: µ and β, which determine the pathogen's unmutated reproduction R 0 , and σ, which controls the rate of mutation. We then follow the spread by measuring the prevalence ρ(t), which monitors the fraction of infected individuals vs. time. We also track the pathogen's evolution via the population averaged fitness F (t) = (1/N ) N i=1 F i (t). Pandemic (Fig. 1a, red) . In our first scenario we set µ = 0.1, β = 8 × 10 −3 and σ = 10 −2 . This captures a pandemic pathogen, which, using k = 15, has R 0 = 1.2 > 1, namely it can spread even without mutation. Indeed ρ(t) rapidly climbs to gain macroscopic coverage, congruent with the prediction of the classic SIS model, but this time constantly growing, due to the gradual, but continuous, increase in fitness F (t).

Mutated (Fig. 1c, green) . Next we reduce the infection rate to β = 1.67 × 10 −3 , an initial reproduction of R 0 = 0.25 < 1. This describes a pathogen whose transmissibility is significantly below the epidemic threshold, and therefore, following the initial outbreak we observe a decline in ρ(t), which by t ∼ 50 almost approaches zero, as the disease seems to be tapering off. In this scenario, however, we set a faster mutation rate σ = 1. As a result, despite the initial remission, at around t ∼ 15, the pathogen undergoes a critical mutation as F (t) crosses the critical F c = 1/R 0 = 4 (grey dashed line) and transitions into the pandemic regime. Consequently, ρ(t) changes course, the disease reemerges and the mutated pathogens successfully spreads.

Lethargic (Fig. 1b, blue) . We now remain in the sub-pandemic regime, with R 0 = 0.25, but with a much slower mutation rate, set again to σ = 10 −2 . As above, ρ(t) declines, however the pathogen evolution is now too slow -it is lethargic, and cannot reach critical fitness on time. Therefore, the disease fails to penetrate the network, lacking the opportunity for the critical mutation to occur.

Taken together, the dynamics of the spread are driven by three parameters: the initial epidemiological characteristics of the pathogen, µ and β, which determine R 0 , and the mutation rate σ, which governs the time-scale for the appearance of the critical mutation. Therefore, to determine the conditions for a mutation-driven contagion we investigate the balance between the decay in ρ(t) vs. the gradual increase in F (t).

To understand the dynamics of the evolving pathogen model, we show in Supplementary Section XXX that at the initial stages of the spread, the prevalence ρ(t) follows

The time-dependent exponential rate ξ(t) is determined by the epidemiological/mutation rates via

whose two terms characterize the pre-mutated vs. post-mutated spread of the pathogen. The first term, linear in t, captures the initial patterns of spread, which are determined by the original pathogen parameters, µ, R 0 . For R 0 < 1 this describes an exponential decay, a là the SIS dynamics in the sub-pandemic regime. At later times, however, as t 3 becomes large, the second term begins to dominate, and the exponential decay is replaced by a rapid proliferation, this time driven by the mutation rate σ. The transition between these two behaviors -decay vs. proliferation -occurs at τ c = 2(1 − R 0 )/3µσ 2 R 2 0 , which quantifies the anticipated time-scale for the appearance of the critical mutation F c in (7). This analysis portrays the mutated phase as a balance between two competing time-scales: on the one hand the exponential decay of the sub-pandemic pathogen, and on the other hand the evolutionary time-scale τ c for the appearance of the critical mutation. For the evolution to win this race the pathogen must not vanish before t = τ c . This imposes the condition (Fig. 2a- 

ensuring that at τ c there are still one or more individuals hosting the pathogen. Indeed, ρ(τ c ) < 1/N indicates that on average, at t = τ c less than a single individual is left in the infected pool. Under this condition, the critical mutation is too late, the spread has already tapered off, and the exponential spread driven by the positive term in (9) is averted.

Taking ρ(τ c ) from (8), we can now use (10) to express the boundary of the mutated phase, predicting the critical mutation rate as

where I 0 = N ρ(t = 0) is the number of individuals infected at t = 0. Equation (11) describes the minimal mutation rate required for the pathogen to evolve a pandemic strain. For R 0 = 1 it predicts σ c = 0, as such pathogen can indeed spread even without mutation. However, as R 0 is decreased, for example under mitigation, the pathogen prevalence rapidly declines, and hence it must evolve at an accelerated rate to reach critical fitness. This is expressed in (11) by an increased σ c , which approaches infinity as R 0 → 0.

To test our predicted phase transition we simulated an array of 1, 050 realizations of Eqs. (1) -(6), representing different epidemiological scenarios. We varied R 0 from 0 to 1.5, i.e. from nontransmissible to highly contagious, and scanned a spectrum of mutation rates from σ = 10 −3 to σ = 10, spanning four orders of magnitude. Simulating each scenario 50 times we observe the probability P for the disease to spread. This is done by tracking the pathogen's longterm prevalence ρ = ρ(t → ∞) and counting the realizations in which ρ → 0 vs. those where ρ > 0. As predicted, we find that the pandemic state, classically observed only at R 0 ≥ 1, now extends to lower R 0 in the presence of sufficiently rapid mutations. This gives rise to the mutated phase (green), in which an initially decaying contagion suddenly turns pandemic. The transition between the lethargic and the mutated states (grey zone) is well-approximated by our theoretical prediction of Eq. (11), depicted by the black solid line.

Equation (11) shows that σ c depends not only on the epidemiological characteristics of the pathogen (µ, R 0 ), but also on the initial condition, here captured by the number of infected individuals I 0 = ρ(t = 0)N . If I 0 is large the critical rate σ c becomes lower, in effect expanding the bounds of the mutated phase. To understand this consider the evolutionary paths followed by the pathogens as they reproduce. These paths represent random trajectories in fitness space, each starting from F i (0) = 1, and with a small probability crossing the critical fitness F c . The more such attempts are made, the higher the chances that at least one of these paths will be successful. Therefore, a higher initial prevalence I 0 of the pathogen increases the probability for the appearance of a critical mutation, enabling a mutated phase under a lower σ. Indeed, in Fig. 2d we find that the phase boundary shifts towards lower σ c as the initial prevalence is increased (grey shaded lines).

Hysteresis. This dependence on I 0 indicates that the transition of Eq. (11) behaves differently if we approach it from the pandemic state or from the healthy state. To observe this let us fix the mutation rate at σ = 0.1 and gradually increase R 0 , seeking the critical point where the system shifts to the mutated phase. This is mapped to a vertical trajectory in the σ, R 0 plane (Fig. 2d, yellow dashed line) . At each value of R 0 we instigate an outbreak with ρ(0) = 0.2, and observe its long-term prevalence ρ. For small R 0 this outbreak decays and the system reverts to the healthy state ρ = 0. However, as we transition into the mutated phase, here predicted at R 0 = R High = 0.6, the pathogen turns pandemic and its prevalence abruptly changes to ρ ≈ 0.85.

To reverse this transition the naïve approach is to push R 0 slightly below this critical point, for instance, by practicing social distancing to reduce transmission. The challenge is that now, moving in the opposite direction -from large to small R 0 -our initial condition is pandemic, with prevalence of order unity (∼ 85%), and hence I 0 ∼ N . Under these conditions, Eq. (11) predicts that, for a fixed σ, the critical R 0 is now lower, at R Low = 0.35. This results in a hysteresis phenomenon, in which criticality occurs at different points depending on the state from which we approach the transition (Fig. 2e) .

We find, therefore, that pathogen evolution fundamentally changes the phase space of epidemic spreading. First it predicts a broad range of conditions -the mutated phase -in which a subpandemic pathogen can gain prevalence. On top of that, it also predicts that this phase exhibits a discontinuous transition, characterized by hysteresis, a phenomenon unobserved in the classic SIS dynamics, yet congruent with with other models 29-35 that incorporate feedback between a pathogen's prevalence (ρ(t)) and its potency (R i (t)). These two observations have direct implications on mitigation:

• Soft mitigation is risky. Most mitigation strategies seek a minimal approach, aiming to drive R 0 just below unity. This is understandable as (i) major restrictions on social interactions are costly and difficult to sustain 36 for extended periods; (ii) having R 0 < 1, even by a small margin, is assumed to naturally suppress the spread, as it leads ρ(t) to decay exponentially towards zero. Our analysis, however, shows that this is insufficient. For R 0 1 we have σ c → 0, indicating that even a relatively stable pathogen, with a low mutation rate, may eventually break through. Our analysis accounts for this effect and through Eq. (11), provide us, given σ, the level of tolerable R 0 that sufficiently mitigates the mutated phase risk.

• The sooner the better. Another common assumption, driven by the classic epidemic phase-diagram, is that the projected state ρ(t → ∞) depends only on R 0 , i.e. the epidemiological parameters. The current state of the spread ρ(0) plays no role. The observed hysteresis, however, shows that successful mitigation strongly depends on the prevalence at the time of instigation. If the pathogen has already gained sufficient ground, we will need to suppress the reproduction number below R Low , namely the lower phase-boundary in Fig. 2e . It is, therefore, crucial to respond early, and instigate our mitigation when ρ(t) is still small, eradicating the pandemic before mutations risk its reemergence.

Bounded fitness. Our mutation process in Eq. (6) allows the pathogen an unbounded random walk in fitness space. In reality, however, there are practical restrictions on fitness, as R i (t) cannot grow ad infinitum. Therefore, we now consider our evolving pathogen model, substituting the mutation in (6) with

in which the pathogen fitness is bounded from above by F max and from below by zero. Setting F max = 20 we now revisit our phase-diagram (Fig. 3a) . For small σ, mutations are slow, and the evolution path is unaffected by the upper bound on F i (t). Therefore, we continue to observe the same transition as in the unbounded model of Fig. Fig1d . As we increase σ, however, we witness a second transition, this time back to the healthy state, indicating that now, mutations are too rapid. This captures the final phase of our evolving pathogen model -the volatile phase:

Volatile (Fig. 3, blue) . When the mutation rate is too high the pathogen fitness becomes unstable. On the one hand it can rapidly reach critical fitness, yet, on the other hand, due to the random nature of its frequent mutations, it fails to sustain this fitness -resulting in an irregular F (t), that fluctuates above and bellow the critical F c (Fig. 3c) .

To gain deeper insight into the volatile phase, consider the natural selection process, here driven by the reproduction benefit of the fitter strains. This process is not instantaneous, and requires several reproduction instances, i.e. generations, to gain a sufficient spreading advantage. With σ to high, natural selection is confounded, the pathogen shown no consistent gain in fitness and, as Fig. 3c indicates, ρ(t) decays exponentially to zero. In Supplementary Section XXX we use a time-scale analysis, similar to the one leading to Eq. (11) , to show that the volatile phase occurs when σ exceeds

This prediction is, indeed, confirmed by our simulated phase diagram (Fig. 3a, black solid line) .

Together, our phase-diagram illustrates the different forces governing the spread of pathogens in the presence of mutations. While classically, spread is prohibited for R 0 < 1, here we observe a new, previously undocumented pandemic phase, in which the disease can successfully permeate despite having an initially low reproduction rate. The conditions for this phase require a balance between three separate time-scales: (i) The time for the initial outbreak ρ(0) to reach near zero prevalence τ r ; (ii) The time for the pathogen to evolve beyond critical fitness τ c ; (iii) The time for the natural selection to lock-in the fitter mutations τ s . Pathogens with small R 0 , we find, can still spread provided that

The l.h.s. of (14) ensures that the pathogen can reach critical fitness before reaching zero prevalence. This gives rise to the first transition in Eq. (11) , between the lethargic and the mutated phases. The r.h.s. of (14) is responsible for the second transition, from mutated to volatile. It ensures that fitter pathogens do not undergo additional mutation before they have time to proliferate via natural selection. Therefore, we observe a Goldilocks zone, in which the mutation rate σ is just right: on the one hand, enabling unfit pathogens to cross the Rubicon towards pandemicity, but on the other hand, avoiding aimless capricious mutations.

Vaccination during an ongoing pandemic is, by nature, a competition between the rate of the vaccine roll-out and the spread of the pathogen. [37] [38] [39] Therefore, naïvely, to win this race all one has to do is disseminate the vaccine as efficiently as possible, aiming to reach the majority of the population before the pathogen does. This, however, ignores the role of mutations, which may gravely impact even the most efficient vaccination campaign. Such mutations may, generally, be less fit epidemiologically, i.e. have a lower F and consequently a lower R i (t). Therefore, absent a vaccine, they will be rapidly overcome by the faster spreading pathogen strains. However, once the vaccine becomes widespread, resistance, even if less contagious, becomes a highly desirable trait, and a resistant mutation, if occurs, will take over the population.

To examine this we consider a pathogen spreading via the susceptible-infected-recovered 18 (SIR) model, and mutating via (6) . However, this time we add two additional components:

• Vaccination. The population is vaccinated at a rate ν, quantifying the percentage of the (susceptible) population that receives the vaccine per unit time (day).

• Resistance. At each time-step, in addition to the mutations considered in Eq. (6), the pathogen may undergo a vaccine-resistant mutation with probability p. This mutation has no bearing on its epidemiological parameters µ, β, hence providing no additional spreading advantage, other than being resistant to the vaccine. We set p = 5×10 −5 , a rare mutation.

In Fig. 4 we track a two vaccination scenarios, capturing slow (ν = 2 × 10 −4 ) and rapid (ν = 2 × 10 −3 ) roll-outs. In the rapid scenario we set a quite optimistic ν, assuming a coverage of ∼ 0.2% of the population per day. Still, in both cases we find that vaccination cannot guarantee the eradication of the pathogen. For example, in our rapid scenario, a resistant pathogen emerges in 44% of realizations, capturing a significant risk for vaccination failure (Fig. 4d) . Hence, we find that rapid vaccination, in and of itself, may be insufficient.

The point is that as the vaccine prevalence rises, so does the selective pressure for resistance. Indeed, during the spread, there may be several instances of resistant mutations. But, in a predominantly unvaccinated population, such mutations rarely take over, since they are not necessarily fitter than all other coexisting strains. However, as the vaccination campaign progresses, and a significant portion of the population is immune to the non-resistant strain, these same mutations provide a fitness boost, and hence, when they occur, they will eventually gain prevalence. This is clearly observed in Fig. 4c ,f, as the serendipitous resistant mutation (black nodes) rapidly penetrates the network, rendering the vaccine moot.

To mitigate this risk we must strongly suppress the probability of the resistance mutation to occur. As we have no control over the microscopic biological processes driving the pathogen mutation, our only means to reduce this risk is by decreasing the prevalence of the pathogen ρ(t). Indeed, as our analysis has shown, critical mutations are more likely to appear if the pathogen is widespread. Therefore, we find that during a vaccination campaign it is crucial to continue the suppression of the pandemic. This, we emphasize is not just to overcome the spread -that can indeed be accomplished by vaccination alone -but mainly to mitigate the risk of vaccine-resistant mutations.

To observe this, in Fig. 4g -l we revisit our vaccination scenarios, but this time reinforced them with standard mitigation efforts, such as social distancing to reduce R 0 below unity. In case of slow vaccination we continue to observe a large portion of failed campaigns (Fig. 4g, 52% ). This is because the slow build-up of immunity afforded the pathogen sufficient time for evolving a critical mutation (red nodes). It then gained prevalence, increasing the chances for the resistance mutation (black nodes). Indeed, in Fig. 4i we observe first the appearance of red nodes, and only then black nodes, indicating that initially a critical mutation a là Eq. (7) boosts the spread, after which resistance is evolved.

Consequently, our only successful vaccination scenario combines both strategies -mitigation together with rapid roll-out. The latter ensures the pathogen has no sufficient time for spreading via critical mutation, while the former suppresses ρ(t), alleviating the risk of a resistance mutation. As a result, the few sporadic instances of resistance mutations fail to penetrate the network thanks to the mitigation (R 0 < 1), while the handful of critical mutations are contained thanks to the rapid vaccination (Fig. 4k,l, red and black) . Hence, with the pathogen prevalence declining, the risk of the resistant/critical mutation is mitigated, and the vaccination campaign emerges successful.

The phase diagram of epidemic spreading is a crucial tool for forecasting and mitigating pandemic risks. First, it identifies the relevant control parameters, such as µ, β and k in our SIS framework, or additional parameters in more complex contagion processes, whose value determines R 0 and hence the expected the patterns of spread. The phase boundaries, then, help us assess the state of the system -healthy or pandemic -and provide guidelines for our response. For example, social distancing to reduce k, therapeutic treatment to increase µ or mask wearing to suppress β, all aimed to navigate the system's location along the pandemic phase-diagram towards the desired healthy state.

The common thread binding all of these strategies is the assumption that the epidemiological control parameters themselves are constant in time, and hence our intervention must just push them beyond the static phase-boundary, form which point on the epidemic will decay towards ρ → 0 spontaneously. This is, indeed, relevant if the temporal evolution of µ, β is slow compared to the epidemic spreading dynamics -as observed in the case of our lethargic phase. However, once the epidemiological parameters can change at sufficient rate, it fundamentally changes the rules of the game. This is because now, not only are the parameters dynamic, but, thanks to natural selection, they also become responsive. If, for instance, we develop drug-based treatment to increase the recovery rate µ, we inevitably also generate selection pressure towards drug persistence. Similarly, if we vaccinate or practice distancing to reduce k, β, we initiate an evolutionary race towards higher transmissibility or vaccine resistance.

The result is a complex interplay between the spreading dynamics (R 0 ), the instantaneous prevalence of the pathogen (ρ(t)), and the dynamic evolution of its parameters (σ), which reshapes the pandemic phase diagram. It not only expands the pandemic risk to a range of R 0 < 1, but also predicts an explosive transition pattern, i.e. the hysteresis of Fig. 2e , that is not observed in standard epidemiological transitions. This altered phase diagram, and its abrupt first-order like transition, we have shown, has crucial implications pertaining to our mitigation strategies. Yet, more broadly, as a physical phenomenon, it offers an interesting mechanism for explosive transitions. Most often, such abrupt phase-shifts are caused by internal suppression rules, that hold back the transition until it breaks through in an explosive fashion. [40] [41] [42] [43] In contrast, here what holds back the transition is the waiting time for the critical mutation. Until its appearance the system behaves in one way (R 0 < 1), but once it occurs, the system suddenly enters the pandemic regime (R 0 > 1). The explosiveness is therefore traced to a local event, whose probability depends on the system's initial parameters (R 0 , I 0 , σ). This local event then changes fundamentally the state of the system -capturing a feedback between the system's phase and its intrinsic control parameters. We believe this describes a unique mechanism, inherent to the basic ingredients of our biological system, reproduction, mutation and selection. We now remain in the sub-pandemic regime R0 < 1, but increased the mutation rate to σ = 1. For small t we observe ρ(t) rapidly decaying (top). However, thanks to the rapid mutations F (t) reaches critical fitness (grey dashed line) within a short time. Following this point the disease reemerges and ρ(t) changes course, turning pandemic. This is observed in the network example through the appearance of sporadic instances of high fitness pathogens (middle, dark red nodes), which then spread to infect the majority of the population. (d) σ, R0 phase diagram. To systematically observe the different phases we varied R0 ∈ (0, 1.5) and σ ∈ (10 −3 , 10), capturing a total of 1, 050 epidemiological scenarios, with different initial µ, β and σ. For each scenario we ran 50 stochastic realizations and measured the probability P to have ρ(t → ∞) > 0, i.e. pandemic. We observe three phases with sharp boundaries between them. First, the pandemic phase (red) for R0 > 1, independent of σ, as predicted by the classic SIS model. In addition to that the subpandemic regime R0 < 1 is split into two phases: Under small σ, P tends to zero (blue) and the pathogen fails to spread, giving rise to the lethargic phase. For large σ, the spreading probability becomes almost certain, as p ∼ 1 (green), and we observe a mutation driven contagion. The gap between these phases (grey) indicates an abrupt transition from P → 0 to P → 1, a dramatic shift occurring within a narrow range of R0, σ values. This grey range is well-approximated by our theoretical prediction (solid black line) as appears in Eq. (11) . All simulations, here and throughout, were done on a random network of N = 5, 000 nodes and k = 15. The disease parameters were set to µ = 0.1 and the infection rate was set variably to β = µR0/k, to obtain the different values of R0.

The mutation rate σ is specified in each figure. In each scenario we set the initial condition to ρ(t) = 0.2. The critical mutation occurs at the minimum point (tc), which is below the unit line. Therefore the epidemic decays prior to the appearance of the critical mutation. Indeed, the stochastic simulation (blue solid line) approaches zero prevalence, never reaching the positive branch of ρ(t).

(b) Setting σ = 0.16 the system is at criticality. ρ(tc) is adjacent to the unit line, and hence we observe critical behavior: some realizations decay (blue), whereas others successfully mutate. (c) Under σ = 0.5, the system is in the mutated phase, ρ(tc) is sufficiently above the unit line and the critical mutation is reached with probability P → 1. (d) The lethargic-mutated phase boundary in Eq. (11) depends on the initial size of the infected population I0. Here we show this boundary for I0 = 10 2 , . . . , 10 8 . (e) The long term prevalence ρ = ρ(t → ∞) vs. R0 under σ = 0.1 (yellow path in panel (d)). Approaching from small R0 we begin with an initial small infection of I0 = 10 2 and observe a transition to the mutates phase at R0 = R High . In the opposite direction, however, as we begin with large R0 we approach the transition from an already pandemic state with I0 ∼ 10 4 . Now the phase boundary traverses through R0 = RLow. We, therefore arrive at a hysteresis phenomenon, in which the critical transition point depends on the current state of the spread. Consequently, preemptive mitigation, done when the spread is still at its embryonic stage (I0 small), is more effective than reactive mitigation, applied when I0 is already large. (12) . We now observe a volatile phase, in which ρ → 0 (blue), when σ is too large. Hence, the mutated phase (green) now only appears in the Goldilocks zone in which the mutation rate in not too high nor too low. The theoretical prediction of (13) is also shown (black solid line on right). (b) ρ vs. R0 under σ = 3 (yellow path in panel (a)). As opposed to the lethargic-mutated phase change, in the shift from volatile to mutated we observe a continuous second order transition. (c) ρ(t) vs. t in the volatile phase decaying, as predicted, to the healthy state ρ = 0. (d) F (t) changes rapidly thanks to the large σ, and crosses the critical Fc (grey dashed line) early on. However the rapid mutations prevent the slower processes of the natural selection from securing a steady increase in F (t). Hence, the achieved fitness cannot be stably sustained for the pathogen to continually spread. (e) Indeed, we observe multiple instances of critical fitness (dark red) that fail to reproduce and dominate the pathogen population.

Patient-derived mutations impact pathogenicity of SARS-CoV-2

SARS-CoV-2 ORF3b is a potent interferon antagonist whose activity is further increased by a naturally occurring elongation variant

The D614G mutation of SARS-CoV-2 spike protein enhances viral infectivity and decreases neutralization sensitivity to individual convalescent sera

Making sense of mutation: what D614G means for the COVID-19 pandemic remains unclear

Models of RNA virus evolution and their roles in vaccine design

Models cannot predict future outbreaks: A/H1N1 virus,the paradigm

Integrating influenza antigenic dynamics with molecular evolution

The evolution of epidemic influenza

Large-scale sequencing of human influenza reveals the dynamic nature of viral genome evolution

Ecology and evolution of the flu

The Evolution of Norovirus, the Gastric Flu

Emergence and clonal expansion of in vitro artemisinin-resistant Plasmodium falciparum kelch13 R561H mutant parasites in Rwanda

A Kelch13-defined endocytosis pathway mediates artemisinin resistance in malaria parasites

Epidemic spreading in scale-free networks

Epidemic spreading in community networks

Epidemic spreading with information-driven vaccination

Epidemic processes in complex networks

Risk factors for the evolutionary emergence of pathogens

Spreading of multiple epidemics with cross immunization

Pathogen Mutation Modeled by Competition Between Site and Bond Percolation

Modeling the evolution of COVID-19

The effects of evolutionary adaptations on spreading processes in complex networks

How will countrybased mitigation measures influence the course of the COVID-19 epidemic? The Lancet

COVID-19: How to win

Alternating quarantine for sustainable epidemic mitigation

Meet the trillions of viruses that make up your virome, The Conversation

Metagenomics and future perspectives in virus discovery

Explosive spreading on complex networks: The role of synergy

Effects of two channels on explosive information spreading

Macroscopic patterns of interacting contagions are indistinguishable from social reinforcement

Epidemic dynamics on an adaptive network

Connectivity disruption sparks explosive epidemic spreading

Epidemic spreading under infectionreduced-recovery

Modelling to contain pandemics

Dynamic patterns of information flow in complex net-works

Spatiotemporal signal propagation in complex networks

Digitizable therapeutics for decentralized mitigation of global pandemics

Explosive transitions in complex networks? structure and dynamics: Percolation and synchronization

Explosive first-order transition to synchrony in networked chaotic oscillators

Explosive synchronization in a general complex network

Explosive synchronization in adaptive and multilayer networks

Author contribution X.Z. developed the concept. X.Z., B.B. and S.B. designed the framework. X.Z. and Z.R. performed the numerical simulations. All authors jointly analyze the results and developed the analytic framework. X.Z., B.B. and S.B. wrote the paper.

All code to study and reproduce the results shown here will be made freely available online upon publication.

Slow vaccination = 2×10 -4 Fast vaccination = 2×10 -4No mitigation 0 = 1.2Mitigation 0 = 0.8 Still, due to the slow roll-out, the pathogen is afforded time to evolve higher fitness, reaching a critical mutation (red nodes), then gaining sufficient prevalence to also develop the rare vaccine resistance. Indeed, we find in panel (i) that first a critical mutation occurs, then only once it spreads, we observe instances of resistance mutations. (k)-(l) Therefore, the only strategy that ensures ∼ 100% success is to combine fast roll-out with mitigation. Now critical mutations, as observed through the sporadic instances of red nodes, fail to spread thanks to the widespread vaccination. Resistance mutations, on the other hand, are rapidly eliminated thanks to our mitigation driving R0 below unity. Together allowing us to eradicate the epidemic before a strain with both R0 > 1 (critical) and vaccine resistance ever evolves.