Evolutionary advantage of small populations on complex fitness landscapes
aa r X i v : . [ q - b i o . P E ] F e b Evolutionary advantage of small populations on complexfitness landscapes
Kavita Jain , Joachim Krug and Su-Chan Park ∗ , Theoretical Sciences Unit and Evolutionary and Organismal Biology Unit, Jawaharlal Nehru Centre for AdvancedScientific Research, Jakkur P.O., Bangalore 560064, India Institut f¨ur Theoretische Physik, Universit¨at zu K¨oln, Z¨ulpicherstr. 77, 50937 K¨oln, Germany Department of Physics, The Catholic University of Korea, Bucheon 420-743, Korea
Email: Kavita Jain - [email protected]; Joachim Krug - [email protected]; Su-Chan Park ∗ - [email protected]; ∗ Corresponding author
Running Title : Advantage of small populations
Contact Information (for all authors)
Kavita Jain postal address : Theoretical Sciences Unit and Evolutionary and Organismal Biology Unit,Jawaharlal Nehru Centre for Advanced Scientific Research, Jakkur P.O., Bangalore 560064, India work telephone number : +91-80-22082948
E-mail :[email protected] Krug postal address : Institute f¨ur Theoretische Physik, Universit¨at zu K¨oln, Z¨ulpicher Str. 77, 50937K¨oln, Germany work telephone number : +49-221-470-2818
E-mail :[email protected] 1u-Chan Park (Corresponding Author) address : Department of Physics, The Catholic University of Korea, 43-1 Yeokgok 2-dong,Wonmi-gu, Bucheon 420-743, Republic of Korea work telephone number : +82-2-2164-4524
E-mail : [email protected] 2 bstract
Recent experimental and theoretical studies have shown that small asexual populations evolving oncomplex fitness landscapes may achieve a higher fitness than large ones due to the increasedheterogeneity of adaptive trajectories. Here we introduce a class of haploid three-locus fitnesslandscapes that allow the investigation of this scenario in a precise and quantitative way. Our mainresult derived analytically shows how the probability of choosing the path of the largest initial fitnessincrease grows with the population size. This makes large populations more likely to get trapped atlocal fitness peaks and implies an advantage of small populations at intermediate time scales. Therange of population sizes where this effect is operative coincides with the onset of clonal interference.Additional studies using ensembles of random fitness landscapes show that the results achieved for aparticular choice of three-locus landscape parameters are robust and also persist as the number ofloci increases. Our study indicates that an advantage for small populations is likely whenever thefitness landscape contains local maxima. The advantage appears at intermediate time scales, whichare long enough for trapping at local fitness maxima to have occurred but too short for peak escapeby the creation of multiple mutants.
KEY WORDS: clonal interference, finite population, fitness landscape, fixation probability,three-locus models 3onsider a population experiencing a recent environmental change. Assuming that the population isill-adapted to the new environment, as is typically the case in the beginning of an evolutionexperiment (Lenski and Travisano 1994), the adaptation to the new environment relies on the supplyof beneficial mutations available to the population. During the early stages of adaptation, it is agood approximation to assume that the supply of beneficial mutations is unlimited. Then as the rateat which the beneficial mutations appear is given by the mutation rate per individual times thepopulation size, a large population is expected to experience more beneficial mutations and henceadapt faster than a small population.This conclusion does not seem to depend on the topography and epistatic interactions in thefitness landscape. If the fitness landscape is non-epistatic in the sense that (beneficial) mutations actmultiplicatively on fitness, a large population adapts faster than a small population, although thisadvantage is strongly reduced in asexuals due to clonal interference (Gerrish and Lenski 1998;de Visser et al. 1999; Wilke 2004; Park and Krug 2007; Park et al. 2010). Even if the fitnesslandscape is highly epistatic such as a maximally rugged ( house-of-cards ) landscape (Kingman1978), a large population still wins on average (Park and Krug 2008). This leaves the question as towhether small populations might be at an advantage in complex fitness landscapes with anintermediate degree of epistasis and ruggedness.In a recent experimental work (Rozen et al. 2008) studying the adaptation dynamics ofpopulations of
E. coli in simple and complex nutrient environments, it was found that smallpopulations could attain higher fitness than large populations in a complex medium which can beexpected, on general grounds, to give rise to a rugged and epistatic fitness landscape. The observedfitness advantage of small populations was associated with a greater heterogeneity in their adaptivetrajectories compared to large populations. Specifically, small populations that eventually reachedthe highest fitness levels were often the ones that initially displayed a rather shallow fitness increase,whereas the fitness of those that gained a large initial advantage tended to level off quickly. In largepopulations trajectories were more uniform and typically showed a rapid initial fitness increasefollowed by a significant slowing down or saturation.These experiments suggest that an adaptive advantage can arise from the higher level ofstochasticity in the incorporation of beneficial mutations displayed by small populations, providedthe topography of the underlying fitness landscape is sufficiently complex. As the detailed structure4f the experimental fitness landscape is unknown and unfeasible to determine, it is useful toinvestigate this mechanism theoretically. In previous work this was done using extensive simulationsfor a class of random fitness landscapes with tunable ruggedness (Handel and Rozen 2009). Themain conclusion from this study was that an advantage of small populations can be observed in asubstantial fraction of random landscapes, but the dependence of the effect on parameters such asthe population size, the mutation rate and the fitness effect of beneficial mutations was not exploredsystematically.In this article, we introduce a minimal, analytically tractable model which captures the dynamicbehavior of the population fitness in the experiments by Rozen et al. (2008). We show that thesimplest fitness landscape that can exhibit a small population advantage is a haploid, diallelicthree-locus landscape where the genotypes of minimal and maximal fitness are separated by threemutational steps. There are then 3! = 6 distinct shortest paths leading from the global fitnessminimum (the wild type) to the global maximum, corresponding to the different orderings in whichthe mutations are introduced into the population (Gokhale et al. 2009). We distinguish between smooth paths along which fitness increases monotonically, and rugged paths containing at least onedeleterious mutation. The three-locus landscape is constructed in such a way that the rugged pathscontain a local fitness maximum, and they confer the greatest initial fitness increase to a populationinitially fixed for the minimal fitness genotype.The existence of rugged paths is the hallmark of sign epistasis , a specific type of geneticinteraction under which a given mutation can be beneficial or deleterious depending on the geneticbackground (Weinreich et al. 2005; Poelwijk et al. 2007). Sign epistasis implies that at least some ofthe mutational pathways leading to the global maximum of the fitness landscape are rugged, andthus inaccessible to an adaptive dynamics that is constrained to increase fitness in each step. Inaddition, sign epistasis is a necessary (but not sufficient) condition for the existence of multiplefitness maxima (Poelwijk et al. 2011). The presence of sign epistasis was established in several recentexperimental studies, where all combinations of a selected set of individually beneficial or deleteriousmutations were constructed and their fitness effects (or some proxy thereof) were measured(Weinreich et al. 2006; Lozovsky et al. 2009; de Visser et al. 2009; Carneiro and Hartl 2010). Thereis also considerable evidence for the existence of multiple fitness maxima from evolution experimentsusing bacterial and viral populations (Korona et al. 1994; Burch and Chao 1999, 2000;5lena and Lenski 2003). It is therefore reasonable to assume that sign epistasis and mulitple fitnessmaxima are present also in the complex environment considered by Rozen et al. (2008).The analysis presented below shows that the speed of adaptation is generally an increasingfunction of population size both along smooth and along rugged paths. However, the probability withwhich a particular type of path is chosen depends on population size in such a way that smallpopulations can be favored at least over a certain range of time scales. In particular, as we shallshow, the probability to choose the rugged path in the three-locus model rises sharply with the onsetof clonal interference, and it approaches unity when the dynamics becomes completely deterministicfor very large populations, because then the mutation with the largest initial fitness increase iscertain to fix (Jain and Krug 2007b). The dynamics of small populations are less predictable, andthey therefore enjoy an advantage by more frequently avoiding getting trapped at the local fitnessmaximum. The main part of the article is devoted to a detailed, quantitative analysis of thisscenario. We then show that the mechanism identified within the specific three-locus model is robustby simulating populations on variants of the house-of-cards fitness landscape with three or more loci,and conclude the paper with a discussion of our key results.
Models
FITNESS LANDSCAPES
In the main part of this work we consider the space of genotypes composed of three loci with twoalleles each, which will be denoted 0 and 1. Each genotype is assigned a fitness W according to W (000) = 1 W (001) = W (010) = 1 + s W (100) = 1 + s W (011) = W (101) = W (110) = (1 + s ) W (111) = (1 + s ) (1)where s > s > s ) < s < (1 + s ) . Thus there is a local fitness maximum atgenotype { } and the global maximum is located at genotype { } .6n addition to the landscape equation (1), we consider two ensembles of random fitnesslandscapes consisting of L -locus genotypes with two alleles at each site. In the first ensemble referredto as unconstrained ensemble, the least fit genotype is assigned the allele 0 at every locus and fitness1 while the rest of the genotypes are given fitnesses W (genotype) = 1 + Sx, (2)where S controls the strength of selection and x is a random number drawn from an exponentialdistribution with mean 1. This is Kingman’s house-of-cards model adapted to a finite number ofdiallelic loci (Kingman 1978; Kauffman and Levin 1987; Jain and Krug 2007b; Park and Krug 2008).It is also instructive to study a constrained version of the above model in which the fittest genotypehas all loci with allele 1 (Kl¨ozer 2008; Carneiro and Hartl 2010). Such landscapes are generated byassigning the maximum value amongst 2 L − POPULATION DYNAMICS
We mainly work with a finite population of size N which evolves according to standardWright-Fisher dynamics in discrete generations. In each generation, an offspring chooses a parentwith a probability proportional to the parent’s fitness and copies the parent’s genotype. Then thepoint mutation process is implemented symmetrically in which 0 ↔ µ . Thisprocess is repeated until all N offspring have been generated. In the actual simulations, we treatedthe population size of each genotype as a random variable which is sampled according to amultinomial distribution; for details see Park and Krug (2007).It is also useful to compare the results obtained for finite N with the predictions of thedeterministic mutation-selection dynamics of “quasispecies” type which applies for infinitepopulations (B¨urger 2000; Jain and Krug 2007a). This is done by iterating the deterministicevolution equations for the frequencies f ( σ, t ) of genotype σ at generation t , which read f ( σ, t + 1) = P σ ′ M ( σ ← σ ′ ) W ( σ ′ ) f ( σ ′ , t ) P σ ′ W ( σ ′ ) f ( σ ′ , t ) . (3)7ere M ( σ ← σ ′ ) = µ d ( σ,σ ′ ) (1 − µ ) L − d ( σ,σ ′ ) is the probability that an L -locus genotype σ ′ mutates togenotype σ at Hamming distance d ( σ, σ ′ ). Results
EVOLUTION TIME SCALES
We begin by a comparison of the time taken to reach the global maximum along the smooth and therugged paths on the fitness landscape defined by equation (1). A population that is initially fixed forthe minimum fitness genotype { } has a choice between three single site mutations. Two of these(to genotypes { } and { } ) lead the population to a smooth path towards the global fitnessmaximum { } , whereas the third leads to the local fitness maximum { } from which thepopulation can escape only by the creation of a double mutant (Iwasa et al. 2004;Weinreich and Chao 2005; Weissman et al. 2009).We estimate the time scales of the relevant evolutionary processes in the strong selection weakmutation (SSWM) regime (Gillespie 1983, 1984; Orr 2002), where N µ ≪ { } . Starting from the wild type, each of the single stepmutants is generated in the population at rate N µ and goes to fixation with probability π ( s ) givenby (Kimura 1962) π ( s ) ≈ − e − s − e − Ns (4)with s = s or s . For N − ≪ s ≪
1, the fixation probability π ( s ) ≈ s . The waiting times for lowfitness ( T ) and high fitness ( T ) mutants that will ultimately fix are therefore T ≈ µN s , T ≈ µN s . (5)Adaptation along one of the smooth paths proceeds by sequentially fixing two additional beneficialmutations with selection coefficient s , and the total evolution time is therefore T smooth ≈ T .By contrast, populations that choose the rugged path need to escape from the local fitness peak { } in order to reach the global maximum. The corresponding escape time can be estimated alongthe lines of Weinreich and Chao (2005). Following these authors we introduce the selection8oefficients s ben = W (111) W (100) − ≈ s − s , s del = W (100) W (101) − ≈ s − s (6)which express the relative fitness advantage of the global maximum compared to the local peak( s ben ) and that of the valley genotypes compared to the local peak ( s del ), respectively. Depending onthe population size, the peak escape can proceed through two distinct pathways. In populationssmaller than a critical size N c (Weinreich and Chao 2005), the two mutations separating thegenotypes { } and { } fix sequentially, while in larger populations they fix simultaneously. Weare interested in population size ≫ N c ≈ ln( s/µ ) /s which can be easily satisfied in the SSWMregime as N s ≫
1. In the simultaneous mode the escape time is given approximately by T esc ≈ s del N µ s ben . (7)Assuming all selection coefficients s , s , s ben , s del to be of a similar magnitude s , we see that T esc T , ∼ sµ ≫ µ ≪ s , which is expected to hold under most conditions. In particular, it is true in theSSWM regime because N µ ≪ N s ≫
1. Equation (8) implies that the evolution time T rugged along a rugged path is dominated by the escape time T esc , and is much larger than T smooth . However,both equation (5) and equation (7) share the same dependence on population size N , so once thetype of evolutionary path is chosen, a large population is always at a relative advantage. MEAN FITNESS EVOLUTION AND HETEROGENEOUS ADAPTIVETRAJECTORIES
Figure 1 shows the evolution of the population fitness obtained from simulations of theWright-Fisher model in the landscape defined by equation (1). Each curve contains data averagedover many stochastic histories for a given value of N , keeping other parameters fixed, and startingwith all individuals at the genotype { } with lowest fitness. At short times the fitness rises morerapidly for larger populations, as expected on the basis of the estimates given in equation (5) for N µ ≪
1, while for
N µ > N . Large populations are also seen to beat an advantage for extremely long times, beyond 10 generations. However, for both parameter sets9isplayed in the figure, the ordering of fitness with increasing population size is reversed in anintermediate time interval, which begins at around 2000 generations.The origin of this reversal is illustrated in Figure 2, which shows individual fitness trajectoriesfor the parameter set of Figure 1 (b) and two different population sizes. Individual fitnesstrajectories display a step-like behavior, which reflect the transitions in the most populated genotype .In particular, populations in which the local peak genotype { } becomes dominant are seen toremain trapped at the local peak for a long time (compare to equation (7)). Although the initial risein fitness is much faster for the large populations than for the small ones, the fraction of trajectoriesthat take the rugged path (and thus get trapped at the local peak) also increases with increasing N ,from 11/20 in the Figure 2(a) ( N = 10 ) to 19/20 in Figure 2(b) ( N = 10 ). As a consequence, thefitness after 10 generations, when averaged over all trajectories, is larger for the small populationsthan for the large ones. Similar to the experimental observations of Rozen et al. (2008) and thesimulations of Handel and Rozen (2009), smaller populations reach a higher fitness level becausetheir adaptive trajectories are more heterogeneous, allowing them to avoid trapping at the localfitness peak in a larger number of trials. PATH PROBABILITY
To quantify the statement that large populations are more likely to take the rugged path, weintroduce the probability P r ( N ) that the rugged path is taken by a population of size N . Insimulations, the probability P r was measured by counting the number of events in which { } becomes the most populated genotype for the first time. In Figure 3 these numerical results arecompared with the analytical expressions (discussed below) and the two are seen to be in very goodagreement. We see that P r generally increases with N (provided µ is not too small) thus supportingour main contention. Before presenting an analytic calculation of P r covering the full range ofpopulation sizes, we discuss the limiting cases of very small and very large populations. SSWM regime:
N µ ≪ N µ ≪
1, the path probability P r is equal to the probability that the first mutant that will befixed in a population initially monomorphic for the genotype { } is the local peak genotype { } .10hat is, P r | Nµ ≪ = π ( s ) π ( s ) + 2 π ( s ) . (9)When selection is weak, in the sense of N s , ≪
1, the fixation probability is given by its neutralvalue π = 1 /N and we obtain P r = 1 /
3. On the other hand, for strong selection (and assuming that s , ≪
1) we have P SSW Mr ≈ s s + 2 s (10)independent of N , which is equal to 0.54 and 0.56, respectively, in the two cases displayed in Figure1. Deterministic quasispecies regime: N → ∞ For very large populations the local peak mutant is always present in the population in considerablenumbers and can therefore be expected to dominate the population with a probability approachingunity (Jain and Krug 2007b). The quantitative analysis of this case is based on the deterministicinfinite population dynamics defined by equation (3). As the initial population is assumed to befixed at the genotype { } , for small mutation rates equation (3) gives f ( σ, t = 1) ∼ µ d ( σ, whichis the same for all genotypes at constant Hamming distance from { } . For t >
1, the genotypicpopulation can be determined by a simple construction described by Jain and Krug (2005) and Jain(2007) in which the population of a genotype increases exponentially with its fitness, starting from f ( σ, f (111 , < f (100 ,
1) but the fitnessvalues W (111) > W (100), it is possible that the genotype { } becomes the most populated onebefore { } thus bypassing the local maximum. The population at sequence { } overtakes that of { } at time t when f (000 , t ) = f (100 , t ), which on using f ( σ, t ) ∝ µ d ( σ, W ( σ ) t gives t = − ln µ/ ln W (100). Similarly the time t at which the population of the global maximumovertakes that of the initial sequence is given by t = − µ/ ln W (111). Thus the condition forbypassing corresponding to t < t reads W (100) < W (111) ⇔ s < s (11)which is ruled out by construction. Thus bypassing cannot occur and we conclude thatlim N →∞ P r ( N ) = 1 . (12)11 lonal interference regime The phenomenon of interest in this paper occurs in the intermediate range of population sizes where P r increases from its small population value in equation (9) to the large population limit in equation(12). This regime is more difficult to analyze because of the presence of multiple competing mutantclones in the population. To find an analytic expression for P r ( N ) in this regime, we first reduce thethree-locus problem into a single locus with three alleles, say A , B , and C with respective fitness 1,1 + s , and 1 + s . The two genotypes { } and { } are lumped into a single allele B . Themutation from A to B occurs with probability 2 µ and that from A to C with µ . No other mutationis possible, which ensures that either B or C will be eventually fixed. It is clear that the fixationprobability of allele B approximates 1 − P r .We now present an approximate calculation of P r for the three-allele model using ideas fromclonal interference theory. At time zero the population is monomorphic for allele A . We would liketo determine the probability that an allele B which originates at some time t > t + t f . It is plausible to assume that the fixation and origination of a mutation are notcorrelated, and hence to treat the two processes separately.Let us first consider the probability p ( t ) for the allele B to originate at time t . An allele B appears in the population at rate 2 N µ and would, in the absence of other mutations, go to fixationwith probability π ( s ). As is customary in the field (Maynard Smith 1971; Gerrish and Lenski 1998;Desai and Fisher 2007), we interpret the fixation probability π as the probability for the mutantpopulation to survive genetic drift and, thus, to reach a size large enough for the further timeevolution to be essentially deterministic. Mutations of type B which reach this level are calledcontending mutations (Gerrish 2001), and they arise at rate 2 N µπ ( s ). To obtain p ( t ) this rate hasto be multiplied with the probability that the contending mutation in question is the first to appearamong all possible contenders for fixation. To estimate the probability that no contending mutation(of any type) has appeared before time t , we use a Poisson approximation in which the probabilityfor the non-occurrence of an event is the negative exponential of the expected number of events.Since the expected number of contending mutations arising up to time t is N µ ( π ( s ) + 2 π ( s )) t , weconclude that p ( t ) = 2 N µπ ( s ) exp( − N µ ( π ( s ) + 2 π ( s )) t ) . (13)12e next determine the probability p that the fixation of the contending mutation of type B isnot impeded by the appearance of a contending mutation of type C at some time larger than t .Such a mutation can only arise from the wildtype population A . According to our assumptions, theevolution of the frequency x of B alleles after time t follows the deterministic, logistic growthequation dxdt = s x (1 − x ) . (14)The expected number of C alleles that arise from the wildtype population until the fixation time t f is therefore Z t f N µ (1 − x ( t )) dt = Z x N µxs dx = − N µs ln( x ) , (15)where x is the initial frequency of the contending B allele. As before, we use a Poissonapproximation to determine the probability that no contending mutation of type C arises until thefixation of the B allele. Since the expected number of such contending mutations is − N µ ln( x ) π ( s ) /s , we have p = exp[ N µ ln( x ) π ( s ) /s ] . (16)To obtain the total probability 1 − P r for the B allele to fix we multiply p by p and integrate overthe initial time t , which gives1 − P r ( N ) = p Z ∞ dt p ( t ) = 2 π ( s ) π ( s ) + 2 π ( s ) e Nµ ln( x ) π ( s ) /s . (17)To complete the analysis we have to determine x . A naive argument would suggest that x = 1 /N because the contending mutation should start from a single mutant. However, this doesnot include the fact that the fixation process conditioned on survival is faster than the logisticgrowth with x = 1 /N (Hermisson and Pennings 2005; Desai and Fisher 2007). A simpleapproximate way to take into account this effect is to let the contending mutant clone start atfrequency x = 1 / ( s N ) ≫ /N (Maynard Smith 1971). Inserting this into equation (17) we obtainthe final result given by P r ( N ) = 1 − π ( s ) π ( s ) + 2 π ( s ) exp( − Q ( N )) , (18)where Q ( N ) = N µ ln(
N s ) π ( s ) /s . (19)13igure 3 shows that equation (18) agrees well with simulations of the three-allele model described atthe beginning of this section. In the strong selection limit, using π ( s ) ≈ s the above expression for P r can be simplified to give P r ( N ) ≈ − (1 − P SSW Mr ) e − Q ( N ) with Q ( N ) ≈ N µ ln(
N s ) s /s . Thusthe path probability is close to the SSWM value when Q ( N ) ≪ Q ( N ). The increase of P r beyond the SSWM value which ultimately gives rise to the reversal in theordering of fitness with increasing population size in Figure 1, takes place when Q ( N ) is of orderunity, N µ ln(
N s ) ∼ O (1) (20)where it is assumed that both selection coefficients have a similar scale s .It is straightforward to generalize the above derivation of P r to an L -locus system where L − s and one mutation confers a higher advantage s > s . This merely increases the mutation rate from allele A to allele B to ( L − µ and leads tothe expression P r ( N ) = 1 − ( L − π ( s ) π ( s ) + ( L − π ( s ) exp( − N µ ln(
N s ) π ( s ) /s ) . (21) HOUSE-OF-CARDS MODELS
At this point the question naturally arises as to how generic our results are with respect to thenumber of loci and the structure of the fitness landscape. To address this question, we simulatedpopulations evolving in two ensembles of random fitness landscapes, the unconstrained andconstrained house-of-cards models. In these models the fitness values of genotypes differing by singlemutational steps are assumed to be uncorrelated. While this is not likely to be the case in realfitness landscapes (Miller et al. 2011), the house-of-cards models constitute the conceptually simplestrealization of a generic, rugged fitness landscape which is essentially parameter-free, apart from theoverall fitness scale S in equation (2).Before discussing the dynamics of adaptation in random landscapes, we may ask how typical thethree-locus landscape equation (1) itself is within the constrained ensemble with L = 3. The maintopographic features of the landscape equation (1) are (i) the existence of a single local maximum, inaddition to the global maximum, and (ii) the existence of 2 rugged and 4 smooth paths from theglobal minimum { } to the global maximum { } . The enumeration of all 6! = 720 possibilities14f ordering the fitness values of the 6 genotypes intermediate between the global minimum and theglobal maximum shows that these two features are shared by a fraction of 11 / ≈ . P ( t, N, N ′ ) of asmall population advantage, defined as the probability that the mean fitness of a population of size N at generation t is larger than that of a population of size N ′ at the same generation, where N < N ′ .To determine P ( t, N, N ′ ) numerically, we first calculated the mean fitness for population size N by averaging over 128 independent evolutionary histories on a single landscape, and then comparedit to that for a different population size N ′ on the same landscape. Finally the outcome of thecomparison was averaged over 10 samples of the random landscape ensemble. In order to avoidspurious contributions from cases where the fitness is in fact independent of population size and anapparent advantage of small populations arises due to fluctuations, we count only instances in whichthe inequality (1 − α ) w ( N ) > w ( N ′ ) is satisfied, where w ( N ) stands for the calculated mean fitnessfor population size N . We choose α = 10 − , which is large enough to remove the error due tofluctuations. Of course, this may also screen out landscapes with very small advantage, but we donot think that this effect is substantial, since a mean fitness advantage of less than 0.1% is negligiblecompared to the scale S = 0 . generations and is most pronounced for population size N ≈ . Note that inFig. 1 a where the selection coefficient is of the order of 0.1 as in Fig. 5, the fitness advantage of a15mall population becomes conspicuous when the population with size 10 is compared to that with10 in the time window between 10 and 10 . That is, the quantitative as well as the qualitativebehavior of the house-of-cards model in the constrained ensemble is similar to the three-locus modelin the previous section. Thus, we may conclude that the landscape equation (1) is quite genericwithin the constrained ensemble, in the sense that between 20% and 30% of all landscapes showsimilar dynamical features. In the unconstrained ensemble the probabilities are reduced by about afactor of 2, but the overall behavior is the same. MORE THAN THREE LOCI
Within the house-of-cards models, it is straightforward to investigate the effect of varying thenumber of loci L . In the simulations reported in this subsection we allow for at most one mutationper individual and generation, and specify the genome-wide mutation probability U rather than themutation probability µ per locus. The general relation between U and µ is (Park and Krug 2008) U = 1 − (1 − µ ) L ≈ µL when µL ≪
1. With increasing L the difference between the constrained andunconstrained ensembles becomes less important, because populations typically do not reach theglobal maximum within the simulation time, and hence its precise location is irrelevant. In thefollowing we therefore restrict the discussion to the standard (unconstrained) model.As an illustration, we present simulation results for L = 20 in Figure 6. In this case theexpected fitness value of the global maximum is ≈ .
44 which is far larger than the mean fitness ofthe population with size 10 within the observation time. This means that most of evolutions up tosize 10 and generation 10 have not arrived at the global maximum, but rather explore thesurroundings of the low fitness starting genotype. As in the three-locus case, the mean fitnessaveraged over all landscapes is monotonic in the population size (Figure 6a), and the plot of P ( t, N, N ′ ) in Figure 6b displays a similar, though more pronounced advantage of small populationsfor N = 10 − . However, for N ≥ a new peak is seen to appear in P ( t, N, N ′ ) at later time,which is not present for L = 3.We may interpret this behavior in terms of the different evolutionary regimes described byJain and Krug (2007b). The disappearance of the first peak in P ( t, N, N ′ ) marks the point wherethe mutation supply rate N U is sufficiently large for the population to easily escape from localfitness maxima by the creation of double mutants. Such a population will however still have16ifficulties to cross wider fitness valleys, and hence it will tend to get trapped at local maxima whichare separated by three or more mutational steps. The mechanism for an advantage of smallpopulations that was found in the three-locus model thus reemerges on a larger scale, giving rise tosecondary peaks in P ( t, N, N ′ ). Discussion
Understanding the effect of population size on the rate of adaptation is a central problem inevolutionary theory, which continues to attract considerable attention (Weinreich and Chao 2005;Gokhale et al. 2009; Weissman et al. 2009; Lynch and Abegg 2010). Motivated by the experimentsof Rozen et al. (2008), the present work has addressed a specific aspect of this general problem. Inthe experiments, several populations of
E. coli consisting of either 5 × or 2 . × individualswere evolved in a complex nutrient medium which can be modeled by a complex fitness landscape.The fitness measurements after 500 generations showed that small populations can achieve higherfitness than large populations.A classic scenario in which a small population can acquire an evolutionary advantage because ofgenetic drift has been put forward in the framework of Wright’s shifting balance theory (Wright1931), referred to as SBT in the following discussion. Apart from the intrinsic shortcomings of theSBT (Coyne et al. 1997), however, there are several reasons why it cannot be directly applied to theexperimental situation of Rozen et al. (2008). First of all, the population in the experiments is notstructured, and it is therefore not possible for different demes to occupy separate fitness peaks asassumed in phase II of the SBT. Second, the number of generations in the experiment is too smallfor the entire population to cross a fitness valley either by the fixation of deleterious mutations(phase I of the SBT) or by the simultaneous fixation of individually deleterious but jointly beneficialmutations (Weinreich and Chao 2005). Hence for a proper explanation of the experimentalobservations it cannot be assumed that the population resides at a local fitness optimum from thebeginning of the process. Rather, the evolutionary trajectories begin in a fitness valley, and thedynamics is determined by the competition between different fitness peaks that are available to thepopulation (Rozen et al. 2008).In the preceding sections we have demonstrated that under these conditions, small populations17ay indeed reach higher fitness levels than large ones because they are more likely to evade trappingat local fitness maxima. Our detailed study of a three-locus model where a single local fitness peakcompetes with the global maximum has shown that the dynamical behavior of the population fitnessis not determined by the time scale to acquire beneficial mutation(s) alone and depends on the pathprobability P r ( N ) also. This probability increases from a constant value 1 / P SSW Mr and finallyapproaches the deterministic value unity as the population size N increases. For the parameterregimes in which P r is constant in N , as the waiting times T rugged and T smooth both decrease with N ,larger populations are at an advantage.However when the probability to take a rugged path increases with N , a larger population mayget trapped at the local fitness maximum thus acquiring lower fitness than a small population. Forthe parameters in Figure 1, the path probability exceeds the SSWM value 0 .
55 for
N > . N = 10 but close to unity for N ≥ and therefore apopulation of size N ∼ is able to acquire a higher fitness than larger populations.A key result of our analysis is that the regime of population sizes in which this mechanismoperates coincides with the onset of clonal interference, which occurs precisely when the criterion inequation (20) is satisfied (Gerrish and Lenski 1998; de Visser et al. 1999; Wilke 2004; Park and Krug2007). In the context of the three-allele model discussed previously, clonal interference implies that ahigh fitness clone (allele C) may arise while the low fitness mutant (allele B) is still on the way tofixation, thus enhancing the probability for C to fix and increasing the probability for the populationto evolve along a rugged path. From the criterion in equation (20), we can determine the mutationprobability for the bacterial populations used in the experiment of Rozen et al. (2008). Since thesmall population advantage is observed around N = 5 × and the characteristic size of selectioncoefficients derived from the fitness trajectories in the experiment is ≃ .
1, the estimated mutationprobability is µ ≈ N ln( N s ) ≃ − , (22)which should be interpreted as a beneficial mutation rate. This estimate is consistent with values forthe beneficial mutation rate in E.coli obtained by other experimental approaches, which range from10 − to 10 − (Hegreness et al. 2006; Perfeito et al. 2007).Our theoretical analysis also shows that the advantage of small populations over large ones istransient and at sufficiently large times, the fitness of the large population exceeds that of the small18opulation. This reversal occurs when the large population escapes the fitness valley at time T esc ∼ ( N µ ) − (see Eq. 7) and should be testable over an experimentally accessible time scale of 10 generations if the population size exceeds (10 µ ) − ∼ where we have used (22). In theexperiments of Rozen et al. (2008), the large populations had a size of N ∼ which is two ordersof magnitude below our prediction and hence the reversal in the ordering of fitness was not observed.It would be interesting to test this effect in experiments using larger populations.In their work, Rozen et al. (2008) attributed the fitness advantage seen for small populations tothe heterogeneity in evolutionary trajectories. This qualitative description can be made precise byconsidering predictability of the first adaptive step along the evolutionary trajectory. Following Orr(2005) and Roy (2009) we define the predictability P to be the probability that two replicatepopulations follow the same evolutionary trajectory, which equals the sum of the squares of theprobability of these trajectories. Specializing to the first adaptive step in our three-locus system, thepossible outcomes { } , { } and { } occur with probability P r , (1 − P r ) / − P r ) / P = P r + (1 − P r ) which increases from P = 1 / P r itself. For the parameter range N µ ≪ , N s ≫ P r is independent of N , the predictability P < ACKNOWLEDGEMENTS
This work was supported by DFG within SFB 680
Molecular Basis of Evolutionary Innovations . Wethank Arjan de Visser, Siegfried Roth and Sijmen Schoustra for helpful discussions and comments,and two anonymous reviewers for their suggestions on the manuscript. K.J. and J.K. acknowledgethe hospitality of KITP, Santa Barbara, and support under NSF grant PHY05-51164 during theinitial stages of this work.
LITERATURE CITED
B¨urger, R. 2000. The Mathematical Theory of Selection, Recombination, and Mutation. Wiley,Chichester.Burch, C. L., and L. Chao. 1999. Evolution by small steps and rugged landscapes in the RNA virus φ
6. Genetics 151:921–927.Burch, C. L., and L. Chao. 2000. Evolvability of an RNA virus is determined by its mutationalneighborhood. Nature 406:625–628.Carneiro, M., and D. L. Hartl. 2010. Adaptive landscapes and protein evolution. Proc. Nat. Acad.Sci. USA 107 (suppl. 1):1747–1751.Coyne, J. A., N. H. Barton, and M. Turelli. 1997. A critique of Sewall Wright’s shifting balancetheory of evolution. Evolution 51:643–671.Desai, M. M., and D. S. Fisher. 2007. Beneficial mutation-selection balance and the effect of linkageon positive selection. Genetics 176:1759–1798.Elena, S. F., and R. E. Lenski. 2003. Evolution experiments with microorganisms: The dynamicsand genetic bases of adaptation. Nature Reviews Genetics 4:457–469.20errish, P. J. 2001. The rhythm of microbial adaptation. Nature 413:299–302.Gerrish, P. J., and R. E. Lenski. 1998. The fate of competing beneficial mutations in an asexualpopulation. Genetica 102/103:127–144.Gillespie, J. H. 1983. Some properties of finite populations experiencing strong selection and weakmutation. Am. Nat. 121:691–708.———. 1984. Molecular evolution over the mutational landscape. Evolution 38:1116–1129.Gokhale, C. S., Y. Iwasa, M. A. Nowak, and A. Traulsen. 2009. The pace of evolution across fitnessvalleys. J. Theor. Biol. 259:613–620.Handel, A., and D. E. Rozen. 2009. The impact of population size on the evolution of asexualmicrobes on smooth versus rugged fitness landscapes. BMC Evolutionary Biology 9:236.Hegreness, M., N. Shoresh, D. Hartl, and R. Kishony. 2006. An equivalence principle for theincorporation of favorable mutations in asexual populations. Science 311:1615–1617.Hermisson, J., and P. S. Pennings. 2005. Soft sweeps: molecular population genetics of adaptationfrom standing genetic variation. Genetics 169:2335–2352.Iwasa, Y., F. Michor, and M. A. Nowak. 2004. Stochastic tunnels in evolutionary dynamics.Genetics 166:1571–1579.Jain, K. 2007. Evolutionary dynamics of the most populated genotype on rugged fitness landscapes.Phys. Rev. E 76:031922.Jain, K., and J. Krug. 2005. Evolutionary trajectories in rugged fitness landscapes. J. Stat. Mech.:Theory Exp. P04008.———. 2007a. Adaptation in simple and complex fitness landscapes. in U. Bastolla, M. Porto, H. E.Roman, and M. Vendruscolo, eds. Structural approaches to sequence evolution: Molecules,networks and populations. Springer, Berlin.———. 2007b. Deterministic and stochastic regimes of asexual evolution on rugged fitnesslandscapes. Genetics 175:1275–1288. 21auffman, S., and S. Levin. 1987. Towards a general theory of adaptive walks on rugged landscapes.J. Theor. Biol. 128:11–45.Kimura, M. 1962. On the probability of fixation of mutant genes in a population. Genetics47:713–719.Kingman, J. F. C. 1978. A simple model for the balance between selection and mutation. J. Appl.Prob. 15:1–12.Kl¨ozer, A. 2008. NK fitness landscapes. Master’s thesis, University of Cologne, Institute forTheoretical Physics.Korona, R., C. H. Nakatsu, L. J. Forney, and R. E. Lenski. 1994. Evidence for multiple adaptivepeaks from populations of bacteria evolving in a structured habitat. Proc. Nat. Acad. Sci. USA91:9037–9041.Lenski, R. E., and M. Travisano. 1994. Dynamics of adaptation and diversification: a10,000-generation experiment with bacterial populations. Proc. Nat. Acad. Sci. USA 91:6808–6814.Lynch, M., and A. Abegg. 2010. The rate of establishment of complex adaptations. Mol. Biol. Evol.27:1404–1414.Lozovsky, E. R., T. Chookajorn, K. M. Brown, M. Imwong, P. J. Shaw, S. Kamchonwongpaisan,D. E. Neafsey, D. M. Weinreich, and D. L. Hartl. 2009. Stepwise acquisition of pyrimethamineresistance in the malaria parasite. Proc. Nat. Acad. Sci. USA 106:12025–12030.Maynard Smith, J. 1971. What use is sex? J. theor. Biol. 30:319–335.Miller, C. R., P. Joyce, and H.A. Wichman. 2011. Mutational effects and population dynamicsduring viral adaptation challenge current models. Genetics 187:185–202,.Orr, H. A. 2002. The population genetics of adaptation: The adaptation of DNA sequences.Evolution 56:1317–1330.———. 2005. The probability of parallel evolution. Evolution 59:216–220.Park, S.-C., and J. Krug. 2007. Clonal interference in large populations. Proc. Nat. Acad. Sci. USA104:18135–18140. 22——. 2008. Evolution in random fitness landscapes: the infinite sites model. J. Stat. Mech.:Theory Exp. P04014.Park, S.-C., D. Simon, and J. Krug. 2010. The speed of evolution in large asexual populations. J.Stat. Phys. 138:381–410.Perfeito, L., L. Fernandes, C. Mota, and I. Gordo. 2007. Adaptive mutations in bacteria: high rateand small effects. Science 317:813–815.Poelwijk, F. J., D. J. Kiviet, D. M. Weinreich, and S. J. Tans. 2007. Empirical fitness landscapesreveal accessible evolutionary paths. Nature 445:383–386.Poelwijk, F. J., S. T˘anase-Nicola, D. J. Kiviet, S. J. Tans. 2011. Reciprocal sign epistasis is anecessary condition for multi-peaked fitness landscapes. J. Theor. Biol. 272:141–144.Roy, S. 2009. Probing evolutionary repeatability: Neutral and double changes and the predictabilityof evolutionary adaptation. PLoS ONE 4:e4500.Rozen, D. E., M. G. J. L. Habets, A. Handel, and J. A. G. M. de Visser. 2008. Heterogeneousadaptive trajectories of small populations on complex fitness landscapes. PLoS ONE 3:e1715.de Visser, J. A. G. M., S.-C. Park, and J. Krug. 2009. Exploring the effect of sex on empirical fitnesslandscapes. Am. Nat. 174:S15–S30.de Visser, J. A. G. M., C. W. Zeyl, P. J. Gerrish, J. L. Blanchard, and R. E. Lenski. 1999.Diminishing returns from mutation supply rate in asexual populations. Science 283:404–406.Weinreich, D. M., and L. Chao. 2005. Rapid evolutionary escape by large populations from localfitness peaks is likely in nature. Evolution 59:1175–1182.Weinreich, D. M., N. F. Delaney, M. A. DePristo, and D. L. Hartl. 2006. Darwinian evolution canfollow only very few mutational paths to fitter proteins. Science 312:111–114.Weinreich, D. M., R. A. Watson, and L. Chao. 2005. Sign epistasis and genetic constraint onevolutionary trajectories. Evolution 59:1165–1174.Weissman, D. B., M. M. Desai, D. S. Fisher, and M. W. Feldman. 2009. The rate at which asexualpopulations cross fitness valleys. Theor. Pop. Biol. 75:286–300.23ilke, C. O. 2004. The speed of adaptation in large asexual populations. Genetics 167:2045–2053.Wright, S. 1931. Evolution in Mendelian populations. Genetics 16:97–159.
FIGURES a b PSfrag replacements generationgeneration m e a nfi t n e ss m e a nfi t n e ss W (111) W (111) Figure 1: Average fitness as a function of time on the fitness landscape defined by equation (1) for µ = 10 − , (a) s = 0 . s = 0 .
25 and (b) s = 0 . s = 0 .
05, and various population sizes indicatedin the figures. As a guide to the eye, the maximum fitness values W (111) for each case are also drawn.The data have been averaged over 10 histories. Fitness increases with population size at short timesand at long times, but in both cases this relationship is reversed for a range of population sizes atintermediate times. 24 b PSfrag replacements generationgeneration m e a nfi t n e ss m e a nfi t n e ss W (111) Figure 2: Population mean fitness as a function of time for 20 histories on the fitness landscapeequation (1). The population size is (a) N = 10 and (b) 10 , the selection coefficients are s = 0 . s = 0 .
05 and the mutation probability is µ = 10 − (same as in Figure 1(b)). The smooth curvesdepict the average over 10 independent runs. a b PSfrag replacements s = 0 . s = 0 . s = 0 . s = 0 . NN P r P r µ = 10 − µ = 10 − µ = 10 − µ = 10 − µ = 10 − µ = 10 − Figure 3: Fixation probability P r of the allele C for the simplified three alleles system obtained usingnumerical simulations for µ = 10 − (triangles), 10 − (circles), and 10 − (squares). Dotted lines showthe analytic prediction of equations (18,19). 25 b [ ¯ w ] PSfrag replacements [ ¯ w ] PSfrag replacements generation generation Figure 4: Mean fitness evolution in random three-locus fitness landscapes. Fitness has been averagedover 10 realizations of the landscape and 128 population histories in each realization, with parameters µ = 10 − and S = 0 .
1, and population sizes as indicated in the figures. Both for (a) the constrainedensemble and (b) the unconstrained ensemble the mean fitness increases monotonically with populationsize for all times. a b P ( t , N , N ) PSfrag replacements P ( t , N , N ) PSfrag replacements generation generation Figure 5: Probability of small population advantage in random three-locus fitness landscapes. Plotsshow P ( t, N, N ′ ) as a function of t , with N ′ = 10 N and N as indicated in the figures, µ = 10 − ,and S = 0 .
1. Part (a) shows the constrained ensemble, part (b) the unconstrained ensembles. Theconstrained ensemble with the maximum fitness at the antipodal point is more likely to allow smallpopulations to have larger mean fitness. 26 b [ ¯ w ] PSfrag replacements P ( t , N , N ) PSfrag replacements generation generation Figure 6: (a) Mean fitness evolution and (b) probability of small population advantage for the uncon-strained ensemble with L = 20 loci. Parameters are U = 10 − and S = 0 ..