Fluctuating Selection in the Moran.

Contrary to classical population genetics theory, experiments demonstrate that fluctuating selection can protect a haploid polymorphism in the absence of frequency dependent effects on fitness. Using forward simulations with the Moran model, we confirm our analytical results showing that a fluctuating selection regime, with a mean selection coefficient of zero, promotes polymorphism. We find that increases in heterozygosity over neutral expectations are especially pronounced when fluctuations are rapid, mutation is weak, the population size is large, and the variance in selection is big. Lowering the frequency of fluctuations makes selection more directional, and so heterozygosity declines. We also show that fluctuating selection raises dn /ds ratios for polymorphism, not only by sweeping selected alleles into the population, but also by purging the neutral variants of selected alleles as they undergo repeated bottlenecks. Our analysis shows that randomly fluctuating selection increases the rate of evolution by increasing the probability of fixation. The impact is especially noticeable when the selection is strong and mutation is weak. Simulations show the increase in the rate of evolution declines as the rate of new mutations entering the population increases, an effect attributable to clonal interference. Intriguingly, fluctuating selection increases the dn /ds ratios for divergence more than for polymorphism, a pattern commonly seen in comparative genomics. Our model, which extends the classical neutral model of molecular evolution by incorporating random fluctuations in selection, accommodates a wide variety of observations, both neutral and selected, with economy.

ALL environments vary. Yet molecular sequence analyses interpret patterns of polymorphism and divergence, assuming that selection is constant and directional. We extend the classic neutral model of molecular evolution to incorporate fluctuating selection. Contrary to classical selection theory, fluctuating selection can promote polymorphism in haploids in the absence of frequency-dependent fitness effects. The conditions for neutral fixation are broadened so that alleles temporarily subject to directional selection might fix as if neutral. Fluctuating selection raises the d/d ratio for divergence more than for polymorphism, a pattern commonly seen in genomic comparisons.

Available field evidence suggests selection fluctuates in direction over time (Dobzhansky 1943; Fisher and Ford 1947; Lynch 1987; Cain ; Cook and Jones 1996; Saccheri ). However, the importance of fluctuating selection in patterning polymorphisms, probabilities of fixation, and evolutionary divergence remains poorly understood. Classical population genetics theory suggests that fluctuating selection promotes polymorphism whenever the geometric mean fitness of the heterozygote is greater than both homozygotes (Kimura 1954; Dempster 1955; Haldane and Jayakar 1963; Gillespie 1972, 1973; Jensen 1973; Karlin and Levikson 1974; Karlin and Liberman 1975; Felsenstein 1976; Maynard Smith 1998; Bell 2008). Fluctuating selection alone cannot promote polymorphisms in haploids simply because the clone with the largest geometric mean fitness inevitably wins the competition. Only if fitness is frequency dependent (Felsenstein 1976; Bell 2008, 2010), as in the lottery model (Chesson and Warner 1981), or if recombination between selected loci occurs (Kirzhner ), can selection promote polymorphism in haploids.

In contrast to classical theory, experiments demonstrate that fluctuating selection can maintain a haploid/clonal polymorphism in the absence of frequency-dependent effects on fitness (Yi and Dean 2013). Key is in recognizing two sources of variability when using the growth rates r and r over a period t to define the absolute fitnesses and One (well-recognized) is variability in the relative growth rate, The other (previously overlooked) is variability in the time available for growth, t.

Consider a serial transfer experiment. Following Yi and Dean (2013), let p and q = 1 − p be the frequencies of A and a immediately following dilution in fresh medium, and let e be the fold increase in population density after growth to carrying capacity in environment i. Then the growth of the mixed culture at carrying capacity is described by . This model [experimentally verified for strains of Escherichia coli competing for limiting glucose (Yi and Dean 2013)] describes exponential growth that ceases abruptly once the carrying capacity is reached. When A is rare the time taken to reach carrying capacity is . When a is rare the time taken is . Placing a bound on population size (i.e., restricting the fold-increase to ) means that the time spent growing to carrying capacity is a dependent variable, one that varies with allele frequency.

This frequency-dependent slippage in “time-for-growth” promotes diversity (Figure 1). Alleles experience more doublings when the least fit is most common. This favors the rare, fitter allele. Alleles experience fewer doublings when the fittest is most common. This mitigates selection against the rare, less fit allele. Fluctuating selection is passively biased in favor of rare alleles.

Figure 1

This serial transfer experiment with fluctuating selection shows two Escherichia coli strains (one resistant to chloramphenicol, Clm, and the other resistant to tetracycline, Tet) coexisting in a stable oscillation. Selection is not frequency dependent because the selection coefficients (the slopes) are constant and the lines are either parallel up or parallel down. On reaching carrying capacity (at the peaks and in the valleys), the culture is immediately diluted into fresh medium and the direction of selection switched. Coexistence is possible because the time spent growing in each antibiotic (horizontal bars) is frequency dependent. It takes longer to reach carrying capacity in tetracycline when Clm is common near the start of the experiment than in the stable oscillation at the end. Similarly, it takes less time to reach carrying capacity in chloramphenicol when Clm is common near the start of the experiment than in the stable oscillation at the end. Coexistence is not possible in the alternative transfer regime where the culture is diluted at a fixed time interval before carrying capacity is reached. Adapted from Yi and Dean (2013).

