Following earlier models (e.g. Stewart and Levin, 1984; Berngruber et al., 2013; Sinha et al., 2017; Wahl et al., 2018), we use ordinary differential equations to describe a well-mixed system consisting of susceptible bacteria, lysogens (i.e. lysogenically infected bacteria), and free phages, but extend this system to include an arbitrium-like signalling peptide (Figure 1A). For simplicity, we consider phages that do not affect the growth of lysogenised host cells; susceptible bacteria and lysogens hence both grow logistically with the same growth rate r and carrying capacity K. Lysogens are spontaneously induced at a low rate α, after which they lyse and release a burst of B free phages per lysing cell. Free phage particles decay at a rate δ and adsorb to bacteria at a rate a. Adsorptions to lysogens result in the decay of the infecting phage, thus describing the well-known effect of superinfection immunity (Hutchison and Sinsheimer, 1971; Susskind et al., 1974; McAllister and Barrett, 1977; Kliem and Dreiseikelmann, 1989; Bondy-Denomy et al., 2016), whereas adsorptions to susceptible bacteria result in infections with success probability b. We consider the lytic cycle to be fast compared to both bacterial growth and the lysogenic cycle (Stewart and Levin, 1984; Berngruber et al., 2013; Sinha et al., 2017; Wahl et al., 2018), so that a lytic infection can be modelled as immediate lysis releasing a burst of B free phages. Since the genes encoding arbitrium production are among the first genes to be expressed when a phage infects a host cell (Erez et al., 2017; Stokar-Avihail et al., 2019), each infection leads to an immediate increase of the arbitrium concentration A by an increment c. The lysis-lysogeny decision is effected by the current arbitrium concentration: a fraction $\phi (A)$ of the infections results in the production of a lysogen, while the remaining fraction $(1-\phi (A))$ results in a lytic infection. Arbitrium does not decay spontaneously in the model (since it is a small peptide, spontaneous extracellular degradation is considered to be negligible), but it is taken up by bacteria at a rate u (e.g. through general bacterial peptide importers such as OPP [Erez et al., 2017]), and then degraded intracellularly, thus reducing the arbitrium concentration A.

Figure 1

Download asset Open asset

Consider competing phage variants i that differ in their (arbitrium-dependent) lysogeny propensity ${\phi }_i(A)$. The population densities (cells or phages per volume unit) of susceptible bacteria S, phage particles $P_i$ and corresponding lysogens $L_i$, and the concentration of arbitrium A can then be described by:

where $N=S+{\sum }_i{L}_i$ is the total density of bacteria.

We study two scenarios for the lysis-lysogeny decision: (i) a baseline scenario in which the arbitrium concentration does not affect the lysis-lysogeny decision; each phage variant has a constant lysogeny propensity ${\varphi }_i$ and (ii) a full scenario in which the arbitrium concentration does affect the lysis-lysogeny decision; each phage variant causes lytic infection when the arbitrium concentration is low, but switches to some lysogeny propensity ${\varphi }_{{\text{max}}_i}$ when the arbitrium concentration exceeds the phage’s response threshold ${\theta }_i$ (Figure 1B; note that we use ${\varphi }_{{\text{max}}_i}$ to denote a constant characteristic of the phage and ${\phi }_i$ to denote the function describing how phage variant i’s lysogeny propensity depends on the arbitrium concentration). In this second scenario, phage variants with constant lysogeny propensity are still included: variants with a response threshold ${\theta }_i=0$ cause lysogenic infections with lysogeny propensity ${\varphi }_{{\text{max}}_i}$ independent of the arbitrium concentration. Scenario (ii) is hence an extension of scenario (i).

On top of the ecological processes described in Equation 1–4, the model also allows for evolution of the phages due to mutations that change the characteristics ${\varphi }_{(\text{max})}$ and θ of phage variants. In Equation 1–4 terms describing these mutations were omitted for readability; they are described in detail in Appendix A1.1. In short, replication of any phage variant was assumed to produce mutants with slightly different characteristics (e.g. a slightly higher or lower lysogeny propensity) with a small probability μ. Under scenario (ii), mutations changing ${\varphi }_{\text{max}}$ and θ are implemented as independent processes.

In natural settings as well as in some laboratory experiments, phages regularly cause large outbreaks in pools of susceptible cells that were previously unavailable to the phages (e.g. when phages are spread to a new area, or when phages are serially passaged in a lab setting). Such outbreaks perturb the phage and cell populations away from their equilibrium. To mimic such repeated perturbations, we expose the system of Equation 1–4 to a phage serial-passaging regime (mimicking the experimental set-up of, for example, Bull et al., 1993; Bull et al., 2004; Bollback and Huelsenbeck, 2007; Betts et al., 2013; Broniewski et al., 2020). We initialise the model with a susceptible bacterial population at carrying capacity ($S=K$ cells per mL) and a small phage population (${\sum }_i{P}_i={10}^6$ phages per mL) and numerically integrate Equation 1–4 for a time of T hours. Then a fraction of the phage population is taken and transferred to a new population of susceptible bacteria at carrying capacity and Equation 1–4 are again integrated for T hours. This cycle is repeated to bring about a long series of epidemics. Throughout the manuscript, a dilution factor of $D=0.01$ is used (i.e. the passaged sample is 1% of the phage population). Passaging does not alter the relative frequency of the different phage variants, thus ensuring that the phage variants that were highly prevalent in the phage population at the end of an episode remain at a high relative frequency at the start of the new episode.

