Open Access

The competition between simple and complex evolutionary trajectories in asexual populations

BMC Evolutionary Biology201515:55

Received: 25 November 2014

Accepted: 11 March 2015

Published: 26 March 2015



On rugged fitness landscapes where sign epistasis is common, adaptation can often involve either individually beneficial “uphill” mutations or more complex mutational trajectories involving fitness valleys or plateaus. The dynamics of the evolutionary process determine the probability that evolution will take any specific path among a variety of competing possible trajectories. Understanding this evolutionary choice is essential if we are to understand the outcomes and predictability of adaptation on rugged landscapes.


We present a simple model to analyze the probability that evolution will eschew immediately uphill paths in favor of crossing fitness valleys or plateaus that lead to higher fitness but less accessible genotypes. We calculate how this probability depends on the population size, mutation rates, and relevant selection pressures, and compare our analytical results to Wright-Fisher simulations.


We find that the probability of valley crossing depends nonmonotonically on population size: intermediate size populations are most likely to follow a “greedy” strategy of acquiring immediately beneficial mutations even if they lead to evolutionary dead ends, while larger and smaller populations are more likely to cross fitness valleys to reach distant advantageous genotypes. We explicitly identify the boundaries between these different regimes in terms of the relevant evolutionary parameters. Above a certain threshold population size, we show that the probability that the population finds the more distant peak depends only on a single simple combination of the relevant parameters.


Epistasis Rugged fitness landscape Fitness valley


In an adapting population, evolution often has the potential to follow many distinct mutational trajectories. In order to predict how the population will adapt, we must understand how evolution chooses among these possibilities. Many experimental and theoretical studies have analyzed this question, focusing primarily on the simple case where epistasis is absent, so that each mutation has some fixed fitness effect [1-6]. This work can explain the probability that a given mutation will fix as a population adapts, as a function of its fitness effect, the population size, mutation rate, distribution of fitness effects of other mutations, and other parameters of the evolutionary process.

However, the fitness effect of a mutation often depends on the genetic background in which it occurs. A particularly interesting form of this phenomenon, sign epistasis, occurs when several mutations are individually neutral or deleterious but their combination is beneficial [7]. Sign epistasis has been observed repeatedly in experiments [8-13], and plays a central role in the evolution of complex phenotypes that involve multiple interacting components. When sign epistasis is present, adaptation can involve passing through genotypes of lower fitness — i.e. a population may have to cross a fitness valley or plateau. Thus the fate of a mutation depends not only on its fitness, but also on its adaptive potential [14].

Several recent theoretical studies have analyzed the evolutionary dynamics of fitness valley crossing [15-20]. This work has focused on calculating the rate at which adapting populations cross a valley or plateau, in the absence of any other possible mutational trajectories. However, individually beneficial mutations may often compete with more complex evolutionary trajectories. We must then ask how likely evolution is to eschew the immediately uphill paths, and instead cross valleys or plateaus to reach better but less accessible genotypes. In other words, when the fitness landscape is rugged, we wish to understand whether evolution will take the more “farsighted” path to reach distant advantageous genotypes, rather than a “greedy” trajectory that fixes immediately beneficial mutations regardless of whether these may lead to evolutionary dead ends.

In this article, we analyze this evolutionary choice between immediately beneficial mutations and more complex mutational trajectories that ultimately lead to higher fitness. We calculate the probability that an adapting population will follow each type of competing trajectory, as a function of the population size, mutation rates, and selection pressures. We focus on asexual populations, where the only way for a population to acquire a complex adaptation is for a single lineage to acquire each mutation in turn. Our analysis is similar in spirit to earlier work which also considered the tradeoff between short-term and long-term fitness advantages [21-24]. However, these earlier studies dealt with competition between different strictly uphill or neutral paths, and considered the case where the less beneficial initial mutation led to better long-term evolutionary opportunities. In contrast, our analysis describes the competition between uphill mutations and more complex trajectories. While these two cases can be qualitatively similar in very small populations, they lead to very different dynamics in larger populations where the sign of the effect of the intermediate mutation can play a crucial role.

Our results show that population size has a crucial impact on how “farsighted” evolution can be. This dependence is not monotonic: evolution at intermediate population sizes is most “greedy”, while both larger and smaller populations are more likely to eschew uphill paths in favor of complex trajectories. In large populations, our results show that a single parameter reliably predicts the extent of this evolutionary “foresight” across a wide range of parameters. Finally, we describe how our analysis can be generalized to predict how evolution will choose among even more complex trajectories, such as broad fitness valleys with multiple intermediate genotypes, and we discuss evolution in genotype spaces with many possible evolutionary paths.


We are interested in how a population makes an evolutionary choice when confronted with multiple possible mutational trajectories. Specifically, we focus on the extent to which adaptation proceeds by crossing fitness valleys rather than acquiring immediately beneficial (uphill) mutations. Of course, the relative frequency of valley crossing will depend on the number of available fitness valleys, their depth, and the fitness advantage of the multiple-mutants, as well as the distribution of fitness effects (DFE) of the uphill mutants. Our goal is to understand how the prevalence of valley crossing depends on these factors.


Throughout most of this article, we consider the simplest context in which we can address this question: the choice between a single uphill path and a single fitness valley. Specifically, we consider a haploid asexual population of constant size N which can either acquire an uphill mutation (u) that confers an immediate fitness advantage s u , or alternatively acquire a deleterious fitness valley intermediate (i) with fitness deficit δ i on which background a double-mutant (v) with fitness s v >s u can arise. This scenario is illustrated in Figure 1. We also consider the case of a fitness plateau, where δ i =0.
Figure 1