Motivated by the experimental results in Figure 1, we re-examined the impact of randomly fluctuating selection on polymorphism and evolution using two continuous time overlapping generation models. The first, of an infinite (though bounded) population, sets the stage for the second, of a finite population undergoing continuous mutation and allelic fixations. Results show that fluctuating selection indeed promotes polymorphism, increases rates of evolution, and raises the ratio of nonsynonymous-to-synonymous substitutions more for divergence than for polymorphism. Kimura’s neutral model of molecular evolution is the limiting case where selection is zeroed.

Methods

Infinite populations

Basic model:

The conditions needed for fluctuating selection to promote polymorphism in a population growing continuously at carrying capacity are identical to those in serial transfer (Dean 2005; Yi and Dean 2013). Consider two clonal (i.e., nonrecombining) populations (densities N and N) competing in a succession of n arbitrary environments (each of time interval t) for a continuously replenished growth-limiting nutrient (concentration R). Death rates are fixed and birth rates vary with resource availability. To make matters more concrete, and without loss of generality, let the populations inhabit a chemostat, a continuous culture device. The differential equations describing the competition are during interval i are:where the per capita growth rates r = f(R) and r = f(R) are increasing monotonic functions of R, the concentration of resource in the growth chamber. R.0 is the concentration of the limiting nutrient in the fresh medium entering the growth chamber. δi is the chemostat dilution rate (the fractional rate at which fresh medium enters the growth chamber and spent medium and cells are washed out of it). Y is the yield coefficient (the number of organisms produced per amount of limiting resource consumed).

During growth at quasi-steady state, dR/dt ≈ 0 and so the mean population growth rate equals the dilution rate, pr + qr = δ. The growth rates vary with allele frequency; when A is rare, and when a is rare, . These changes are driven by changes in the quasi-steady state concentration of the limiting nutrient, R. The total population density remains constant because R << R0.0 and so N + N ≈ YR.0 (Dean 2005).

Conditions for coexistence:

For a rare competitor to invade a resident population at equilibrium with its environment requires that . Hence, mutual invasion, and thus polymorphism, is ensured whenever the weighted arithmetic mean relative growth rates of A and a (respectively) are both greater than one:The weights, δ, are proportional to the number of population doublings, δLog2e, needed to maintain population density (equivalent to the DLog2e population doublings needed to reach carrying capacity in the serial transfer experiment). The weights vary with changes in the chemostat dilution rate and in the time spent in each environment (the fold-dilution and/or the carrying capacity in serial transfer experiments).

Taking reciprocals reveals that selective sweeps by common alleles are prevented when their weighted harmonic mean relative growth rates are <1. For A and a (respectively):Hence, polymorphism is ensured whenever an allele’s weighted arithmetic mean relative growth rate is >1 and its weighted harmonic mean relative growth rate is <1 (Dean 2005; Yi and Dean 2013).

Relative fitness in discrete and continuous time:

Population geneticists often define relative fitness in discrete time as and relative fitness in continuous time as . However, the suggestion that our model of selection is inherently frequency-dependent because r and r vary with allele frequency is wrong. Halving the growth rate of a necessarily doubles its generation time. A better way to write relative fitness is where r = Loge2 corresponds to one a generation (doubling), and t is the time needed to complete it. When quasi-steady-state growth rates are proportional to the concentration of the limiting resource, r = α and r = α, and both relative growth rate r/r = α/α and relative fitness (α/α − 1)Loge2 are unaffected by the changes in allele frequency. The conditions for a protected polymorphism simplify toTemporal fluctuations in fitness can protect a polymorphism in the absence of frequency dependent effects on relative fitness (Figure 1).

It is the number of allelic doublings (i.e., generations) per environment that varies in a frequency-dependent manner. To invade, A must grow faster than the resident a growing at equilibrium with the chemostat dilution rate, . As A approaches fixation, its growth rate must slow to the dilution rate while a is washed out of the growth chamber, . Hence, the number of doublings experienced by A is frequency dependent; A experiences doublings when rare but only doublings when common. Both alleles experience more doublings in environments where the least fit allele is most common.

Nonoverlapping generations:

Although we shall focus on continuous time overlapping generation models in this paper, it is worth pausing for a moment to consider the discrete time nonoverlapping generation version of our serial transfer model. The exponential growth of two competitors, can be rewritten as where W and W are the absolute fitnesses and g is the frequency-dependent number of generations to carrying capacity during period i. This reformulation admits two geometric progressions (say 1, 2, 4, 8, 16… and 1, 3, 9, 27, 81…). Adults need not survive from generation to generation. Overlapping generations are not needed for coexistence. This contrasts with the lottery model of Chesson and Warner (1981), where overlapping generations are essential to coexistence.

Frequency-dependent selection:

Classical negative frequency-dependent selection, with and , increases the likelihood of invasion by rare alleles and reduces the likelihood of fixation by common alleles. Positive frequency-dependent selection, with and , opposes the intrinsic bias in favor of rare alleles and can prevent a polymorphism arising. In all that follows, we assume relative growth rates are not frequency dependent (i.e., ), which is common during starvation competition for a single limiting resource.

Probability of a polymorphism:

Many environments vary cyclically and dramatically: tides, diurnal cycles, the seasons, etc. Imagine an infinite population with two alleles subject to fluctuating selection in an environment that repeatedly cycles through an arbitrary number of seasons equal in length. Simulations (Figure 2) show that the variability in growth parameters α and α needed to protect a polymorphism becomes progressively smaller, and the correlation between them becomes progressively less important, as the number of periods per cycle is increased. At 10 periods per cycle, coexistence is likely only when α and α are highly variable and negatively correlated (Figure 2, n = 10). Conversely, the more similar the responses to environmental change, the less likely the two alleles can coexist. Nevertheless, a protected polymorphism is still possible with correlation coefficients as high as 0.9, given α and α are sufficiently variable. By 104 seasons per cycle, a protected polymorphism is very likely, unless α and α vary little and are tightly correlated (Figure 2, n = 10,000).

Figure 2

The probability that a polymorphism is protected in an infinite population increases with the number of seasons per cycle from n = 10 to n = 10,000, and as the growth parameters become more variable (SD σ) and as trade-offs become more severe (correlation coefficient ρ < 0). Simulations were implemented in Mathematica. For each season i, pairs of growth rates (α and α) were drawn randomly from a bivariate normal distribution, with means μ = μ = 1, a common SD (0 ≤ σ ≤ 0.35 and the distribution truncated at zero), and a given correlation coefficient (−1 ≤ ρ ≤ 1). The probability of a protected polymorphism was determined as the proportion of 104 replicates in which both arithmetic mean relative growth rates were >1.

Random environments:

Environmental variables can vary erratically without cycling. Let the α and α for each period i be drawn randomly from a bivariate distribution, with means Eα ≠ Eα, SDs Var(α) ≠ Var(α), and a correlation coefficient −1 < ρ < 1. The expected relative growth rates of A and a when each is rare are, approximately (Kendall and Stuart 1979),assuming that the doublings in each period (DLog2e in a serial transfer experiment and δLog2e in a chemostat) are uncorrelated with the growth parameters α and α. The model suggests a rare less fit allele can still increase in frequency if the superior competitor’s growth rate is sufficiently variable. For example, A will increase in frequency when Eα/Eα < 1 if a’s variance is sufficiently large to make E(α/α) > 1. Displacing a is not possible, however, because Eα/Eα > 1; any variability in A’s growth rate further inflates E(α/α). The result is a balanced polymorphism in which the high fitness of a is offset by its increased sensitivity to environmental variability. Coexistence is not possible when the difference between the scaled variance and covariance [Var(α)/E2α − Cov(α, α)/EαEα] is less than the difference in the relative growth rates (Eα/Eα − 1). Polymorphism is possible only when a considerable overlap in the fitness distributions exists.

If the expected growth parameters are equal (Eα = Eα = 1) and identically [Var(α) = Var(α)] though perhaps not independently [Cov(α, α) ≠ 0] distributed, thenThe slightest variability in clonal growth rates () guarantees a protected polymorphism except when they are perfectly positively correlated (ρ = 1). This result stands in marked contrast to all analyses from classical population genetics, which predict that clonal polymorphisms cannot be protected by fluctuating selection (Kimura 1954; Dempster 1955; Haldane and Jayakar 1963; Gillespie 1972, 1973; Jensen 1973; Karlin and Levikson 1974; Karlin and Liberman 1975; Felsenstein 1976; Chesson and Warner 1981; Maynard Smith 1998; Bell 2008, 2010).

Finite populations without mutation

Polymorphisms cannot be protected in finite populations because rare alleles risk stochastic loss. Instead, fluctuating selection promotes polymorphism by retarding, rather than preventing, the loss of alleles. Weak selection fails to retard the loss of alleles to drift. Strong selection combined with infrequent environmental fluctuations drives the directional loss of genetic diversity. The capacity for fluctuating selection to promote polymorphism therefore depends on the interaction between drift, selection, and the frequency of environmental changes.

The Moran model:

We explored the interaction between drift, selection, and the frequency of environmental changes with a Moran model (Moran 1962) that we embedded in a resource, thereby forming an obvious stochastic analog to competition in chemostats. In both, reproduction and death are continuous and the population finite and fixed in size (N). The probability of a birth is proportional to the resource concentration. The probability of a death is determined by a fixed density independent death rate, d, which is equivalent to the chemostat dilution rate, δ.

In the Moran model, each birth is coupled to a death and so the population size remains constant. The transition state probabilities for A, frequency j/N, arewhere the growth rates r and r fluctuate over periods i = 1, 2… The density independent death rate does not appear in Equation 9a because it always cancels [e.g., d/(d) = j/N].

Just as in chemostats, the absolute growth rates vary with allele frequency. When A is rare, (as before, is the resource abundance needed to maintain a pure population of a at its equilibrium density), and when a is rare, ( is the resource abundance needed to maintain pure population of A at its equilibrium density). Again, proportional changes in growth rates render the relative growth rates, r/r = α/α, independent of allele frequency. Equation 9a can be rewritten aswhich shows that the probability that the next event, either the birth or the death of an A allele, is unaffected by changes in resource abundance. In both chemostats and the Moran model, the numbers of births per hr are determined by the death rates (with the allelic growth rates forced to vary proportionally to the allele frequencies). The time taken to produce N births is t = Log2/δ hrs in chemostats and t = Loge2/d in the Moran model. In consequence, an environment that changes every xt hrs also changes every xN birth–death events.

The distribution of A alleles (v = 0 to N) following each birth–death event is given by the vector v + 1, = Mv (i = 1, 2 …) where M is the square (N + 1) tridiagonal transition matrix for period i [with probabilities of: no change in number (k = j) on the diagonal, unitary increases above (k = j + 1), and unitary decreases below (k = j − 1)] and v is the vector of frequencies at time t (measured in birth-death events). In a cyclical environment with two alternating periods the mean number of generations, , until loss of one or other allele is (i = 1, 2), where t = N is one generation.