In this set-up, only phages are passaged from one epidemic episode to the next. To assess the robustness of simulation results to changes in this protocol, a second set-up was considered in which a fraction of the full sample (susceptible cells, lysogens, phages, and arbitrium) was passaged. We furthermore tested how the results are affected by variation in the bacterial carrying capacity. For this, at the start of each episode a carrying capacity value was sampled from a gamma distribution with mean K. We control the level of noise through the variance of this gamma distribution.

In total, the model (Equation 1–4) has nine parameters (excluding the phage characteristics ${\varphi }_i$, ${\varphi }_{{\text{max}}_i}$, and ${\theta }_i$, which vary between phage variants present in any given simulation). As far as we are aware, none of these have been estimated for phages of the SPBeta group, but many have been measured for other phages, most of which infect E. coli (Table 1, estimates taken from Little et al., 1999; De Paepe and Taddei, 2006; Wang, 2006; Shao and Wang, 2008; Zong et al., 2010; Berngruber et al., 2013). To reduce the number of parameters in our analysis, we nondimensionalised the equations to obtain five scaled parameter values (Appendix A1.3) and used the literature estimates to derive default values for these scaled parameters (Table 1). To account for the uncertainty in these estimates, we performed parameter sweeps consisting of 500 simulations with parameter values randomly sampled from broad parameter ranges (Table 1). To ensure that low values of the parameters were well-represented, parameter values were sampled log-uniformly.

Numerical integration was performed in Matlab R2017b, using the default built-in ODE-solver ode45. Scripts are available from https://github.com/hiljedoekes/PhageCom.

Next to numerical integration results, we analytically found expressions for the model equilibria and derived expressions for the evolutionarily stable strategy (ESS) under serial passaging in both scenarios (excluding and including arbitrium communication). Detailed derivations are provided in Appendix A3.

A common approach to analysing ODE-models such as Equation 1–4 is to characterise the model’s equilibrium states (Stewart and Levin, 1984; Wahl et al., 2018; Cortes et al., 2019). Such an analysis is provided in Appendix A2. However, we will here argue that to understand the evolution of arbitrium communication, and the lysis-lysogeny decision in general, considering the equilibrium states is insufficient.

Firstly, the function of the arbitrium system is to allow phages to respond to changes in the density of susceptible cells and phages as reflected in the arbitrium concentration. But when the system approaches an equilibrium state, the densities of susceptible cells and phages become constant, and so does the arbitrium concentration. Equilibrium conditions hence defeat the purpose of small-molecule communication such as the arbitrium system. Evolution of small-molecule communication must be driven by dynamical ecological processes, and hence can only be studied in populations that are regularly perturbed away from their ecological steady state.

Secondly, under equilibrium conditions natural selection can act on the lysis-lysogeny decision only if infections still take place, and hence lysis-lysogeny decisions are still taken. We argue that this is unlikely. If the phage population is viable (i.e. if the parameter values are such that the phages proliferate when introduced into a fully susceptible host population), the model converges to one of two qualitatively different equilibria, depending on parameter conditions (Appendix A2): either (i) susceptible host cells, lysogens and free phages all coexist, or (ii) all susceptible host cells have been infected so that only lysogens and free phages remain. The evolution of a constant lysogeny propensity in a host-phage population with a stable equilibrium of type (i) was recently addressed by Wahl et al., who show that under these conditions selection always favours phage variants with high lysogeny propensity (i.e. $\varphi =1$) (Wahl et al., 2018). However, only a narrow sliver of parameter conditions permits a stable equilibrium of type (i) (Sinha et al., 2017; Cortes et al., 2019), and when we estimated reasonable parameter conditions based on a variety of well-studied phages, we found that they typically lead to a stable equilibrium of type (ii) (Appendix A2, parameter estimates based on Little et al., 1999; De Paepe and Taddei, 2006; Wang, 2006; Shao and Wang, 2008; Zong et al., 2010; Berngruber et al., 2013). This is because phage infections tend to be highly effective: their large burst size and consequent high infectivity cause temperate phages to completely deplete susceptible host cell populations, replacing them with lysogens that are immune to superinfection and hence have a strong competitive advantage over the susceptible cells (Bossi et al., 2003; Gama et al., 2013). Under common parameter conditions, after a short epidemic the susceptible cells population is hence depleted and no more infections take place, causing the competition between different phage variants to cease (see Figure 2A for example dynamics). Then there is no long-term selection on the lysis-lysogeny decision, and studying its evolution in this state is pointless.

Figure 2 with 2 supplements see all

Download asset Open asset

We therefore consider a scenario in which the phage and cell populations are regularly perturbed away from equilibrium. To do so, we simulate serial-passaging experiments by periodically transferring a small fraction of the phages to a new population of susceptible host cells at carrying capacity, thus simulating cycles of repeated outbreaks (see Materials and methods).

To form a baseline expectation of the evolution of the lysis-lysogeny decision under the serial-passaging regime, we first considered a population of phage variants that do not engage in arbitrium communication, but do differ in their constant lysogeny propensity ${\varphi }_i$. Under typical parameter conditions (default values in Table 1), each passaging episode starts with an epidemic in which the susceptible cell population is depleted, followed by a period in which the bacterial population is made up of lysogens only (Figure 2A, dynamics shown for a passaging episode length $T=12$ h). The composition of the phage and lysogen populations initially changes over subsequent passaging episodes (Figure 2A), but eventually an evolutionarily steady state is reached in which one phage variant dominates the phage population ($\varphi =0.04$; Figure 2B), confirming that the lysis-lysogeny decision is indeed under selection.