Model and characteristic trajectories. (a) The model to study fitness valley crossing prevalence. The population starts as wild type (w), and then acquires uphill mutations (u) at rate μ u that confer an immediate fitness advantage s u , and acquires deleterious fitness valley intermediates (i) at rate μ i with fitness deficit δ i on which background double-mutants (v) with fitness s v >s u arise at rate μ v . (b)-(e) The four main forms of fitness valley crossing. (b) Small populations are characterized by low genetic diversity and strong genetic drift, leading sequential fixation of intermediates to dominate the dynamics. (c) For larger populations, genetic diversity is maintained longer, and double mutants will tend to arise on transient single-mutant backgrounds, in a process known as stochastic tunneling. (d) If the drift time is small compared to the maximal rate of change in background fitness, we can approximate the drift time of the intermediate by its expectation, dramatically simplifying the mathematical analysis. (e) For very large populations, we can treat single-mutants deterministically, in a process dubbed semi-deterministic tunneling.

Because we are interested in the evolutionary choice between competing mutational trajectories, we assume that these two trajectories are mutually exclusive (i.e. the mixed genotypes ui and uv are strongly deleterious), so that only one genotype (either u or v) can eventually fix in the population. As a measure of evolutionary foresight, we analyze the probability that the double-mutant v fixes as a function of the relevant mutation rates, selection coefficients, and population size. In some situations, we could imagine that after either genotype u or v fixes, another set of competing potential trajectories become available. In this case, our analysis predicts the long-term relative ratio of fixed uphill versus valley-crossing mutations. In the Discussion, we consider how this model can be extended to the situation where there are many different competing uphill paths and valleys, and to broader fitness valleys involving multiple intermediate genotypes.


In addition, we compare our analytical predictions for valley crossing probability to Wright-Fisher simulations. Each simulated population was evolved until either the uphill genotype or valley-crossing genotype fixed. Valley crossing probabilities were then inferred from the number of trials in which the valley-crossing genotype fixed, out of 1000 trials per parameter set.


In the absence of the uphill genotype, fitness valley crossing can be modeled as a homogeneous Poisson process with rates as calculated by [17]. In small populations, the primary role of the uphill genotype is to introduce an effective time limit on this process: once an uphill mutation destined to survive drift first occurs, it very quickly fixes, leading to the extinction of the wild-type. The probability of valley-crossing can thus be calculated as the probability that the intermediate i fixes before the uphill genotype u. An example of this is shown in Figure 1b.

In larger populations, the dynamics are more complex, as illustrated in Figure 1c. Rather than leading to a single cutoff time for valley-crossing to occur, the single-mutant occurs and gradually increases in frequency. This leads to a decline in the size of the wild-type background on which intermediate and valley-crossing mutants can arise, and a corresponding increase in the mean fitness of the population (Figure 1c). These effects gradually reduce the rate at which intermediates are produced, and make these intermediates effectively more deleterious relative to the mean fitness. These factors reduce the rate of the valley-crossing process. Thus valley-crossing becomes an inhomogenous Poisson process, with a rate that depends on the random appearance time T u of the uphill mutant.

In general, these effects of interference and tunneling are complex. However, the analysis becomes simpler in two specific regimes. When the expected drift time of the intermediate genotype is short, we can neglect the changing background fitness due to the uphill mutant during this drift time (Figure 1d). Alternatively, for very large populations (N μ>1), the Poisson process approximation breaks down and both uphill and intermediate mutations can be treated deterministically (Figure 1e), and only the valley-crossing genotype must be treated stochastically.

These various regimes are illustrated in Figure 2. We now analyze each in turn, assuming weak selection (s j <1) for all genotypes throughout. Taken together, this provides a complete picture of the probability that evolution will eschew the immediately uphill path in favor of the more complex adaptation.
Figure 2

Regimes of valley crossing. Phase plot summarizing the different regimes of fitness valley crossing. (a) The regime boundaries in terms of population size N and uphill fitness s u . (b) The boundaries in terms of population size N and intermediate deleteriousness δ i .

Small populations

When the population size is small enough that the probability of stochastic tunneling is very low, the population is generally clonal or nearly clonal, and moves in Markovian jumps between neighboring genotypes. The transition between genotypes i and j occurs at rate
$$ r_{i j} = N \mu_{i j} \pi_{i j}, $$
where π i j is the probability that a single j mutant will give rise to a lineage that fixes, given by the standard formula,
$$ \pi_{i j} = \frac{1-e^{-2(s_{j}-s_{i})}}{1-e^{-2N(s_{j}-s_{i})}}. $$
We refer to this as the sequential fixation regime. Because we are considering neutral and weakly deleterious intermediates, we account for the possibility of back-mutation to the wild type if the intermediate fixes. Therefore, the process can be modeled as an absorbing states Markov chain, where the wild type and intermediates act as transient states, and the uphill and double-mutant genotypes act as absorbing states. From elementary Markov chain theory, we find
$$ {\fontsize{8.8pt}{9.6pt}\selectfont{\begin{aligned} {}{P_{\text{cross}}} \,=\, \frac{r_{w i} r_{i v}}{r_{w i} r_{i v} \!+ {r_{w u}} ({r_{i w}} \!+ r_{i v})}\,=\, \left[\!1 \,+\, \left(\! \frac{\pi_{wu}}{\pi_{wi}} \!\right)\!\left(\! \frac{{\mu_{u}}}{{\mu_{i}}}\! \right) \!\left(\!1\,+\,\frac{{\mu_{i}} \pi_{iw}}{{\mu_{v}} \pi_{iv}}\! \right)\!\right]^{\!-1}\!. \end{aligned}}} $$