Simulations were implemented in Mathematica. The mean time to fixation (the loss of either allele), for a haploid population of N = 50 individuals, initially with j = 25 A alleles and assuming evolution by random genetic drift alone (α = α), is = 34,164 generations. This estimate is close to the expected time to fixation of = N(1 − p)/p Loge(1 − p) = 50 Loge2 = 34,657 generations. We attribute the small discrepancy to the theory assuming a much larger population size than 50 (Moran 1962).

Symmetrically fluctuating fecundities

Median time to fixation:

What effect does fluctuating selection have on the rate of loss of polymorphism in a finite population? We computed the median time to fixation starting with equal allele frequencies (j/N = 0.5) and with symmetric relative fitnesses (α0.1/α0.1 = α0.2/α0.2 ≠ 1), which would guarantee a protected polymorphism if the population were infinite. Rapid switches in the direction of selection (every generation, τs = 1) retard the loss of polymorphism (Figure 3).

Figure 3

The capacity of fluctuating selection to promote polymorphism in a finite population depends not only on population size, but also on the frequency of environmental change. In the symmetric fitness model (α0.1/α0.1 = α0.2/α0.2 ≠ 1), fluctuations in every generation (τ = 1) shepherd alleles close to the nodal frequency of 0.5, thereby retarding the loss to drift. As fluctuations become less frequent (every 10 or 20 generations) the amplitude of the oscillations in allele frequency increases and, consequently, the risk to loss by drift increases. The extreme limit is directional selection where polymorphism is no longer promoted. Simulations were implemented in Mathematica.

Switches every τ = 10 generations produce more nuanced results. Weak selection does little to promote polymorphism because it fails to shepherd allele frequencies close to the nodal frequency of 0.5. Strong selection accelerates the stochastic loss of alleles by driving their frequencies close to the boundaries of 0 and 1 at the end of each period. With τ = 10, the loss of diversity is slowest in large populations with intermediate selection intensities (s ≈ 0.65 in Figure 3). Switches every τ = 20 generations further reduce both the time to fixation and the range in fitness that promotes polymorphism.

How does variability in the frequency of switching affect the rate of loss of polymorphism? Variability in the frequency of switching has little impact on the rate of loss of polymorphism when selection is moderate in populations with N > 100 (compare switching every one and 10 generations for s < 0.1 in Figure 3). A single period >20 generations is sufficient to purge polymorphism in populations with N < 103 when selection is strong (s > 0.5). Rare extended periods of environmental stasis purge variation, especially when selection is strong.

Probabilities of fixation:

What effect does fluctuating selection have on an allele’s probability of ultimate fixation? When the selection coefficient s = |α/α − 1| < 1/N, alleles behave as selectively neutral (Moran 1962) regardless of whether selection fluctuates. With strong symmetric fluctuating selection (α0.1/α0.1 = α0.2/α0.2 >> 1) the rare currently beneficial allele will be driven back down in frequency when the selection reverses. However, it is not expected to return to its starting frequency because more doublings occur on the way up, when the least fit allele is most common, and fewer doublings occur on the way down, when the fittest allele is most common. This bias ensures that rare alleles are more likely to fix than if selectively neutral (Figure 4).

Figure 4

Probability of ultimate fixation (Up0) of an allele, initial frequency p0, subject to symmetric fluctuating selection (α0.1/α0.1 = α0.2/α0.2 = 1 + s) in a population (N = 103) with selection fluctuating each generation (τ = 1). Rare alleles are more likely to fix than neutral alleles (for p0 < 0.5, as the surface steepens relative to the neutral edge at s = 0 as s increases). The probability that a rare allele reaches fixation is far lower when selection fluctuates than when selection is directional. Simulations were implemented in Mathematica.

Before considering what happens with fluctuating selection, let us first remind ourselves of fixation by a new allele subject to directional selection. The probability of ultimate fixation of a new allele, frequency 1/N with constant relative fitness 1 + s, is Up0 ≈ s for Ns > 1 (Moran 1962). Most beneficial alleles are lost when rare. Any beneficial allele that reaches an appreciable frequency will march, almost deterministically, toward fixation. For example, U0.05 = 0.99 for s = 0.01 and N = 10,000. In effect, Up0 ≈ s is the probability that a beneficial allele, though still rare, attains sufficient numbers (∼p0Ns > 4.6 for Up0 = 0.99) to have escaped loss to drift.

Now consider a new allele (p0 = 1/N) with symmetric fitnesses (α0.1/α0.1 = α0.2/α0.2 = 1 + s ≠ 1) in a finite population (size N) with rapid switches in the direction of selection (τs < 1). Under these circumstances the mean selection coefficient of a rare allele is ≈ s2/2. Hence, the probability that the allele attains a frequency sufficient to escape initial loss to drift is ∼s2/2. The allele now enters a “zone of attraction,” where the impact of drift is sufficiently minimal that selection trains the allele to a stable oscillation around the nodal frequency of 0.5. This is reflected in Figure 4 as the inflection in the Up0 surface at p0 = 0.5 when selection is very strong (s = 0.25). Escape from the oscillation is random and does not depend on the initial allele frequency, hence the inflection. The probability of fixation from the stable oscillation is 1/2 and so the probability of ultimate fixation of a new allele with symmetric fecundities [1 + s and (1 + s)−1] is Up0 ≈ 1/2 s2/2. When s2/4 < 1/N, the allele will fix as if selectively neutral.