The distribution of phage variants at evolutionarily steady state depends on the time between passages, T (Figure 2C). If this time is short ($T\le 5$ h), the phage variant with $\varphi =0$ dominates at evolutionarily steady state. This is an intuitive result: under these conditions phages are mostly exposed to environments with a high density of susceptible cells, in which a lytic strategy is favourable. Surprisingly, however, if the time between passages is sufficiently long ($T>5$ h), the viral population at evolutionarily steady state always centres around the same phage variant, independent of T ($\varphi =0.04$; Figure 2C). This result can be explained by considering the dynamics within a passaging episode (see Figure 2A): Once the susceptible cell population has collapsed, free phages no longer cause new infections and are hence ‘dead ends’. New phage particles are then formed by reactivation of lysogens only, so that the distribution of variants among the free phages comes to reflect the relative variant frequencies in the lysogen population. Hence, when the time between passages is sufficiently long, the phage type that is most frequent in the sample that is eventually passaged is the one that is most frequent in the population of lysogens (Figure 2—figure supplement 1). Under default parameter conditions, these are the phages with a low lysogeny propensity of $\varphi =0.04$.

Note that although a single phage variant clearly dominates the population, some diversity is maintained (Figure 2B–C). This is due to a mutation-selection balance: because mutants with slightly different $\varphi$-values continuously arise from the dominant phage variant and selection against these mutants is weak (due to their similarity to the dominant phage variant), the balance between influx of mutant variants by mutation and their efflux by selection results in the long-term presence of these mutants in the population. Such a population consisting of a dominant variant and its close mutants is called a quasi-species (Eigen, 1971).

Next, we assessed the robustness of the results to changes in the serial passaging protocol. In the standard protocol, only phages are passaged between episodes. If instead the passaged sample consists of the full system (susceptible cells, lysogens, and phages), almost identical results are obtained (Figure 2—figure supplement 2). This is again explained by realising that, as long as the time between passages is sufficiently long, the distribution of variants in the free phages is equal to the distribution in the lysogens. Since lysogens need to be induced to contribute to a new outbreak and the induction rate α is the same for all phage variants, the contribution of passaged lysogens to the new outbreak does not alter the relative frequency of phage variants.

To examine how these results depend on the model parameters, we determined which phage variant was most abundant at evolutionarily steady state for 500 randomly chosen parameter sets (see Table 1 for parameter ranges), always using a long time between passages ($T=24$ h). The selected $\varphi$-values for all parameter settings lie between $\varphi =0$ and $\varphi =0.12$ (y-axis of Figure 2D). We can hence conclude that selection favours phages with low but usually non-zero lysogeny propensities. These phages employ a bet-hedging strategy: throughout the epidemic they ‘invest’ a small part of their infection events in the production of lysogens, such that they are maximally represented in the eventual lysogen population.

To better understand how the lysogeny propensity $\varphi$ that is selected depends on parameter values, we derived an analytical approximation for the evolutionarily stable strategy (ESS) under the serial-passaging regime if the time between passages is sufficiently long (Appendix A3.1–2). Because the phage dynamics during an epidemic affect the dynamics of the susceptible cells and vice versa, phage fitness is frequency dependent and the ESS is not found by a simple optimisation procedure, but by identifying the particular $\varphi$-value, denoted ${\varphi }^{*}$, that maximises phage fitness given that this strategy ${\varphi }^{*}$ itself shapes the dynamics of the epidemic (Box 1). We find that the ESS can be approximated by the surprisingly simple expression

where P₀ is the density of phages at the start of a passaging episode. This approximation corresponds well with the results of the parameter sweep (Figure 2D), indicating that it indeed captures the most important factors shaping the evolution of the lysogeny propensity $\varphi$.

Box 1.

Lysogeny propensity of the evolutionarily stable strategy (ESS).

An evolutionarily stable strategy (ESS) is a strategy that cannot be invaded by any other strategy. In the context of the lysogeny propensity $\varphi$, it is the value ${\varphi }^{*}$ such that a population currently dominated by a phage with $\varphi ={\varphi }^{*}$ cannot be invaded by any phage variant with a different $\varphi$-value. A phage variant with $\varphi ={\varphi }_i$ invading in a resident population with the same $\varphi ={\varphi }_i$ always grows exactly like the resident. If this is the best possible invader, any other phage variant must perform worse than the resident and cannot invade. Hence, the ESS is the optimal response to itself. However, we still have to define what it means to be the ‘‘best possible invader’’ under the serial-passaging regime. Note that if the time between passages is sufficiently long, phages are selected on their ability to produce lysogens during the active epidemic (see Main Text). The optimal invader is hence the phage variant that, when introduced at a very low frequency, produces the most lysogens per capita between time $t=0$ and the time that the susceptible cell population collapses, $T_{\mathrm{E}}$. The $\varphi$-value of the optimal invader depends on $T_{\mathrm{E}}$ (red line in plot): if the epidemic phase is short, lysogens have to be produced quickly and a high $\varphi$-value is optimal, while if the epidemic lasts longer, phages can profit more from lytic replication and a lower $\varphi$-value is optimal. In turn, however, the duration of the epidemic $T_{\mathrm{E}}$ depends on the lysogeny propensity $\varphi$ of the resident phage population (blue line in plot): phages with a lower value of $\varphi$ replicate more rapidly and hence cause an earlier collapse of the susceptible population. The ESS is the value ${\varphi }^{*}$ that is optimal given the collapse time $T_{\mathrm{E}}({\varphi }^{*})$ that results when ${\varphi }^{*}$ itself is the resident strategy. Graphically, this value can be identified as the intersection of $T_{\mathrm{E}}(\varphi )$ and ${\varphi }_{\mathrm{opt}}(T_{\mathrm{E}})$ (the red and blue lines).