As the population size increases, π wu →2s u and π wi →0, so \(\frac {\pi _{\textit {wu}}}{\pi _{\textit {wi}}} \rightarrow \infty \), and P cross→0. Thus we find that within the sequential fixation regime, larger population sizes are less likely to cross fitness valleys.

Stochastic tunneling

For large populations, the probability that deleterious intermediates will fix declines drastically, and successful double mutants will instead arise on the unfixed single-mutant background in a process known as stochastic tunneling [16]. This transition occurs when
$$ N > \frac{1}{2 \delta_{i}} \log \left[ 1 + \frac{\exp (2 \delta_{i}) - 1}{p_{v}} \right], $$
where p v is the probability that the intermediate lineage survives drift long enough to give rise to an ultimately successful double mutant lineage (we will explicitly calculate this probability below). We can then model the appearance of an intermediate mutant lineage destined to give rise to a double mutant lineage as a Poisson process. The rate λ v at which these lineages appear is given by the rate at which intermediate mutations arise times the probability of success of the lineage, integrated over the drift time t d after appearance of the single-mutant intermediate:
$$ {\lambda_{v}}=N_{\text{wt}} {\mu_{i}} \int_{0}^{\infty} \frac{\partial {p_{i}}(t_{d})}{\partial t_{d}} dt_{d}. $$
Here N wt is the wild-type population size, and p i (t d ) is the cumulative probability that a single-mutant lineage will give rise to a successful double-mutant lineage by time t d after it appears. This probability is given by [17]:
$$ {\fontsize{8.3pt}{9.6pt}\selectfont{\begin{aligned} {}{p_{i}} (t_{d})\! =\! 2 \left(\! \frac{(a_{+}\! -1)(1-a_{-})(1-\exp[-(1-{\delta_{i,\text{eff}}})(a_{+}-a_{-})t_{d}]\!)}{a_{+} -1+(1-a_{-})(1-\exp[-(1-{\delta_{i,\text{eff}}})(a_{+}-a_{-})t_{d}]\!)} \!\right), \end{aligned}}} $$
$$ \begin{aligned} {}a_{\pm} = \frac{2 -{\delta_{i,\text{eff}}} - {\mu_{v}} {s_{v,\text{eff}}} \pm \sqrt{({\delta_{i,\text{eff}}}+{\mu_{v}} {s_{v,\text{eff}}})^{2}+4{\mu_{v}} {s_{v,\text{eff}}}}}{2(1-{\delta_{i,\text{eff}}})}, \end{aligned} $$

and δ i,eff and s v,eff are the fitnesses of the intermediate and valley-crossing genotypes relative to the (time-dependent) mean population fitness. These effective fitnesses and N wt depend on the background at time t d in a way we must now consider. The background in turn is determined by the frequency f u (T v +t d ) of the uphill genotype at time t d after the appearance of the first single-mutant intermediate lineage destined for success at time T v (Figure 1c). Thus the uphill genotype frequency sets the fitness background on which valley-crossing probabilities are determined.

To calculate the relevant effective parameters, we note that the appearance of uphill lineages destined for success can be modeled as a Poisson process. Moreover, because the valley-crossing genotypes make up a tiny fraction of the population (unless the double-mutant has already established), we can treat the genetic background on which these uphill lineages appear as essentially fixed. Therefore, the first uphill lineage destined to survive genetic drift will appear at time T u , distributed exponentially with rate
$$ {\lambda_{u}}=N {\mu_{u}} \pi_{wu} \approx N {\mu_{u}} (2 {s_{u}}). $$
Once a successful uphill lineage appears, we assume it establishes in time τ est=γ e /(2s u ), where γ e ≈.577 is the Euler-Mascheroni constant [25], and then sweeps deterministically according to
$$ {f_{u}}(\hat{t}) = \frac{1-\exp\left[-({\mu_{u}}+{s_{u}})\hat{t}\,\right]} {1+({s_{u}}/{\mu_{u}})\exp\left[-({\mu_{u}}+{s_{u}})\hat{t}\,\right]}, $$

where \(\hat {t}\equiv t-{T_{u}}-\tau _{\text {est}}\) is the time after establishment of the uphill mutant.

Conditioning on the appearance time T u , we can thus work out our effective parameters
$$\begin{array}{@{}rcl@{}} N_{\text{wt}}(t \, | \, {T_{u}}) & = & N(1-{f_{u}}(t-{T_{u}}-\tau_{\text{est}})) \end{array} $$
$$\begin{array}{@{}rcl@{}} \delta_{i,\text{eff}}(t+t_{d} \, | \, {T_{u}}) & = & -\delta_{i}\!+{f_{u}}(t\!+t_{d}\!-{T_{u}}\!-\tau_{\text{est}}){s_{u}} \end{array} $$
$$\begin{array}{@{}rcl@{}} s_{v,\text{eff}}(t+t_{d} \, | \, {T_{u}}) & = & s_{v}\!-{f_{u}}(t\!+t_{d}\!-{T_{u}}-\tau_{\text{est}}) {s_{u}}. \end{array} $$