Simulations (Figure 5A) confirm that for sufficiently strong selection, the probability of ultimate fixation is proportional to the square of the selection coefficient, Up0 ≈ s2/4. Alleles fix as if selectively neutral (or nearly so) for Ns2 < 4 (found graphically as the intersection of lines Up0 = 1/N and Up0 = s2/4 on the log-log plots in Figure 5A). These results contrast with the classic analytical result for fluctuating selection in the Fisher–Wright model. There, τs < 1, rather than Ns2 < 4, dictates that alleles fix as if selectively neutral (Takahata ).

Figure 5

Forward simulations of the impact of random genetic drift, fluctuating selection, and the frequency of environmental change on the probability of ultimate fixation of a new allele, p0 = 1/N, (A) with symmetric fitnesses (α0.1/α0.1 = α0.2/α0.2 = 1 + s ≠ 1) and (B) with random growth parameters (α.i and α.i) drawn from a normal distribution N(1, σ). Infrequent environmental changes (solid lines, τ = 105) impose directional selection with the probability of fixation for a selected allele effectively halved. Rapidly fluctuating selection (dashed lines, τ = 1) expands the region where alleles fix as if selectively neutral. Simulations were implemented in C.

How does the probability of ultimate fixation depend on the frequency of environmental switching? Infrequent switching imposes directional selection for many generations. A new allele either fixes with probability s if it is sufficiently beneficial (probability 0.5), or is purged with probability 1 if it is deleterious (probability 0.5). The chance that a new mutant allele ultimately sweeps through a population inhabiting an infrequently switching environment is therefore Up0 = s/2 and the boundary for effective neutrality is Ns < 2 (confirmed in simulations as the intersection of lines Up0 = 1/N and Up0 = s/2 on the log-log plots in Figure 5A).

We explored the impact of τ on the probability of ultimate fixation in a population sized N = 103. Probabilities of ultimate fixation rise steeply between τ = 10 and τ = 100 generations for s > 0.04 (Figure 6A) and more gradually between τ = 100 and τ = 103 generations for s ≤ 0.02. Weakly selected alleles (s < 0.04) fix as if selectively neutral in rapidly changing environments but can become exposed to selection as τ is increased. Strongly selected alleles (s > 0.02 in these simulations) are more likely to fix in rapidly changing environments (τs < 1) than in the Fisher–Wright model (Takahata ). Simulations show Up0 ≈ 0.01 > 1/N = 1/103 with s = 0.2 and τ = 0.1 even though τs = 0.02 << 1 (Figure 6A). Variations in switching rates above τ = 103 and below τ = 10 have little effect on fixation probabilities in these simulations.

Figure 6

Forward simulations of the impact of random genetic drift, fluctuating selection, and period length (τ) on the probability of ultimate fixation of a new allele, p0 = 1/N in a population sized N = 103. (A) Simulations with symmetric fitnesses (α0.1/α0.1 = α0.2/α0.2 = 1 + s ≠ 1). (B) Simulations used random fitnesses drawn from a normal distribution, N(1, σ). Simulations were implemented in C.

Randomly fluctuating fecundities

We explored the balance between drift and selection for a biallelic Moran model using computer simulations with the following loop encoded in C:

Randomly choose an allele destined to die (based on allele frequency alone).

Randomly choose an allele destined to reproduce (based on allele frequency weighted by growth rate).

Increase the “birth” allele frequency by one.

Reduce the “death” allele by one.

Return to step 1.

At each environmental change, new allelic growth rates were drawn randomly from a normal distribution, N(1, σ), truncated at zero and with 0 ≤ σ ≤ 0.25. In this model, there is no recombination and the current growth rates are independent of the previous growth rates.

Probability of fixation:

The fitness SD in our model plays a role similar to the selection coefficient in the symmetric fitness model. With rapid changes (τs ≤ 1), the mean fitness of a rare allele is (Equation 8, a and b), making the probability of ultimate fixation of a new allele . New alleles fix as if neutral when . Forward simulations with ρ = 0 confirm that lines Up0 = 1/N and Up0 = σ2/2 intersect at on the log-log plot in Figure 5B.

In an infrequently switching environment selection is effectively directional and so the probability that a new allele reaches fixation is with weak selection. The selection coefficient, s = α − α, is distributed as (Kendall and Stuart 1979). A new allele is purged with probability 1 if s < 0 (probability 0.5). The probability of ultimate fixation when s > 0 is therefore approximately half the mean of the folded distribution (Kendall and Stuart 1979), . Simulations with ρ = 0 confirm that lines Up0 = 1/N and intersect at approximately on the log-log plot in Figure 5B. Variability in the frequency of environmental shifts outside the range 0.1 < τs < 102 has little impact on the probabilities of fixation (corresponding to the range 1 < τ < 103 in Figure 6B).

Finite populations with mutation

We explored the balance between mutation, drift, and selection for a multiallelic Moran model using by modifying the computer simulation loop:

Randomly choose an allele destined to die (based on allele frequency alone).

Randomly choose an allele destined to reproduce (based on allele frequency weighted by growth rate).

Choose a random number to determine if a mutation occurs.

If no mutation occurs, increase the “birth” allele frequency by one else create a new allele with a new randomly chosen growth rate.

Reduce the “death” allele by one.

Return to step 1.

As before, new allelic growth rates were drawn randomly from a normal distribution, N(1, σ), at each environmental change. There is no recombination, current growth rates are independent of previous growth rates, and the growth rates of the mutants and their parents are uncorrelated.

Randomly fluctuating selection, like symmetric fluctuating selection, ensures that rare alleles are more likely to fix than if selectively neutral. This increases the rate of evolution (Figure 7). The coalescent is also affected. Unlike a typical neutral coalescent, many mutations exist only in descendants, the original alleles on which the mutations first arose having gone extinct. Genealogies can appear as ancient balanced polymorphisms, even though they are relatively young (Figure 7).

Figure 7

Examples of genealogies for alleles that are selectively neutral [α ∼ N(1, 0), top] and subject to fluctuating selection [α ∼ N(1, 0.25), bottom] each generation (τ = 1) in a population sized N = 104 with mutation rate u = 10−4. With random fluctuating selection, the genealogy has many more mutations per branch (dots), giving it the appearance of an ancient balanced polymorphism.

Allelic heterozygosity:

Simulations show that rapid random fluctuations in selection (τs < 1) elevate polymorphism. The impact is especially dramatic, with allelic heterozygosity (H = 1 − ∑p2) raised ≥1000-fold above neutral expectations, when the number of new mutant alleles produced each generation (Nu) is low and the selection (Nσ2) is strong (Figure 8A). Random fluctuations in growth rates have little impact on heterozygosity when Nσ2 < 1. The ability of randomly fluctuating selection to increase the absolute number of alleles in the population is relatively modest (Figure 8B). The curvature of the surface suggests that many more alleles might accumulate in larger populations (N > 104). Simulating such large populations is impractical, even on supercomputers, as each individual birth and death must be counted in Moran models. Heterozygosities are raised by flattening the Ewens (1972) sampling distribution (Figure 9) in small and moderately sized populations (N < 104).

Figure 8

The impact of fluctuating selection on polymorphism. (A) Heterozygosity as a function of population size (N), variance in growth rate (σ2), and mutation rate (u) with fitnesses drawn from a normal distribution [N(1, σ)] each generation (τ = 1). The broad red line denotes the neutral expectation. Heterozygosities appear unaffected by fluctuating selection in populations to the right of the thin red line denoting N σ2 = 1. (B) The number of alleles is less affected by fluctuating selection than is heterozygosity in both small (N = 100, red) and moderately large (N = 10,000, tan) populations. The curvature of the surface suggests that large populations (N > 10,000 require huge amounts of CPU time and so were not simulated) may be able to accumulate many more alleles. Parameters were varied as follows: 102 ≤ N ≥ 104, 10−8 ≤ u ≥ 10−1, 10−3 ≤ Nu ≥ 1, τ = 1.

Figure 9

The impact of fluctuating selection on the site frequency spectrum in a sample of 100 alleles drawn from a population of N = 104 with u = 10−4, and σ = 0 (red line, Ewens distribution) and σ = 0.25.

Mixed model: selection and neutrality:

We modified our multiallelic Moran simulation to include neutral mutations as follows:

Randomly choose an allele destined to die (based on allele frequency alone).

Randomly choose an allele destined to reproduce (based on allele frequency weighted by growth rate).

Choose a random number to determine if a mutation occurs.

If no mutation occurs, increase the “birth” allele frequency by one else create a new allele and choose a random number to determine if the mutation is neutral or confers a new randomly chosen growth rate.

Reduce the “death” allele by one.

Return to step 1.

As before, new allelic growth rates were drawn randomly from a normal distribution, N(1, σ), at each environmental change. There is no recombination and current growth rates are independent of previous growth rates. New alleles carrying neutral mutations are assigned the parental fitness in all environments. The growth rates of new alleles with selectable mutations are uncorrelated with those of their parents.

A minority of mutations subject to randomly fluctuating selection can increase the rate of evolution dramatically when the supply of new mutants is limited (small Nu) and selection is strong (large N). Under these circumstances populations rarely have more than two alleles segregating at a time and so the rate of evolution is simply the product of the number of new alleles appearing each generation, Nu, multiplied by the probability of ultimate fixation (1/N for neutral alleles and σ2/2 for selected alleles). The rate of evolution relative to neutral expectations (with r mutations selected and 1 − x neutral) is ∼[(d + d)/u] = (xNu/2 + (1 − x)u)/u = xN/2 + 1 − x. Increases in the rate of evolution are entirely attributable to random fluctuations in selection.

The relative rate of evolution [(d + d)/u] declines as the supply of new mutants (Nu) increases (Figure 10A), an effect attributable to the tendency of clonal interference to slow rates of fixation (it takes much longer to fix an allele when multiple alleles of similar fitness compete). Increases in either population size or the proportion of new mutants subject to randomly fluctuating selection also increase [(d + d)/u]. However, as both intensify clonal interference (both increase the number of segregating alleles in the population), the impact on [(d + d)/u] is smaller than with increases in alone.

Figure 10

Simulations with 10% of new mutants subject to selection. Impact of fluctuating selection on (A) the rate of evolution and (B) the ratio of the d/d ratios of divergence over polymorphism.

The nonsynonymous-to-synonymous substitution ratio for divergence (d/d.div) follows a similar pattern, with clonal interference reducing d relative to d as the supply of new mutants (Nu) increases. Changes in the nonsynonymous-to-synonymous ratio for polymorphism (d/d.poly) are less dramatic. In consequence, the d/d ratio for divergence is raised relative to the d/d ratio for polymorphism (Figure 10B).

Very strong fluctuating selection reduces polymorphism by driving alleles to fixation before the environment changes. Weaker selection, sufficiently strong to maintain selected alleles in a rapidly changing environment, also reduces polymorphism (Figure 11). The phenomenon arises because oscillations in frequency allowing drift to purge alleles of their selectively equivalent variants. The net effect is for fluctuating selection to purge neutral variation.

Figure 11

Strong fluctuating selection can either increase or decrease heterozygosity, depending on the number of new alleles (Nu) entering the population each generation. When Nu is small, fluctuating selection inflates heterozygosity by evening the distribution of selected alleles compared to the neutral expectation. When Nu is large, oscillations in the frequencies of selected alleles cause them to bottleneck each generation, purging them of their neutral variants. Selection varied from σ = 0.1 to 0.25 with Nu varying between 0.1 and 10 using N = 100 (with u varied between 10−3 and 10−1), N = 103 (with u varied between 10−4 and 10−2), and N = 20,000 (with u varied between 10−5 and 10−3), with τ = 1.

Data availability

The authors state that all data necessary for confirming the conclusions presented in the article are represented fully within the article.

Results and Discussion

We explored two models of fluctuating selection in a haploid species. The first, a continuous time model of an infinite population with two alleles subject to periodic changes in fitness, was explored analytically. The second, a probabilistic model of a finite population subject to periodic changes in selection, was explored analytically and by simulation using the Moran framework. Both confirm that, contrary to classical population genetics theory (Kimura 1954; Dempster 1955; Haldane and Jayakar 1963; Gillespie 1972, 1973; Jensen 1973; Karlin and Levikson 1974; Karlin and Liberman 1975; Takahata ; Felsenstein 1976; Maynard Smith 1998; Bell 2008, Huerta-Sanchez ; Uecker and Hermisson 2011; Waxman 2011; Gossmann ; Cvijovic ), fluctuating selection can promote polymorphism in haploids in the absence of frequency-dependent effects on fitness.

The conditions for coexistence in a temporally variable environment (Equations 4 and 5) are similar to those derived by Levene (1953) for coexistence in a spatial model with random dispersal of alleles into habitats at each generation. A common criticism of the Levene model as a mechanism for maintaining polymorphism is that it is robust only with large fitness differences (Hoekstra ). This criticism applies to our model when cycling between two habitats. However, the likelihood of coexistence increases with the number of habitats (Figure 2) to the point that, in a randomly changing environment, a protected polymorphism is virtually guaranteed when the mean fitnesses are equal (Equation 8).

The capacity of fluctuating selection to promote polymorphism is huge. The time to fixation for either of two alleles, with initial frequencies of 0.5 in a population sized N = 5000, is = 5000 Loge2 = 34,657 generations if they are neutral, and = 1021 generations if their fitnesses alternate between 1.4 and 1.4−1 each generation (Figure 3), which is 2000 times the age of the universe, assuming one generation per second. While a cyclical environment with symmetric fitnesses represents an extreme, it illustrates the huge potential of fluctuating selection to maintain polymorphism.

Simulations (Figure 8) show that rapidly fluctuating selection has the capacity to elevate polymorphism, especially when selection is strong and mutation is weak. At least in moderately sized populations, this is achieved by flattening the Ewens sampling distribution rather than by increasing the number of alleles in the population. Indeed, rapidly fluctuating selection can purge neutral variants of selected alleles through repeated bottlenecks (Figure 11). This suggests that high frequency nonsynonymous substitutions flanked by divergent regions depauperate in silent replacements might reflect the action of fluctuating selection (Tiffen ).

Fluctuating selection must be strong to maintain a balanced polymorphism. The symmetric fluctuating selection model needs and the random selection model needs . Although the intensity of selection declines as the square root of the population size, we still need or in a population sized N = 106. The fluctuations in selection must be even larger if the alleles differ in their mean growth rates.

To what extent does fluctuating selection promote diversity in natural populations and communities? Classic work on Drosophila (Dobzhansky 1943), Panaxia dominula (Fisher and Ford 1947; Cook and Jones 1996), Cepaea nemoralis (Cain ), Daphnia (Lynch 1987), and Bison betularia (Saccheri ) suggests that selection can be sufficiently strong, sufficiently variable, and the fluctuations sufficiently rapid to promote polymorphisms. A recent genome-wide survey of a North American population of Drosophila melanogaster (Bergland ) illustrates the huge capacity of fluctuating selection to maintain polymorphisms for many generations. This study identified hundreds of seasonally oscillating single nucleotide polymorphisms, with estimated selection coefficients lying between 0.05 and 0.5 per generation. Some polymorphisms predate the ∼5 million year old divergence with D. simulans. If these results are typical then many polymorphisms might be maintained by fluctuating selection.

What is the probability of establishing a balanced polymorphism in a seasonal environment? The probability is simply the probability that a new allele escapes initial loss to drift; for the symmetric fluctuating selection model and for the random selection model. Assuming genic selection, the probability that a new allele escapes drift to establish a balanced fluctuating polymorphism lies between 0.052/2 = 0.00125 and 0.52 = 0.25. With a mutation rate of 10−9 per site per generation, 1 out of 20 mutations being viable, an average gene with 103 bases, and a genome with 2 × 104 genes, we might expect between 1.25 × 10−6 N and 2.5 × 10−4 N balanced fluctuating polymorphisms to arise per generation. That hundreds of single nucleotide polymorphisms oscillate seasonally in a North American population of D. melanogaster (Bergland ) is by no means excessive, given an effective size of N ≈ 104.

Many mutations of small functional effect likely experience fluctuating selection as fitness optima wobble with changes in the environment. For these alleles, the criterion to fix as if selectively neutral is less stringent than if selection were constant and directional (e.g., rather than in the classic neutral model). A new allele subject to a selection coefficient of s = 0.01 that changes direction each generation in a population sized N = 103 has the same probability of fixing as a neutral allele. This result holds true even when switches in the direction of selection are as infrequent as one in 102 generations (Figure 6A). Now a selection coefficient of 0.01 is readily measured in a competition experiment over a period far greater than 100 generations (Lunzer ). An allele, demonstrated to have changed frequency by natural selection under current conditions, may nevertheless fix as if selectively neutral, were the selection to change direction from time to time. Consequently, conclusions drawn from ecological genetics might occasionally be in conflict with the conclusions drawn from molecular evolution.

Huerta-Sanchez and Gossmann explored the impact of fluctuating selection on molecular evolution using the Fisher–Wright model. They concluded that randomly fluctuating selection increases the probability of ultimate fixation and raises the d/d ratio for divergence over polymorphism in a manner similar to that expected under positive directional selection. We confirmed their results with our Moran simulations (Figure 10). However, our model extends their results in several ways. We relaxed the assumption that selection is strong and mutation is weak to allow many alleles to simultaneously cosegregate in the population. We find that increasing the mutation rate (u) increases polymorphism and thereby intensifies clonal interference, causing the relative rates of fixation and evolution [(d + d)/u] to slow and the d/d ratio for divergence to decline (Figure 10A). We assumed complete linkage so that the d/d ratio for polymorphism receives a boost when Nu > 1 from the large amplitudes experienced by selected alleles, which causes them to lose their neutral variants to drift whenever they become rare (Figure 11). Neither mechanism operates in the models of Huerta-Sanchez and Gossmann because: (1) clonal interference is not possible in a two-allele model, (2) fluctuating selection cannot purge neutral variation unlinked to selected sites, and (3) fluctuating selection reduces polymorphism in the Fisher–Wright model.

In their analyses of fluctuating selection using the Fisher–Wright model, Karlin and Levikson (1974), Takahata , Gillespie (1991), and Cvijovic showed that the frequency of temporal fluctuations in fitness influences the outcome of selection. They found that a new allele would fix as if selectively neutral for τs < 1. For example, simulations show the probability of fixing a new allele in a population sized N = 104 with s = ±0.1 and τ = 1, is U0 ≈ 0.0001 (i.e., a neutral outcome for τs < 1). In the Moran model with rapidly fluctuating symmetric selection, these exact same parameters produce U0 ≈ 0.12/4 = 0.0025 (i.e., a nonneutral outcome for τs < 1, Figure 6A).

Why do the models differ? Allele frequencies are expected to return to their starting values at the end of each cycle in a Fisher–Wright model with symmetric fitnesses (selection is not biased in favor of either allele). Rapid fluctuations in fitness produce such small displacements in allele frequencies that alleles drift as if selectively neutral. Fluctuating selection assumes a different architecture in our model. Rare alleles are favored because both competitors experience more generations when the least fit is most common (which benefits the fitter rare allele) and fewer generations when the most fit is most common (which mitigates selection against the same rare allele when it is less fit). Selection is not expected to return the allele frequencies to their starting values at the end of each cycle even when the changes in fitness are exactly symmetric. Instead, the rare allele is expected to increase in frequency. Selection must be weak () and environmental fluctuations rapid (τs < 1) to overcome this intrinsic bias (Figure 5 and Figure 6).

In 1924, Haldane started his exploration of the genetic basis of evolution using a discrete time nonoverlapping generation model (Haldane 1924). To accommodate sexual reproduction in diploids, he abandoned certain key demographic and ecological processes (e.g., continuous reproduction in bounded populations). His simplification was brilliant, and perfectly suited for his goal of studying genes in populations (Haldane 1932). There can be no greater praise for this pioneering enzymologist (Haldane 1930) than to note that his is the approach we all use to study evolutionary genetics today. Yet, abandoning continuous reproduction in bounded populations came with a hidden cost. Haldane’s simplification, so ideal for studying directional selection in constant environments, forces the generations to march in lock step with seasonal changes. His model does not accommodate frequency-dependent changes in the number of allelic doublings (generations) per season. Both in our experiments (Yi and Dean 2013) and in our models, this frequency-dependent slippage in generations promotes diversity.

In our random selection model (equivalent to a diploid model with genic selection) evolutionary outcomes are the emergent properties of three inescapable processes: mutation, geometric growth, and resource depletion (to place a bound on population size). The model brings selection and neutrality under the same umbrella to reveal a middle ground, where weakly selected alleles contribute both to polymorphism and fix as if selectively neutral. Unlike the classical population genetic models of fluctuating selection, our model has direct empirical support (Figure 1 shows selection can protect a polymorphism in a haploid species) and predicts the higher d/d ratios for divergence than for polymorphism commonly observed among closely related taxa (Hughes 2008; Wagner 2008). We suggest that our random selection model provides a realistic alternative to the Fisher–Wright model for describing fluctuating selection in finite populations.