Box 1—figure 1

Download asset Open asset

Equation 5 shows that the ESS depends on the initial phage density in a passaging episode, P₀, relative to the burst size B and maximal host-cell density K, and the effective burst size $bB$, which represents the expected number of progeny phages per phage that adsorbs to a susceptible bacterium. The ESS ${\varphi }^{*}$ decreases with the dilution factor of the phages upon passage (i.e. with lower P₀). On the other hand, ${\varphi }^{*}$ increases with the effective burst size $bB$ (note that ${(bB)}^{-1}$ decreases when $(bB)$ increases). Both effects can be intuitively understood by considering how these factors affect the duration of the epidemic, $T_{\mathrm{E}}$. If the phage density is low at the start of a passaging episode or if the phages have a small effective burst size, it takes a while before the phage population has grown sufficiently to cause the susceptible population to collapse. Since a lytic strategy is favoured early in the epidemic, when the susceptible cell density is still high, a longer epidemic favours phages with lower values of $\varphi$ (see the red line in the figure in Box 1). On the other hand, if the initial phage density is high or if the phages have a high effective burst size, the susceptible cell population collapses quickly, phages have a much shorter window of opportunity for lysogen production and hence phages with higher $\varphi$-values are favoured.

Next, we included the possibility of arbitrium communication and let phage variants be characterised by two properties: their arbitrium response threshold, ${\theta }_i$, and their lysogeny propensity when the arbitrium concentration exceeds their response threshold, ${\varphi }_{{\text{max}}_i}$ (see Figure 1B). We then again considered the dynamics of our model under a serial-passaging regime.

In Figure 3A, example dynamics are shown for three competing phage variants, all with ${\varphi }_{\text{max}}=1$ but with different response thresholds ${\theta }_i$. The arbitrium concentration increases over the course of the epidemic, and the phage variants switch from lytic infection to lysogen production at different times because of their different response thresholds.

Figure 3 with 2 supplements see all

Download asset Open asset

Note that the maximum arbitrium concentration obtained during a passaging episode is approximately $A=cK$ (Figure 3A). This is because during the epidemic the dynamics of the susceptible cell density are mostly determined by infection events and not so much by the (slower) bacterial growth. Since the arbitrium concentration increases by an increment c every time a susceptible cell is infected, the infection of all initial susceptible cells will result in an arbitrium concentration of $A=cK$ (assuming that the degradation of arbitrium is also slow and can be ignored during the growth phase of the epidemic). The arbitrium concentration during the early epidemic then is a direct reflection of the fraction of susceptible cells that have so far been infected.

To study the evolution of arbitrium communication, we again considered the distribution of phage variants at evolutionary steady state for varying values of the time between passages, T. Similar to the results shown in Figure 2C, we find two main regimes (Figure 3B): if the time between passages is short ($T<4$ h, illustrated by $T=2$ h in Figure 3B), selection favours phage variants that only cause lytic infections (${\varphi }_{\text{max}}=0$); if the time between passages is sufficiently long ($T\ge 5$ h, illustrated by $T=12$ h in Figure 3B), the phage population is dominated by variants with ${\varphi }_{\text{max}}=1$ and $\theta \approx 0.65cK$. For $4\le T<5$ h, we see a transition between these two regimes (Figure 3—figure supplement 1). If the time between passages is sufficiently long ($T>5$ h), phage variants are hence selected that switch from a completely lytic to a completely lysogenic strategy when the arbitrium concentration exceeds a certain threshold.

In the simulations of Figure 3B, phage variants could have emerged that use the bet-hedging strategy found in the absence of communication (in phage variants with $\theta =0$, the lysogeny propensity is always ${\varphi }_{\max }$, independently of the arbitrium concentration), but this did not happen. We can hence conclude that any bet-hedging phage variants were outcompeted by variants that do use arbitrium communication. To underscore this conclusion, we simulated a competition experiment between the bet-hedging phage variant that was selected in the absence of communication and the communicating variant selected when arbitrium dynamics were included (Figure 3C). The communicating phage quickly invades on a population of bet-hedging phages and takes over, confirming that communication is indeed favoured over bet-hedging.

If a full sample (susceptible cells, lysogens, phages, and arbitrium) is passaged instead of phages only, again almost identical results are found (Figure 3—figure supplement 2). As was the case for the simulations in which arbitrium was absent, passaged lysogens do not alter the distribution of phage variants in the new outbreak. The passaged arbitrium does not significantly affect the outbreak dynamics either, because its concentration after dilution is much lower than the response threshold θ of the phage variants that are selected.

To study how the evolution of phage communication depends on phage and bacterial characteristics, 500 simulations were performed with randomly sampled sets of parameter values (Table 1), using a long time between serial passages ($T=24$ h). For each simulation, we determined which phage variant was most prevalent at evolutionary steady state. Although we varied the parameter values over several orders of magnitude, the most prevalent phage variant had a lysogeny propensity of ${\varphi }_{\text{max}}=1$ and a response threshold of $\theta =0.5cK$ or $\theta =0.6cK$ in almost all simulations (Figure 4A). Hence, over a broad range of parameter values, phages are selected that use the arbitrium system to switch from a fully lytic to a fully lysogenic strategy (i.e. ${\varphi }_{\text{max}}=1$). This suggests that over the course of an epidemic, there is an initial phase during which the lytic strategy is a ‘better’ choice (i.e. produces the most progeny on the long run), while later in the epidemic the production of lysogens is favoured and residual lytic infections that would results from a lysogeny propensity $\phi <1$ are selected against.

Figure 4

Download asset Open asset

To better understand the intriguing consistency in θ-values found in the parameter sweep, we used a similar approach as before to analytically derive an approximation for the response threshold ${\theta }^{*}$ of the evolutionarily stable strategy under the condition that the time between passages is long (Appendix A3.3). Again, we find a surprisingly simple expression for the ESS:

Note that the expression in Equation 6 again depends on the effective burst size $bB$, which is an indicator of the phage’s infectivity. The evolutionarily stable response threshold ${\theta }^{*}$ declines as the effective burst size increases, converging to a value of ${\theta }^{*}=\frac{1}{2}cK$ for highly infective phages (Figure 4B, green line). The same result was found for simulations of the competition between phage variants with different θ-values under different effective burst sizes (Figure 4B, blue dots). We see that Equation 6 provides a good prediction for the response threshold value that is selected over evolutionary time, especially for phages with high effective burst size (Figure 4B).

For phages with a very small effective burst size, the response threshold selected in the simulations tends to be lower than the analytical approximation. This is due to a violation of one of the simplifying assumptions made to arrive at the analytical approximation of Equation 6, namely that during the active epidemic the dynamics of the arbitrium concentration are dominated by its production through infections and arbitrium uptake and degradation by susceptible cells can be ignored. While this is a reasonable assumption in a fast progressing epidemic, it breaks down if the dynamics of the epidemic are slow, which is exactly the case if the effective burst size $bB$ is small. Under these conditions, the uptake and degradation of arbitrium by susceptible cells cause the arbitrium concentration to be lower than assumed in the analytical derivation. Consequently, the actual selected response thresholds (which are essentially arbitrium concentration values) are lower than the analytically predicted values.

The result in Equation 6 can be further understood biologically. Remember that the arbitrium concentration during the epidemic varies between $A=0$ and $A=cK$, and is a reflection of the fraction of susceptible cells that have so far been infected. It makes sense that the evolutionarily stable response threshold causes phages to switch infection strategy somewhere in the middle of the epidemic: if a phage variant switches to the lysogenic strategy too early, its free phage population does not expand enough to compete with phages that switch later, but if it switches too late, the susceptible-cell density has decreased to such a degree that the phage has missed the window of opportunity for lysogen production. The ESS results from a balance between the fast production of phage progeny during the initial lytic cycles and the eventual production of sufficient lysogens. For phages with a high effective burst size, this balance occurs around the time that half of the available susceptible cells have been infected. Phages with lower effective burst size are, however, predicted to switch later, because these phages need to invest a larger portion of the available susceptible cells in the production of free phages to produce a sufficient pool of phages that can later form lysogens. Note, however, that the range of biologically reasonable effective burst sizes includes many high values (range of x-axis in Figure 4B, Table 1), that is, many real-life phages have high infectivity. Hence, for natural phages in general, we predict that if they evolve an arbitrium-like communication system, communication will be used to switch from causing mostly lytic to mostly lysogenic infections when in an outbreak approximately half of the pool of susceptible bacteria has been infected.

So far, we have considered the evolution of arbitrium communication under highly predictable settings, with each outbreak taking place in a population of bacteria with the same initial density (i.e. the bacterial carrying capacity was constant). As argued above, in such a set-up the arbitrium concentration provides information on the density of susceptible cells still available for infection, and the phages use this to inform their lysis-lysogeny decision. While the bacterial carrying capacity can be kept constant in lab experiments, it is far from obvious that this would be the case in natural environments. This warrants the question of how robust the results are to variation in the bacterial carrying capacity.

We therefore performed simulations in which the carrying capacity varies from outbreak to outbreak. For a long time between passages (T = 24 h), at the start of each passaging episode a random carrying capacity was drawn from a gamma distribution with mean K and a pre-set variance that differs from simulation to simulation. We use the coefficient of variation (CV), which is defined as the standard deviation relative to the mean, to describe the level of noise.

Figure 5 summarises the results of these additional simulations. Surprisingly, a communication strategy with ${\varphi }_{\text{max}}=1$ and $\theta \approx 0.5cK$ is selected for a large range of carrying capacity noise up to CV ≤0.35 (illustrated by CV = 0.22 in Figure 5A; see Figure 5—figure supplement 1 for full data). In other words, even if the carrying capacity varies with a standard deviation up to one third of its mean value, the communication strategy described in the previous section is still selected.

Figure 5 with 1 supplement see all

Download asset Open asset

As the coefficient of variation increases even further, the arbitrium response threshold value θ of the selected phages decreases, and so does the lysis-lysogeny propensity that is used at high arbitrium concentration ${\varphi }_{\text{max}}$ (Figure 5B and C). These results make sense: if the carrying capacity strongly varies between passaging episodes, the phages regularly cause outbreaks in bacterial populations with low density. Phages with a response threshold value larger than the bacterial carrying capacity do not produce any lysogens during such an outbreak, which is disastrous for their long-term fitness. Hence, lower response thresholds are selected. The corresponding lower ${\varphi }_{\text{max}}$ values likely evolve to compensate for the earlier switch to lysogen production caused by the lower θ-values. In highly variable conditions, phages are hence selected to switch from a lytic strategy very early in the epidemic to a bet-hedging strategy later.

While we find much lower response threshold values when the variation in bacterial carrying capacity is high, these threshold values do remain clearly larger than zero (Figure 5C). This is true even if the carrying capacity is exponentially distributed (CV = 1; see Figure 5—figure supplement 1). Hence, even under very high variation of bacterial density a form of arbitrium communication (in which phages use the arbitrium signal to switch from a lytic to a bet-hedging strategy) is still favoured over completely bet-hedging strategies.

In this section, we find the dynamical equilibria existing in the model, and derive parameter conditions for their stability. This analysis provides us with a baseline expectation of the densities of phages, lysogens, and susceptible cells that the model converges to after sufficient time.

Equilibria of the model are found by equating Equation A.7–A.10 to zero and solving for the model variables. By the definition of equilibrium, the arbitrium concentration at equilibrium is constant:

where the bar indicates equilibrium values. Because the equilibrium arbitrium concentration is constant, at equilibrium the different phage variants are characterised by their lysogeny propensity ${\phi }_i(\overline{A})$ only, irrespective of whether arbitrium communication is included in the model or not. Hence, expressions for the model equilibria are the same in the absence and presence of arbitrium communication. Below, we derive these expressions, and determine under what conditions the different model equilibria are stable.

In the model, four qualitatively different types of equilibria can occur.

Firstly, there is a trivial equilibrium at $\overline{S}=0$, and ${\overline{L}}_i=0$, ${\overline{P}}_i=0$ for all phage variants i. This equilibrium is unstable as long as the bacterial logistic growth rate satisfies $r>0$.

Secondly, there is an equilibrium in which the susceptible population is at carrying capacity and the infection is absent: $\overline{S}=1$, and ${\overline{L}}_i=0$, ${\overline{P}}_i=0$ for all phage variants i. This equilibrium is stable if no phage-lysogen pairs $P_i$-$L_i$ can invade on the susceptible population. To derive stability conditions, we approximate the dynamics of $P_i$ and $L_i$ in the vicinity of the equilibrium by the linearised equations

where ${\overline{\varphi }}_i={\phi }_i(0)$ is the lysogeny propensity of phage type i at the equilibrium. No phage-lysogen pair can invade (i.e. the equilibrium is stable) precisely if the real parts of both eigenvalues of the Jacobian matrix J are negative for all i. The eigenvalues of the Jacobian matrix J are given by

The real parts of ${\lambda }_{+/-}$ are both negative precisely if $\mathrm{\Delta }<{\mathrm{\Gamma }}^2$ and $\mathrm{\Gamma }<0$. From $\mathrm{\Gamma }<0$, we find

while $\mathrm{\Delta }<{\mathrm{\Gamma }}^2$ yields

Since all parameters are non-negative and $0\le {\overline{\varphi }}_i\le 1$, the condition in Equation A.16 is more stringent than the condition in Equation A.15. Note also that the condition in Equation A.16 does not depend on the lysogeny propensity of the invading phages, $\overline{\varphi }$. Hence, the equilibrium $\overline{S}=1$, ${\overline{L}}_i=0$, and ${\overline{P}}_i=0$ for all phage variants i is stable exactly if the condition in Equation A.16 is satisfied. This condition makes sense: phages cannot spread in a susceptible cell population at carrying capacity if their infection rate $Ba\overline{S}$ is smaller than the decay rate of phage particles $\delta +a\overline{S}$.

Thirdly, there is a class of equilibria in which $\overline{S}=0$, and ${\overline{L}}_i>0$, ${\overline{P}}_i>0$ for some i. In these equilibria, the total densities of lysogens $\overline{L}$ and of free phages $\overline{P}$ are given by

In the absence of susceptible cells at equilibrium, no infections can take place and hence all phage variants behave identically (since phages vary only in ${\phi }_i(A)$, which occurs exclusively in the infection terms). This is reflected in the equations for the different lysogen variants, which for $\overline{S}=0$ and $\overline{L}=1-\alpha$ (Equation A.17) reduce to

Hence, any combination of ${\overline{L}}_i$ values with ${\sum }_i{\overline{L}}_i=\overline{L}=1-\alpha$ and corresponding ${\overline{P}}_i$-values,

is an equilibrium. Analogous to the reasoning above, such an equilibrium is stable if the susceptible cells cannot invade the phage-lysogen population at equilibrium. The linearised equation for the dynamics of S near the equilibrium of Equation A.17 reads

The right-hand side of this equation is negative (i.e. susceptible cells cannot invade) precisely if $\delta +a(1-\alpha )<Ba(1-\alpha )$, or summarised

(Note that to arrive at condition Equation A.20 we assume $\delta +a(1-\alpha )>0$. This assumption is justified because the spontaneous induction rate of lysogens α is small, and hence $1-\alpha >0$ (see Table 1).)

Lastly, there can be an equilibrium in which the susceptible cells, lysogens and phages all coexist. Expressions for $\overline{S}$, ${\overline{P}}_i$, and ${\overline{L}}_i$ at this equilibrium are bulky and not directly insightful. This type of equilibrium was however extensively analysed in recent work by Wahl et al. for phages with a constant lysogeny propensity (i.e. the restricted model where arbitrium communication is excluded) (Wahl et al., 2018). Remember that in an equilibrium state, phage variants are characterised by their lysogeny propensity ${\phi }_i(\overline{A})$ only (i.e. differences in response threshold θ are relevant only if they are reflected in differences in ${\phi }_i(\overline{A})$; phage variants i and j with ${\theta }_i\ne {\theta }_j$ but ${\phi }_i(\overline{A})={\phi }_j(\overline{A})$ can for all practical purposes be considered the same), and hence the results found by Wahl et al. can be extended to the model analysed here. Phage variants with different lysogeny propensity ${\phi }_i(\overline{A})$ can be seen as consumers that compete for a single resource, namely susceptible cells to infect. As long as $\overline{S}>0$, we hence expect competitive exclusion, and the phage population at equilibrium will be dominated by phages with a single lysogeny propensity value $\overline{\varphi }\equiv \phi (\overline{A})$ (when mutations are ignored, these phages will be the only ones present; otherwise we find a quasispecies). Furthermore, Wahl et al. show that if susceptible host cells coexist with a resident lysogen-phage population with some lysogeny propensity ${\overline{\varphi }}_{\text{r}}$, a phage-lysogen pair of a variant with higher lysogeny propensity ${\overline{\varphi }}_{\text{i}}>{\overline{\varphi }}_{\text{r}}$ can always invade on this equilibrium. Hence, in this equilibrium, the dominant phage will be the one with the highest equilibrium lysogeny propensity, $\overline{\varphi }=1$ .

As has been demonstrated previously (Stewart and Levin, 1984; Wahl et al., 2018), the ‘coexistence equilibrium’ is stable only if phages and lysogens can invade on a susceptible population at carrying capacity (i.e. the condition in Equation A.16 is violated), and susceptible cells can invade on the phage-lysogen population in equilibrium (i.e. the condition in Equation A.20 is violated). Hence, susceptible cells, lysogens and phages all coexist precisely if

If the effective burst size $B>1$ (a necessary condition for the phage to be viable), this can be rewritten as

Because the spontaneous induction rate of lysogens, α, is small ($\alpha <0.01$, see Table 1), the condition in Equation A.22 is very specific. Susceptible cells, lysogens and phages coexist only if the exponential growth rate of a lytic phage spreading in a susceptible population at carrying capacity, $(B-1)a-\delta$, is positive but very small, that is, if the epidemic is viable but only barely so. In reality, however, most phage epidemics are characterised by a high infectivity, mainly because of a large burst size (De Paepe and Taddei, 2006). Therefore, the condition in Equation A.22 is rarely satisfied, and for most phages we should instead expect to converge to equilibria of the third type ($\overline{S}=0,\overline{L}>0,\overline{P}>0$).

This observation has consequences for the selection pressures on phage variants over the course of a typical epidemic. As soon as the pool of susceptible host cells is depleted, competition between the different phage variants vanishes and the relative frequency of the variants freezes (see Figure 2—figure supplement 1 for an illustration). Under these conditions, no infections take place and hence there is no selection on the lysis-lysogeny decision. We conclude that the evolution of the lysis-lysogeny decision of typical phages requires regular perturbations away from equilibrium conditions.

In the main text, we present simulation results of the evolution of the lysis-lysogeny decision of phages exposed to a serial-passaging regime, both when arbitrium communication is excluded from the model (Figure 2), and when it is included (Figures 3 and 4). These simulations show that, after many passaging episodes, a single variant dominates the phage population (accompanied by its quasispecies, see Figure 2C, Figure 3B). We therefore assume that, within the parameter range of interest, a pure Evolutionarily Stable Strategy (ESS) exists. In this section, we derive analytical expressions for the ESS. We first provide a definition of the ESS under a serial-passaging regime, and give a general description of how this ESS can be found (section A3.1). Then, we apply this general approach to derive the ESS of phages that differ in a constant lysogeny propensity, $\varphi$ (absence of arbitrium communication, section A3.2), and the ESS of communicating phages that differ in their response threshold θ (section A3.3).

Consider a population of phages under a serial-passaging regime with long time T between passages. At the start of each passaging episode, a fraction D of the phages is taken from the end of the previous episode, and is added to a ‘fresh’ population of susceptible bacteria at carrying capacity. This procedure is repeated over many episodes. Within each episode, the dynamics of the susceptible bacteria, lysogens, phages, and arbitrium are described by Equation A.7–A.12.

An evolutionarily stable strategy (ESS) is defined as a strategy that cannot be invaded by any other strategy that is initially rare. To find the ESS, we therefore consider a scenario where an invader phage variant attempts to invade a resident phage variant. Below, we specify what it means for a phage variant to be able to invade in a resident phage population under the imposed serial-passaging regime.

Envision a resident population consisting of an isogenic phage population that has gone through many passaging episodes. Over time, the dynamics within these episodes have converged to a repeatable trajectory characterised by $P_{\mathrm{r}}(t)$, the resident phage density over time, $S(t)$, the density of susceptible bacteria over time, and $L_{\mathrm{r}}(t)$, the density of lysogens over time. At the start of one episode, now suddenly introduce a second phage with its own (possibly different) strategy. Crucially, the initial density of this invading phage $P_{\mathrm{i}}(0)$ is infinitesimally small. Consequently, during the first passaging episode the dynamics of the resident phage and the bacteria ($P_{\mathrm{r}}(t)$, $S(t)$, and $L_{\mathrm{r}}(t)$) are not affected by the invader. The invader is able to invade precisely if at the end of the first episode its frequency has increased relative to that of the resident, that is, if $P_{\mathrm{i}}(T)/P_{\mathrm{r}}(T)>P_{\mathrm{i}}(0)/P_{\mathrm{r}}(0)$.

Note that the relative frequency of an invader with exactly the same strategy as the resident does not change during an episode. Suppose that such an invader is the best-performing invader under the environment set up by the resident; then this implies that no invader can increase in frequency over a passaging episode, and therefore the resident strategy must be an ESS. In other words, in a resident phage population consisting of the ESS only, the ESS itself is the optimal strategy for an invading phage variant, that is, the ESS is the optimal response to itself.

What does it take to be the best-performing invader? To answer this question, we consider the dynamics within a single passaging episode in more detail.

If the time between passages T is long (see below for exact condition), and the parameter conditions are such that the system converges to an equilibrium with $S=0,P>0,L>0$ (typical parameter values, see Appendix 2), the dynamics within an episode can be divided in three distinct phases (see Figure 2—figure supplement 1):

1. Epidemic phase. A substantial population of susceptible bacteria ($S>0$) supports an ongoing epidemic. Free phages and lysogens are formed through infection of susceptible bacteria.

2. Transition phase. The population of susceptible cells has collapsed ($S\approx 0$). The lysogen population expands to fill up the space left behind by lysed cells. Free phages particles can no longer infect susceptible cells, and disappear through decay and adsorption to lysogens.

3. Equilibrium phase. The composition of the population is well-characterised by an equilibrium of ‘type 3’ (see Equation A.17). There is a small but consistent population of free phages that originates from lysogens through spontaneous induction. The distribution of phage variants in the free phage population is a direct representation of the relative frequency of the variants in the lysogen population (Equation A.18).

Let $T_{\mathrm{E}}$ be the time at which the susceptible population collapses, that is, the end of the epidemic phase. If the time between passages T is sufficiently larger than $T_{\mathrm{E}}$, the passage takes place during the equilibrium phase. The relative frequency of phage variants in the passaged sample then directly reflects the relative frequency of the corresponding lysogens. Since lysogens are only differentially formed through infection dynamics (and not through lysogen replication, which happens at the same rate for all lysogen variants), the relative frequency of the different lysogens is established during the epidemic phase and does not change afterwards. The performance of an invading phage can hence be measured by its lysogen production between $t=0$ and $t=T_{\mathrm{E}}$.

The dynamics of $S(t)$, and consequently $T_{\mathrm{E}}$, are determined by the resident phage: highly virulent resident phages (that cause many lytic infections, for instance because of a low $\varphi$-value) deplete the susceptible cell population faster than less virulent residents. The optimal invader under a certain resident is the phage variant that produces the most lysogens during the limited window of opportunity that it is offered by the environment set up by the resident. Since the ESS is the optimal response to itself, it is the strategy that, as an invader, produces the most lysogens during an epidemic phase set up by itself. We will use this reasoning to identify the ESS.

Consider the restricted model in which arbitrium communication is excluded and phages are characterised by a fixed lysogeny propensity $\varphi$. To find the ESS under this scenario, we take the following steps:

1. Derive how the duration of the epidemic, $T_{\mathrm{E}}$, depends on the $\varphi$-value of a resident phage population.

2. Find the optimal $\varphi$ given a fixed value of $T_{\mathrm{E}}$, that is, the $\varphi$-value that yields a maximal lysogen density at time $T_{\mathrm{E}}$.

3. Combine 1. and 2. to find the ESS: the $\varphi$-value ${\varphi }^{*}$ that maximises its lysogen density at time $T_{\mathrm{E}}({\varphi }^{*})$, the duration of the epidemic as dictated by its own $\varphi$-value, ${\varphi }^{*}$.

This approach is summarised in Box 1. Below, we provide the full analysis.

Next, consider a population of phages that do engage in arbitrium communication, again under a serial-passaging regime with sufficiently long time between the passages. Below, we use an approach similar to section A3.2, but more general, to derive the evolutionarily stable arbitrium response threshold, ${\theta }^{*}$. We take the following steps:

1. Describe the dynamics of an invading phage and its corresponding lysogens in an environment dictated by a resident phage.

2. Find the optimal invader response threshold under a fixed resident response threshold, that is, find the θ-value that maximises the invader’s lysogen production at time $T_{\mathrm{E}}$ when the dynamics of susceptible cells (and hence $T_{\mathrm{E}}$) are determined by a fixed resident phage.

3. Determine the ESS, ${\theta }^{*}$, as the optimal response to itself: the optimal invader response threshold (as found in step 2) if that same response threshold is the resident strategy.

We found that the results below are best understood in terms of the non-scaled model; in particular the (non-scaled) burst size of the phages turns out to be an important parameter. Therefore, the derivations below are presented for the dimensionalised equations Equation A.1–A.6.