These effective parameters encompass the two main effects of the sweeping uphill mutation on the valley crossing probability: the first represents the declining wild-type background on which new mutations can arise, and the remaining two represent the decreasing relative fitness of the valley-crossing lineage.

We are interested in the probability that a double-mutant lineage destined for success appears before the uphill genotype fixes. Integrating over all possible appearance times T u , this is given by:

$$\begin{array}{@{}rcl@{}} {P_{\text{cross}}} &\! =\! & \int_{0}^{\infty} dt_{u} \left({\lambda_{u}} e^{-{\lambda_{u}} t_{u}} \right) \times\\ & & \left[ \!1 - \exp\left[\!-\int_{0}^{\infty} dt \left(N_{\text{wt}}(t,{t_{u}}) {\mu_{i}} \int_{0}^{\infty} dt_{d} \frac{\partial {p_{i}}(t+t_{d},t_{u})}{\partial t_{d}} \right)\right]\right]. \end{array} $$

This integral is a complete solution for the probability of valley-crossing in the stochastic tunneling regime, provided that the population size is small enough that the Poisson process approximation above holds. Although it does not have a simple closed-form solution, we can easily evaluate the integral numerically. Alternatively, there is a simple and relevant parameter regime in which background fitnesses change slowly. We now consider this case, and show that it allows us to evaluate our expression for the valley-crossing probability explicitly. In a later section below, we turn to the alternative case where the Poisson process approximation breaks down, and we can instead treat all single-mutants deterministically; the valley-crossing probability also simplifies considerably in this very large population regime.

Slowly changing background fitness

One of the main complications of Equation (13) is the integral over possible drift times t d , which reflects the increasing effective deleteriousness of the intermediate genotype as the uphill genotype sweeps to fixation and increases the mean fitness of the background population. However, when s u is small or δ i is large, this background fitness changes slowly compared to the intermediate drift time. In this case, we can treat the background during intermediate drift as effectively constant (Figure 1d). This eliminates the need for an integral over t d , since the time-dependent probability of success of a single-mutant at time t is fully determined by f u (t+t d )≈f u (t).

We can further simplify the analysis if we treat the probability of crossing the valley as a function of two probabilities: the probability P v,1 that the first successful valley-crossing lineage appears before the first successful uphill lineage establishes, and the probability P v,2 that a successful valley-crossing lineage appears after the uphill mutant establishes:
$$ {P_{\text{cross}}} = P_{v,1} + \left(1-P_{v,1} \right) P_{v,2}. $$
The calculation of P v,1 takes place on a purely wild-type background, so we can use
$$ {\lambda_{v}} = N{\mu_{i}} p_{v}, $$
where p v , the probability that the intermediate lineage survives drift long enough to give rise to an ultimately successful double mutant lineage, is given by [17]:
$$ {}p_{v} = -{\delta_{i}}+\sqrt{{\delta_{i}}^{2}+4{\mu_{v}} {s_{v}}} \approx \left\{ \begin{array}{ll} 2\sqrt{{\mu_{v}} {s_{v}}} & \text{if}~~ {\delta_{i}} \ll 2\sqrt{{\mu_{v}} {s_{v}}} \\ 2{\mu_{v}} {s_{v}}/{\delta_{i}} & \text{if}~~ {\delta_{i}} \gg 2\sqrt{{\mu_{v}} {s_{v}}}. \end{array} \right. $$
Meanwhile, λ u is unchanged from the original analysis. P v,1 is determined by a race between these two exponential random variables. Using basic properties of the exponential, the probability that the double-mutant appears before the uphill genotype establishes is therefore
$$\begin{array}{*{20}l} P_{v,1} &= \left\{ \begin{array}{ll} \left(\frac{{\lambda_{v}}}{{\lambda_{v}}+{\lambda_{u}}}\right) e^{-{\lambda_{u}} \Delta \tau} & \text{if}~~ \Delta \tau>0\\ 1-\left(\frac{{\lambda_{u}}}{{\lambda_{v}}+{\lambda_{u}}}\right) e^{-{\lambda_{v}} \Delta \tau} & \text{if}~~ \Delta \tau<0, \end{array} \right. \end{array} $$
where Δ τ represents the difference between the mean drift time τ drift of the valley-crossing lineage and the mean establishment time τ est of the uphill lineage,
$$ \Delta \tau = \tau_{\text{drift}} - \tau_{\text{est}}. $$
Here τ est is as given above, and from [17], we can approximate the drift time as
$$\begin{array}{*{20}l} \tau_{\text{drift}} &= \left\{ \begin{array}{ll} \log 2/\sqrt{{\mu_{v}} {s_{v}}} &\text{if}~~ {\delta_{i}} \ll 2\sqrt{{\mu_{v}} {s_{v}}}\\ 1/{\delta_{i}} &\text{if}~~ {\delta_{i}} \gg 2\sqrt{{\mu_{v}} {s_{v}}}. \end{array} \right. \end{array} $$
If the successful uphill lineage establishes, the first successful valley-crossing lineage still has a chance to appear and outcompete it, albeit on a declining wild-type background. Thus for P v,2 we get a similar integral as in the original analysis. However, since we are approximating p i as constant, the rate of successful lineage generation simplifies to
$$ {\fontsize{8.7pt}{9.6pt}\selectfont{\begin{aligned} {}\hat{\lambda}_{v}(t) &= {\mu_{i}} \, N_{\text{wt}} \; p_{i} (\delta_{i,\text{eff}},s_{v,\text{eff}})\\ &= {\mu_{i}} \, N(1\!-{f_{u}}) \left(\!-({\delta_{i}}\,+\,{s_{u}} {f_{u}})\,+\,\sqrt{({\delta_{i}}\,+\,{s_{u}} {f_{u}})^{2}\,+\,4{\mu_{v}} ({s_{v}}\,-\,{s_{u}} {f_{u}})}\!\right)\!. \end{aligned}}} $$
Integrating and assuming mutation rates are small compared to selection pressures, we find:
$$ \int_{0}^{\infty} dt \hat{\lambda}_{v}(t) \approx \left(\frac{\log{N{s_{u}}}}{{s_{u}}}\right) N {\mu_{i}} {p_{i}} = \left(\frac{\log{N{s_{u}}}}{{s_{u}}}\right) {\lambda_{v}}. $$
Combining these results, we find
$$\begin{array}{*{20}l} {P_{\text{cross}}} &= P_{v,1} + \left(1-P_{v,1}\right) \left(1-e^{-\left(\frac{\log{N{s_{u}}}}{{s_{u}}}\right) {\lambda_{v}}}\right) \notag \\ &= 1 - \left(1-P_{v,1}\right) e^{-\left(\frac{\log{N{s_{u}}}}{{s_{u}}}\right) {\lambda_{v}}}. \end{array} $$
We expect this result to be valid provided that δ i is effectively constant over the expected drift time. This will hold when
$$\begin{array}{*{20}l} {s_{u}} &\ll \left\{ \begin{array}{ll} 2 \sqrt{\frac{2}{\log 2}} \sqrt{{\mu_{v}} {s_{v}}} & \text{if}~~ {\delta_{i}} \ll 2\sqrt{{\mu_{v}} {s_{v}}}\\ 2 {\delta_{i}} & \text{if}~~ {\delta_{i}} \gg 2\sqrt{{\mu_{v}} {s_{v}}} \end{array} \right.. \end{array} $$

Semi-deterministic tunneling

We now consider the case where N μ i >1 and N μ u >1, and hence the Poisson process approximation used to derive Equation (13) breaks down. Fortunately, in this regime the number of single-mutant intermediates and uphill mutants in the population are well approximated by their deterministic expectation (Figure 1e). Thus the only random variable is the appearance time of the first successful double-mutant lineage, which occurs with rate
$$ {\lambda_{v}}(t) = N {f_{i}} {\mu_{v}} \pi_{iv} \approx 2 N {f_{i}} {\mu_{v}} ({s_{v}} - {s_{u}} {f_{u}}). $$
Because intermediates never make up a large portion of the population, f u is unaffected by f i , and hence is still given by Equation (9). We can then approximate the frequency f i of the single-mutant intermediates using mutation-selection balance with a declining wild-type population:
$$ {f_{i}} = {f_{i}}^{*} (1 - {f_{u}}), $$
where f i gives the independent deterministic dynamics (mutation-selection balance) of single mutants on the wild-type background, and (1−f u ) is the size of the wild-type background. It is useful to transform this into an integral in the frequency domain:
$$ \int_{0}^{\infty} {\lambda_{v}}(t) dt = {\int_{0}^{1}} {\lambda_{v}}({f_{u}}) \left(\frac{\partial {f_{u}}}{\partial t}\right)^{-1} d{f_{u}}. $$
We note from equation 9 that
$$ \frac{\partial {f_{u}}}{\partial t} = {\mu_{u}} (1-{f_{u}}) + {s_{u}} {f_{u}}(1-{f_{u}}), $$
$$ t({f_{u}})=\frac{1}{{\mu_{u}}+{s_{u}}} \log \left(\frac{1+({s_{u}}/{\mu_{u}}){f_{u}}}{1-{f_{u}}}\right). $$
We further approximate
$$ {f_{i}}^{*} = \left\{ \begin{array}{ll} {\mu_{i}} t &\text{if}~~ {\delta_{i}} \ll 2\sqrt{{\mu_{v}} {s_{v}}}\\ \frac{{\mu_{i}}}{{\delta_{i}}} (1-\exp(-{\delta_{i}} t)) &\text{if}~~ {\delta_{i}} \gg 2\sqrt{{\mu_{v}} {s_{v}}}. \end{array} \right. $$
Combining these expressions and assuming μ u s u , we find
$$ \int_{0}^{\infty} {\lambda_{v}}(t) dt = 2N{\mu_{i}} {\mu_{v}} {s_{v}} \gamma /{s_{u}}, $$
$$ {\fontsize{8.1pt}{9.6pt}\selectfont{{}\gamma\! =\! \left\{\!\! \begin{array}{ll} {s_{u}}^{\!\!-1} \left(\tfrac{1}{2} \log ({s_{u}}/{\mu_{u}})^{2} \!- \!({s_{u}}/{s_{v}}\!\right) \log ({s_{u}}/\!{\mu_{u}}) \,+\, \pi^{2}/6) &\text{if}~~ {\delta_{i}} \!\ll 2\sqrt{{\mu_{u}} {s_{v}}}\\ {\delta_{i}}^{\!-1} \left(\log ({s_{u}}/{\mu_{u}})\! - \!({s_{u}}/{\delta_{i}}) (1\,-\,({s_{u}}/\!{\mu_{u}})^{-{\delta_{i}}/{s_{u}}}) \right) &\text{if}~~ {\delta_{i}} \!\gg 2\sqrt{{\mu_{u}} {s_{v}}}. \end{array} \right.}} $$
This gives
$$ \begin{aligned} {P_{\text{cross}}} = 1&-\exp \left[-\int {\lambda_{v}}(t) dt\right] =1\\&- \exp \left[\frac{-2N{\mu_{i}} {\mu_{v}} {s_{v}} \gamma}{{s_{u}}}\right] = 1- \exp \left[-2N{\Gamma}\right], \end{aligned} $$
where we have defined the useful quantity
$$ {\Gamma} \equiv \frac{{\mu_{i}} {\mu_{v}} {s_{v}} \gamma}{{s_{u}}} \,. $$

Thus we see that in very large populations, the probability of valley-crossing depends in a simple way on the single composite parameter Γ. The form of this composite parameter depends crucially on whether the fitness cost of the intermediate genotype is large or small, as defined in Equation (31).


Our results have shown that the size of a population strongly influences how “farsighted” it can be. In small populations, genetic drift is strong relative to selection, so the evolutionary dynamics proceeds by sequential fixation. Since fixation of a deleterious intermediate becomes less likely in larger populations, this means that increasing the population size initially decreases the relative influence of fitness-valley crossing. However, as the population size increases further, beneficial mutations take longer to fix, maintaining diversity in the population and allowing double mutants to stochastically tunnel on the declining wild-type background [15-17]. Together, these effects lead to a non-monotonic relationship between population size and the probability that evolution will favor the complex adaptation over the directly uphill path, as illustrated in Figure 3a. This nonmonotonic dependence on population size is similar in spirit to the results of earlier work analyzing evolution on epistatic landscapes in the absence of fitness valleys [21-23].
Figure 3

Valley crossing probability. (a) Simulation results for μ u =5×10−6, μ i =μ v =5×10−5, δ i =0, and s v =.07. The black vertical dashed lines indicate the boundaries between sequential fixation, stochastic tunneling, and semi-deterministic tunneling. Markers represent inferred valley-crossing probabilities from 1000 simulations per point. Lines represent theoretical predictions in each regime: in the stochastic tunneling regime, dashed lines represent the slowly changing background fitness approximation, and solid lines represent the full integral solution Equation (13). The color of the line indicates the uphill fitness s u . (b) Crossing probability for populations in the semi-deterministic regime across a wide range of parameters, plotted against the predictive parameter N Γ. Filled markers represent deleterious intermediates (\({\delta _{i}} = 10\sqrt {{\mu _{v}} {s_{v}}}\)), while open markers represent neutral intermediates (δ i =0).

It is interesting to note that P cross does not immediately begin to rise with the onset of tunneling. Instead, the dependence is more complex, as a consequence of the tradeoff between increasing mutation rates and fixation times. Nevertheless, for populations in which the transition to valley-crossing behavior occurs in the semi-deterministic regime, we can derive a simple expression for the threshold size at which the population will tend to cross valleys with probability P cross. A straightforward inversion of (32) gives
$$ N= \frac{-\log \left[1-{P_{\text{cross}}}\right]}{2 {\Gamma}}, $$

valid for large population sizes in the semi-deterministic regime. Thus in this regime the threshold size above which a population exhibits a given degree of foresight (i.e. has a particular P cross) depends only on Γ. To illustrate this, in Figure 3b we show P cross as a function of N Γ for a variety of simulations across different values of μ u ,μ i ,μ v ,s u ,δ i and s v . It is clear that even across a wide parameter range, N Γ is a reliable predictor of valley crossing probability in the semi-deterministic limit.

Multiple intermediates and evolutionary predictability

Recently, Szendro et al. [26] simulated evolution across a wide range of population sizes on an experimentally-derived epistatic fitness landscape, finding that “evolutionary entropy” (i.e. unpredictability over all possible outcomes, given by \(S = - \sum p_{j} \log p_{j}\)) varied nonmonotonically with population size. Specifically, these authors found that entropy initially decreased with N above a characteristic population size N1/μ, before increasing again above a second characteristic size N1/μ 2; they argued that these points were related to the supply rate of single and double mutants respectively. Our analysis is consistent with these results. For example, the increase in entropy at N1/μ 2 found by [26] corresponds to our result that valley-crossing begins to significantly influence evolution when N1/Γ: this is approximately proportional to 1/μ 2, albeit with an additional log dependence on μ from the γ factor that would be harder to observe experimentally.

We find related behavior if we extend our analysis to valleys with more intermediates. As a simple example, we consider a fitness landscape where we add deleterious intermediates u 0 and v 0 (each occurring at rate μ 0 and leading to fitness cost δ i ) to the uphill and valley-crossing branches respectively, so that we now have competition between a single-intermediate valley and a two-intermediate valley. In a large enough population, mutation-selection balance will ensure that
$$ N_{u_{0}} = N_{v_{0}} = \frac{N \mu_{0}}{{\delta_{i}}}. $$
If we assume that these sub-populations are large enough that double-mutants behave deterministically, then we find the crossing probability obeys
$$ {}-\log \left[1-{P_{\text{cross}}}\right] = 2 \left(\frac{\mu_{0}}{{\delta_{i}}}\right) N {\Gamma} = 2 \left(\frac{{s_{v}} \gamma}{{s_{u}} {\delta_{i}}} \right) N \mu_{0} {\mu_{i}} {\mu_{v}}. $$

Thus our analysis of valleys with multiple intermediates suggests that the N μ 2 entropy peak is not unique: as the population grows larger, there should be entropy peaks corresponding to foresight across valleys with increasing numbers of intermediates. The emergence of such second peaks has been observed in simulations [23], and our model offers a quantitative outline of where such peaks should occur given the relevant evolutionary parameters. In the example above, for instance, there would be an entropy peak approximately proportional to N μ 3, and in general, we could expect entropy peaks at points approximately proportional to N μ n for (n−1)-intermediate valleys. In practice, however, the semi-deterministic approximation will break down for any sizable number of intermediates, unless the population size is unrealistically large.

Many paths

Throughout this paper, we have assumed the presence of a single uphill mutation and fitness valley. We now consider how our analysis can be extended to predict how evolution chooses among many such possible mutational trajectories.

In small populations that are in the sequential fixation regime, we simply add additional transient transition matrix elements representing different mutations, with the uphill mutations transitioning to the uphill absorbing state, and similarly for the valley-crossing mutations. When stochastic tunneling is important, we must instead add the rates of single valley-crossing mutants to get a total rate of
$$ {\Lambda_{v}} = \sum_{v} {\lambda_{v}}, $$
and similarly find the total rate of uphill mutants that are destined to survive drift,
$$ {}{\Lambda_{u}} = \sum_{u} {\lambda_{u}} = \sum_{u} N {\mu_{u}} \pi({s_{u}}) = \int N U_{b} \pi({s_{u}}) \rho({s_{u}}) d {s_{u}}, $$
where in the last equality we have replaced a discrete collection of uphill mutations with a continuous distribution of uphill fitness effects ρ(s), and a discrete collection of mutation rates μ u with a total beneficial mutation rate U b . This is valid as long as Λ u 1. Once an uphill mutation destined to survive drift occurs, the probability that it has fitness s is given by the ratio between its partial rate and the total rate; formally, the probability density is given by:
$$ f(s) = N {\mu_{u}} {\pi(s)} \rho(s) / {\Lambda_{u}}. $$
Using these expressions, we can integrate our results from the analysis over all possible trajectories. However, we note that if there are a large number of weakly beneficial mutations, it is possible the first successful lineage to appear will be outcompeted by a stronger uphill mutation that arises later but fixes first. Our analysis applies provided that we consider only uphill mutations that reach a significant portion k of the population before a new, more fit uphill mutant is expected to be produced: i.e.
$$ {\Lambda_{u}} \tau_{k} = \left(\int_{s_{\text{cutoff}}}^{\infty} N \mu_{s} {\pi(s)} \rho(s) ds\right)(\tau_{k}) <1, $$

where τ k is the expected time for a single-mutant destined for success to make up frequency k of the population. This is consistent with our intuition that as the population size grows larger, we increasingly expect the mutations of largest effect to dominate the dynamics.


Using a simple three-locus fitness landscape model (Figure 1a), we identified several regimes of valley-crossing with qualitatively different behavior (Figure 2). By examining the behavior in each of these regimes in turn, we found that the probability of valley-crossing has a complex, non-monotonic dependence on population size (Figure 3a), and identified a parameter Γ that reliably predicts the population size at which valley crossing becomes preferred (Figure 3b). Finally, we showed how these results can be extended to fitness valleys with more intermediates, and to fitness landscapes with many possible evolutionary trajectories, as is the case in most naturally occurring populations.



Distribution of fitness effects



We thank Benjamin Good, Sergey Kryazhimskiy, Elizabeth Jerison, and other members of the Desai lab for useful discussions, and Katya Kosheleva for help with the figures. This work was supported by the Harvard HCRP, Harvard Herchel Smith, and Harvard PRISE programs (I.E.O.), and the James S. McDonnell Foundation, the Alfred P. Sloan Foundation, the Harvard Milton Fund, grant PHY 1313638 from the NSF, and grant GM104239 from the NIH (M.M.D.). Simulations in this paper were performed on the Odyssey cluster supported by the Research Computing Group at Harvard University.

Authors’ Affiliations

Department of Organismic and Evolutionary Biology, Department of Physics, and FAS Center for Systems Biology, Harvard University


  1. Perfeito L, Fernandes L, Mota C, Gordo I. Adaptive mutations in bacteria: High rate and small effects. Science. 2007; 317(5839):813–5.View ArticlePubMedGoogle Scholar
  2. Tenaillon O, Rodríguez-Verdugo A, Gaut RL, McDonald P, Bennett AF, Long AD, et al.The molecular diversity of adaptive convergence. Science. 2012; 335(6067):457–61.View ArticlePubMedGoogle Scholar
  3. Lang GI, Rice DP, Hickman MJ, Sodergren E, Weinstock GM, Botstein D, Desai MM. Pervasive genetic hitchhiking and clonal interference in forty evolving yeast populations. Nature. 2013; 500(7464):571–4.View ArticlePubMed CentralPubMedGoogle Scholar
  4. Gerrish P, Lenski R. The fate of competing beneficial mutations in an asexual population. Genetica. 1998; 102:127–44.View ArticlePubMedGoogle Scholar
  5. Schiffels S, Szollosi GJ, Mustonen V, Lassig M. Emergent neutrality in adaptive asexual evolution. Genetics. 2011; 189(4):1361–75.View ArticlePubMed CentralPubMedGoogle Scholar
  6. Good BH, Rouzine IM, Balick DJ, Hallatschek O, Desai MM. Distribution of fixed beneficial mutations and the rate of adaptation in asexual populations. Proc Nat Acad Sci. 2012; 109:4950–5.View ArticlePubMed CentralPubMedGoogle Scholar
  7. Weinreich DM, Watson RA, Chao L. Perspective: Sign epistasis and genetic constraint on evolutionary trajectories. Evolution. 2005; 59(6):1165–74.View ArticlePubMedGoogle Scholar
  8. Weinreich DM, Delaney NF, DePristo MA, Hartl DL. Darwinian evolution can follow only very few mutational paths to fitter proteins. Science. 2006; 312(5770):111–4.View ArticlePubMedGoogle Scholar
  9. Poelwijk FJ, Kiviet DJ, Weinreich DM, Tans SJ. Empirical fitness landscapes reveal accessible evolutionary paths. Nature. 2007; 445:383–6.View ArticlePubMedGoogle Scholar
  10. Kvitek DJ, Sherlock G. Reciprocal sign epistasis between frequently experimentally evolved adaptive mutations causes a rugged fitness landscape. PLoS Genet. 2011; 7(4):1002056.View ArticleGoogle Scholar
  11. Silva RF, Mendonça SC, Carvalho LM, Reis AM, Gordo I, Trindade S, et al.Pervasive sign epistasis between conjugative plasmids and drug-resistance chromosomal mutations. PLoS Genet. 2011; 7(7):1002181.View ArticleGoogle Scholar
  12. Dawid A, Kiviet DJ, Kogenaru M, de Vos M, Tans SJ. Multiple peaks and reciprocal sign epistasis in an empirically determined genotype-phenotype landscape. Chaos: Interdisciplinary J Nonlinear Sci. 2010; 20(2):026105.View ArticleGoogle Scholar
  13. Salverda ML, Dellus E, Gorter FA, Debets AJ, Van Der Oost J, Hoekstra RF, et al.Initial mutations direct alternative pathways of protein evolution. PLoS Genet. 2011; 7(3):1001321.View ArticleGoogle Scholar
  14. Woods RJ, Barrick JE, Cooper TF, Shrestha U, Kauth MR, Lenski RE. Second-order selection for evolvability in a large escherichia coli population. Science. 2011; 331(6023):1433–6.View ArticlePubMed CentralPubMedGoogle Scholar
  15. Weinreich DM, Chao L. Rapid evolutionary escape by large populations from local fitness peaks is likely in nature. Evolution. 2005; 59:1175–82.View ArticlePubMedGoogle Scholar
  16. Iwasa Y, Michor F, Nowak MA. Stochastic tunnels in evolutionary dynamics. Genetics. 2004; 166(3):1571–9.View ArticlePubMed CentralPubMedGoogle Scholar
  17. Weissman DB, Desai MM, Fisher DS, Feldman MW. The rate at which asexual populations cross fitness valleys. Theor Population Biol. 2009; 75(4):286–300.View ArticleGoogle Scholar
  18. Gokhale CS, Iwasa Y, Nowak MA, Traulsen A. The pace of evolution across fitness valleys. J Theor Biol. 2009; 259(3):613–20.View ArticlePubMed CentralPubMedGoogle Scholar
  19. Altland A, Fischer A, Krug J, Szendro IG. Rare events in population genetics: stochastic tunneling in a two-locus model with recombination. Phys Rev Lett. 2011; 106(8):088101.View ArticlePubMedGoogle Scholar
  20. Weissman DB, Feldman MW, Fisher DS. The rate of fitness-valley crossing in sexual populations. Genetics. 2010; 186(4):1389–410.View ArticlePubMed CentralPubMedGoogle Scholar
  21. Rozen DE, Habets MG, Handel A, de Visser JAG. Heterogeneous adaptive trajectories of small populations on complex fitness landscapes. PLoS One. 2008; 3(3):1715.View ArticleGoogle Scholar
  22. Handel A, Rozen DE. The impact of population size on the evolution of asexual microbes on smooth versus rugged fitness landscapes. BMC Evolutionary Biol. 2009; 9(1):236.View ArticleGoogle Scholar
  23. Jain K, Krug J, Park S-C. Evolutionary advantage of small populations on complex fitness landscapes. Evolution. 2011; 65(7):1945–55.View ArticlePubMedGoogle Scholar
  24. Van Nimwegen E, Crutchfield JP. Metastable evolutionary dynamics: crossing fitness barriers or escaping via neutral paths?. Bull Math Biol. 2000; 62(5):799–848.View ArticlePubMedGoogle Scholar
  25. Desai MM, Fisher DS. Beneficial mutation–selection balance and the effect of linkage on positive selection. Genetics. 2007; 176(3):1759–8.View ArticlePubMed CentralPubMedGoogle Scholar
  26. Szendro IG, Franke J, de Visser JAG, Krug J. Predictability of evolution depends nonmonotonically on population size. Proc Nat Acad Sci. 2013; 110(2):571–6.View ArticlePubMed CentralPubMedGoogle Scholar


© Ochs and Desai; licensee BioMed Central. 2015

This is an Open Access article distributed under the terms of the Creative Commons Attribution License (, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly credited. The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated.