Fixation probability of rare nonmutator and evolution of mutation rates

Abstract Although mutations drive the evolutionary process, the rates at which the mutations occur are themselves subject to evolutionary forces. Our purpose here is to understand the role of selection and random genetic drift in the evolution of mutation rates, and we address this question in asexual populations at mutation‐selection equilibrium neglecting selective sweeps. Using a multitype branching process, we calculate the fixation probability of a rare nonmutator in a large asexual population of mutators and find that a nonmutator is more likely to fix when the deleterious mutation rate of the mutator population is high. Compensatory mutations in the mutator population are found to decrease the fixation probability of a nonmutator when the selection coefficient is large. But, surprisingly, the fixation probability changes nonmonotonically with increasing compensatory mutation rate when the selection is mild. Using these results for the fixation probability and a drift‐barrier argument, we find a novel relationship between the mutation rates and the population size. We also discuss the time to fix the nonmutator in an adapted population of asexual mutators, and compare our results with experiments.


Introduction
Because most mutations are deleterious, the mutation rate can not be too high; in fact, in an infinitely large population, for a broad class of fitness functions, an error threshold has been shown to exist above which the deleterious effects of mutation cannot be compensated by selection (Eigen 1971;Jain and Krug 2007). The mutation rate is not zero either (Baer et al. 2007), and it has been argued that the stochastic fluctuations in a finite population limit the evolution of mutation rates below a certain level since in small enough populations, the advantage gained by lowering the mutation rate cannot compensate the effect of random genetic drift (Lynch 2010). Empirical data for organisms with widely different effective population size show a negative correlation between the deleterious mutation rate and the population size (Sung et al. 2012), and some quantitative insight into this relationship has been obtained by treating all deleterious mutations to be lethal (Lynch 2011). However, this is clearly an extreme scenario, and it is important to ask how the deleterious mutation rate evolves when mutations are only weakly deleterious.
Many theoretical and experimental investigations have also shown that in an adapting asexual population, a mutator allele causing a higher mutation rate than that of the nonmutator can get fixed [see a recent review by Raynes and Sniegowski (2014)]. As the mutators produce not only deleterious but also beneficial mutations at a higher rate than the nonmutators, the mutator allele can hitchhike to fixation with favorable mutations (Smith and Haigh 1974;Taddei et al. 1997). However, once the population has reached a high fitness level, high mutation rates are detrimental because most mutations will now be deleterious, and in such a situation, the mutation rate is expected to decrease (Liberman and Feldman 1986). Indeed, in some experiments (Tr€ obner and Piechocki 1984;Notley-McRobb et al. 2002;McDonald et al. 2012;Turrientes et al. 2013;Wielgoss et al. 2013), the mutation rate of an adapted population carrying a mutator allele has been seen to decrease and the time to fixation has been measured, but a theoretical understanding of this time scale is missing.
To address the issues discussed above, we study the fate of a rare nonmutator in a large asexual population of mutators using a multitype branching process (Patwa and Wahl 2008). An important difference between the previous works on mutator hitchhiking (Taddei et al. 1997;Andre and Godelle 2006;Wylie et al. 2009;Desai and Fisher 2011) and our study is that here the mutator population is assumed to be at mutation-selection equilibrium and is therefore not under positive selection. However, compensatory mutations that alleviate the effect of deleterious mutations are included in our model. We find that when only deleterious mutations are present, a nonmutator can get fixed with a probability that increases with the deleterious mutation rate of the mutator. Compensatory mutations in the mutator population are expected to decrease the fixation probability of the nonmutator, and we find that this intuition is indeed correct when deleterious mutations in the mutator are effectively lethal. But, surprisingly, when the deleterious mutations are mildly harmful, the fixation probability is found to initially increase and then decrease as the rate of compensatory mutations increases. Our study thus identifies the conditions under which the spread of nonmutators is suppressed in the absence of positive selection, and complements earlier works in which a mutator hitchhikes with beneficial mutations to fixation (Taddei et al. 1997;Andre and Godelle 2006;Wylie et al. 2009;Desai and Fisher 2011).
Using our results for the fixation probability and a drift-barrier argument which states that the advantage offered by a decrease in the deleterious mutation rate is limited by random genetic drift in a finite population (Lynch 2010), we find that the deleterious mutation rate decreases with increasing population size in accordance with experimental data (Sung et al. 2012). However, unlike previous theoretical work that treats the deleterious mutations to be effectively lethal (Lynch 2011), here we consider both strongly and weakly deleterious mutations, and not only reproduce the result in Lynch (2011), but also find a new scaling law in the latter case. We also use the results for the fixation probability to find the time to lower the mutation rate in an adapted population of mutators and compare our theoretical results with recent experiments (McDonald et al. 2012;Wielgoss et al. 2013).

Model and Methods
We consider an asexual population in which the fitness of an individual with k deleterious mutations is given by WðkÞ ¼ ð1 À sÞ k , where the selection coefficient 0 < s < 1. A deleterious mutation is allowed to occur at a rate U d and a beneficial one at a rate U b \U d . We are interested in the fate of a nonmutator that arises in this population and whose total mutation rate is smaller than that of the mutator. In a sufficiently large population of mutators in which stochastic fluctuations due to genetic drift may be ignored, this can be addressed using a branching process (Patwa and Wahl 2008), as described below.
The fixation probability p(k, t) of a single copy of a nonmutator allele with fitness W(k) present at generation t changes according to (Johnson and Barton 2002) 1 À pðk; tÞ ¼ exp À WðkÞ WðtÞ where WðtÞ ¼ P 1 k¼0 WðkÞpðk; tÞ is the average fitness of the mutator population and p(k, t) is the mutator frequency. The above equation expresses the fact that a single copy of the rare allele in the fitness class k whose offspring distribution is Poisson with mean W(k)/W(t) will be lost eventually if each of its offspring, which may undergo mutations with probability M(k?k 0 ), do not survive. Here we consider strong mutators whose mutation rate is much higher than that of the nonmutator (Sniegowski et al. 1997;Oliver et al. 2000) and therefore neglect the mutation rate of the latter in most of the following discussion (however, see Fig. 1). We also assume that the mutator population is at mutation-  Figure 1. Dependence of the fixation probability obtained using a multitype branching process on the deleterious mutation rate U d for two values of the selection coefficient s and compensatory mutation rate U b ¼ 0. The points are obtained by numerically solving (2) when the mutation rate of the nonmutator is zero (s,□), and the stationary state solution of (1) when the nonmutator's mutation rate is 50 times lower than that of the mutator (+, 9). The lines show the analytical result (6). selection equilibrium as is likely to be the case in large populations that have been evolving for a long time in a constant environment. As a result, the probability p(k, t) becomes time-independent. These considerations lead to a relatively simpler, but still highly nonlinear equation given by The above expression, of course, reduces to the wellknown single locus equation (Fisher 1922;Haldane 1927) when the nonmutator can be present in only one genetic background, but here we are dealing with a multitype branching process because a nonmutator can arise in any fitness class.
The total fixation probability is obtained on summing over all genetic backgrounds (Johnson and Barton 2002), where the probability that a nonmutator arises in a background of k deleterious mutations is given by the mutator frequency p(k) in the stationary state.
Although the steady-state frequency p(k) in the absence of compensatory mutations that mitigate the effect of deleterious mutations is known exactly (Kimura and Maruyama 1966;Haigh 1978), the corresponding solution with nonzero U b is not known. We therefore compute the mutator frequency numerically for nonzero U b using (A1) given in Appendix 1, and use these results in (2) to find the fixation probability for arbitrary U b . To make analytical progress, we use a perturbation theory in which the effect of the small dimensionless parameter U b =s can be studied by expanding the quantities of interest in a power series in U b =s, and write The terms p 0 ðkÞ and p 0 ðkÞ corresponding to n = 0 in the above expansion give the results in the absence of compensatory mutations, and in Appendix 1, we calculate the stationary state fraction p(k) to linear order in U b =s.

Fixation probability
In the absence of compensatory mutations We first consider the case when U b ¼ 0. Taking the logarithm on both sides of (2), and expanding the left hand side (LHS) up to p 2 0 ðkÞ, we find that either p 0 ðkÞ ¼ 0, or where the average fitness W 0 ¼ e ÀU d and the average number of deleterious mutations k 0 ¼ U d =s (Kimura and Maruyama 1966;Haigh 1978). The last expression on the right hand side (RHS) of (5) is obtained by expanding the exponentials as the parameters U d and s are small. As the fixation probability must not be negative, the expression (5) is valid when k\b k 0 c, and the solution p 0 ðkÞ ¼ 0 holds otherwise. Here ⌊x⌋ denotes the largest integer less than or equal to x. More generally, a nonmutator can get fixed if its fitness WðkÞ % e Àsk is larger than the average fitness e Às k of the mutator population, or k\b kc, k being the average number of deleterious mutations (Johnson and Barton 2002). Equation (5) shows that the fixation probability p 0 ðkÞ decreases as the number of deleterious mutations increase, as one would intuitively expect. However, the probability p 0 ðkÞ that a nonmutator would arise in a background with k\ k 0 deleterious mutations increases. On summing over the backgrounds in which a nonmutator can arise, as explained in Appendix 2, we find that the total fixation probability falls in two distinct regimes defined by whether U d is below or above s: For k 0 ( 1, as a mutation is costly, it can be treated as effectively lethal (Johnson 1999). In this situation, the advantage conferred by the nonmutator is simply given by 1 À e ÀU d % U d and the classical result for the single locus problem gives the fixation probability to be 2U d (Fisher 1922;Haldane 1927). For k 0 ) 1, the total fixation probability apparently receives contribution from k 0 genetic backgrounds, but merely ffiffiffiffi ffi k 0 p genetic backgrounds are actually relevant because the Poisson-distributed frequency p 0 ðkÞ has a substantial weight for fitness classes that lie within a width ffiffiffiffi ffi k 0 p of the mean (also, see Appendix 2). Equation (6) shows that for fixed s, the nonmutator is more likely to be fixed when U d is large. But, for a given U d , the fixation probability initially increases with the selection coefficient and then saturates to 2U d . In Figure 1, the analytical results above are compared with those obtained by numerically iterating (2) and (1) when the mutation rate of the nonmutator is zero and U d =50, respectively, and we see a good agreement in both cases.

Including compensatory mutations
We now study how compensatory mutations in the mutator population affect the fixation probability of the nonmutator. Figure 2 shows that when k 0 ( 1, the fixation probability decreases with U b , but for k 0 ) 1, it changes nonmonotonically: it first increases and then decreases with increasing U b . To understand this behavior, consider the change dp tot ¼ p tot À p 0 in the fixation probability due to compensatory mutations which is simply given by dp tot ¼ X b kc k¼0 p 0 ðkÞdpðkÞ þ p 0 ðkÞdpðkÞ þ dpðkÞdpðkÞ: (7) When U b is nonzero, the change in the fixation probability dpðkÞ ¼ pðkÞ À p 0 ðkÞ and the mutator frequency dpðkÞ ¼ pðkÞ À p 0 ðkÞ behave in a qualitatively different manner. With increasing U b , the average fitness of the mutator population increases which, by virtue of (2), decreases the fixation probability of the nonmutator, i.e., dp(k) < 0. However, as the frequency of individuals with less deleterious mutations increases when U b is nonzero, the change in the mutator fraction dp(k) > 0. Thus, the change in the total fixation probability given by (7) receives both positive and negative contributions, and it is not obvious which one of these factors would have a larger effect.
To address this question, we calculate the fixation probability for small U b =s as described below. Substituting (4) in the expression (7) for dp tot , and neglecting terms of order ðU b =sÞ 2 and higher, we find that dp tot % ðU b =sÞp 1 , where The contribution p 1 ðkÞ is calculated in Appendix 3, and we find that p 1 ðkÞ % À2s k 0 ð1 À p 0 ðkÞÞ; k\b k 0 c; which is negative, as expected. An expression for the fraction p 1 ðkÞ is obtained in Appendix 1, and its behavior is shown in Figure 3 for small and large k 0 . For small k 0 , the frequency p 0 ðkÞ is close to one in the zeroth fitness class and zero elsewhere. But the correction p 1 ðkÞ is negligible in all the fitness classes. For large k 0 , the contribution p 1 ðkÞ is significantly different from zero in many fitness classes and can be approximated by Thus, as claimed above, the fraction p 1 ðkÞ is positive for k\ k 0 and negative for k [ k 0 (also, see Fig. 3).
When U d ( s, as already mentioned, the fraction p 1 ðkÞ is negligible in all the fitness classes and p 0 ð0Þ % 1. Using these results in (8) and (9), we get p 1 ¼ À2s k 0 , and thus dp tot This reduction in the fixation probability of the nonmutator when U b is nonzero is expected as the effect of compensatory mutation is to restore the mutators that have suffered lethal mutation to the zeroth mutation class, thus enabling them to offer competition to the nonmutators. When U d ) s, as shown in Appendix 3, we can obtain a quantitative estimate of the initial increase in dp tot by calculating the sum on the RHS of (8) to obtain (A14), and thence  (11) and (12). The broken curve for U b =s [ 0:1 is a linear fit, 0:1 À 0:24U b =s, to the numerical data. For U d =s ¼ 0:1, the ratio U b =s is also below 0.1 as U b is assumed to be smaller than U d .
Thus, we find that for small U b , the increase of the mutator frequency in fitness classes with fewer deleterious mutations dominates the increase in the mutator fitness resulting in positive dp tot . However, for large U b , the net change in the fixation probability is negative because the last term in the summand of (7), which is also negative, enters the picture. As the maximum in dp tot occurs at large U b =s, the perturbation theory described here can not capture the eventual decrease in this parameter regime. A quantitative comparison of the results obtained by numerically solving (2) and (A1) for arbitrary U b with the analytical results (11) and (12) for small U b =s is shown in Figure 2, and we observe a good match when U b =s is small. For large U b =s and U d =s, a fit to the numerical data shows that the fixation probability decreases linearly with U b .

Evolution of mutation rates in finite populations
The drift-barrier hypothesis states that in a finite population, the beneficial effect of lower deleterious mutation rate can be outweighed by the stochastic effects of random genetic drift which limits the evolution of mutation rates (Lynch 2010). In a finite population of size N, a mutation that decreases the deleterious mutation rate confers an indirect selective advantage and will spread through the population. However, as U d decreases, the fixation probability of such a mutant decreases until it reaches its neutral value p neu ¼ 1=N. Here we have calculated the fixation probability p 0 neglecting stochastic fluctuations. The full fixation probability Π that includes the neutral and the large population limit may be obtained as follows. The fixation time for a mutator in a finite population of nonmutators when all mutations are deleterious has been calculated using a diffusion theory by Jain and Nagar (2013), and shown to increase exponentially with the population size. The fixation probability $ e À2NS is thus exponentially small in the population size (Kimura 1980;Assaf and Mobilia 2011), where we have identified the rate of decrease of fixation probability with a selection coefficient 2S. This effective selection coefficient is found to match exactly with the result (6) for the fixation probability p 0 obtained here using a branching process. Although this is not a rigorous proof, these observations strongly suggest that the fixation probability of a nonmutator in a finite population of size N is of the classical form (Kimura 1962) where S ¼ p 0 =2. We also mention that the probability 2S depends on the difference in the deleterious mutation rate of the mutator and the nonmutator when the mutation rate of the nonmutator is nonzero (Jain and Nagar 2013), and has also been shown to be insensitive to the distribution of selective effects (Desai and Fisher 2011). Thus, according to (13), a crossover between positive selection and neutral regime occurs when p 0 $ N À1 and gives a lower bound on the mutation rates. We recall that the fixation probability p 0 in (6) shows a transition when U d $ s, and at this mutation rate, the fixation probability p 0 $ s. This translates into a change in the behavior of U d when Ns crosses one, and we have Thus, in the weak selection regime (Ns ( 1), the deleterious mutation rate depends on the selection coefficient and decreases faster than when the selection is strong. Figure 4 shows the preliminary results of our numerical simulations for a finite size population of mutators with mutation rate U d in which nonmutators with mutation rate U d =2 can arise with a certain probability. This population of nonmutators and mutators evolves via standard Wright-Fisher dynamics, and the time to fix the nonmutators is measured (Jain and Nagar 2013). For a fixed N, the fixation time is found to increase as the mutation rate of the mutator is decreased until a minimum mutation rate is reached below which the fixation time remains constant. This lower bound, shown in Figure 4, exhibits different scaling behavior in the weak and strong selection regimes, in accordance with (14).

Fixation probability
A rare mutator arising in a population of nonmutators carries a higher load of deleterious mutations but offers indirect benefit by producing more beneficial mutations. The fixation probability of a rare mutator in a finite nonmutator population has been studied by Andre and Godelle (2006) and Wylie et al. (2009) analytically, and found to vary nonmonotonically with the mutation rate of the mutator. It has been shown that the fixation probability is of the classical form (13) where the effective selection coefficient S when scaled by the selective advantage s increases (decreases) when the ratio of mutation rate to selection coefficient is below (above) one. Here, we studied a situation in which a nonmutator appears in a mutator population and is beneficial as it produces fewer deleterious mutations, and calculated its fixation probability p tot using a branching process. The mutator population is assumed to be at mutation-selection balance, and therefore, by definition, selective sweeps resulting in the spread of favorable mutations are neglected. However, it is interesting to note that the scaled fixation probability of the nonmutator obtained here also changes its behavior when the deleterious mutation rate is of the order of the selection coefficient, see (6). Our work significantly extends the previous result of Lynch (2011) as the deleterious effect of mutations is allowed to be mild here, and therefore, we are dealing with a truly multilocus problem.
Compensatory mutations that alleviate the effect of deleterious mutations are found to have a surprising effect on the fixation probability of the nonmutator. Although they improve the fitness of the mutator population, it also means that the nonmutator can arise in a better genetic background where it has a better chance of fixation. Thus, compensatory mutations affect both the resident mutator population and the invading nonmutator allele in a positive manner. The effect of these two factors on the fixation probability of the nonmutator is, however, opposite and can result in an unexpected increase in the fixation probability of the nonmutator when compensatory mutations are present. Here we have shown analytically that this scenario is realized when the mutations are weakly deleterious and the compensatory mutation rate is small, as illustrated in Figure 2. The increase in the fixation probability due to compensatory mutations can be quite high, but we do not have analytical estimates for this. An exact solution of (A1) would, of course, pave the way for a better analytical understanding but is currently not available.

Fixation time
In a maladapted asexual population, the mutators can sweep the population as they facilitate rapid adaptation (Raynes and Sniegowski 2014). But as the population adapts and the supply of beneficial mutations diminishes, mutators have a detrimental effect on the population fitness and a mutation that lowers the mutation rate is favored. In bacteria Escherichia coli, several genes (such as mut T and mut Y) are involved in avoiding or repairing the errors that occur during the replication process, and defects in these genes can lead to the mutator phenotype (Miller 1996). But compensatory mutations in the defective error-repair machinery can reduce the mutation rate, at least, partially (Wielgoss et al. 2013). We therefore model this situation by assigning a probability b with which mutators can convert into nonmutators due to a mutation in the proofreading or error-repair region. In E. coli, the conversion probability f from nonmutator to mutators has been estimated to be $ 10 À6 per bacterium per generation (Boe et al. 2000). But the probability b for the reverse mutation is not known, although one expects b < f, possibly because it is a gain-of-function mutation (Wielgoss et al. 2013).
When the rate Nb at which the nonmutators are produced from the mutators is small enough that the new alleles behave independently, the time taken to fix the nonmutator population is given by T ¼ ðNbp tot Þ À1 . In a long-term evolution experiment on E. coli, Wielgoss et al. (2013) found the mutation rate to decrease by about a factor two in a nearly adapted mutator population with a mutation rate 150 times that of the wild type in two lineages. As the population size in Lenski's experiments has been estimated to be about 10 7 (Wahl et al. 2002), the product Nb can be at most ten which is not too large. We first note that in the experiment of Wielgoss et al. (2013), the fixation time was longer in the lineage in which the mutation rate decreased by a smaller amount, in accordance with (6). To make a quantitative comparison, we consider the ratio of the times for the two lineages, as T depends strongly on the probability b which is not known experimentally. Using the data in Table 2 of Wielgoss et al. (2013), we find the ratio of fixation time in mutT mutY-L background to that in mutT mutY-E background to be 9209/5157%1.8. The theoretical formula (6), on replacing U d by the difference between the mutation rate of the nonmutator and mutator, yields 1.5 (1.2) when mutations are assumed to be strongly (weakly) deleterious and the selection coefficient same in both lineages. As (6) is obtained assuming that the mutators are strong whereas the mutation rates decreased merely by a factor two in the experiment, a more careful examination is needed. Solving (1) numerically in the stationary state, we find that the ratio is unaffected when the mutations are strongly deleterious. But using the mutation rates in Table 2 of Wielgoss et al. (2013) and s $ 0.01 yield the ratio to be about 4.5. Although the theoretical conclusions (1.5 À 4.5) are in reasonable agreement with experiments, the above analysis suggests that the reversion probability b may not be too small (i.e., Nb [ $ 1), and a more sophisticated theory that takes care of the interference between the nonmutators (Gerrish and Lenski 1998) may be required to obtain a closer match. We close this discussion by noting that in an experiment on Saccharomyces cerevisiae in which the adapted population reduced its genomewide mutation rate by almost a factor four in two of the experimental lines (McDonald et al. 2012), the fixation time seems to increase with the mutation rate, in contradiction with the experiment of Wielgoss et al. (2013) and the theory presented here.

Evolution of mutation rates
Experiments show that the mutation rate decays as N À0:7 for prokaryotes and N À0:9 for eukaryotes (Sung et al. 2012). The population size and deleterious mutation rates are negatively correlated as deleterious mutations can get fixed in small populations due to stochastic fluctuations, but not in large populations where the genetic drift is ineffective (Lynch 2010). Here, we have shown that a reciprocal relationship between the population size and mutation rate holds for large populations, but for small populations, the deleterious mutation rate decreases much faster, see Figure 4. This is in contrast to experimental results mentioned above where the data has been fitted assuming a single scaling law. In view of our theoretical results discussed above, a more careful analysis of experimental data is required.
While the evolution of deleterious mutation rate has received much attention, to the best of our knowledge, analogous theoretical predictions for the beneficial mutation rate are not available. As large populations experience clonal interference (Gerrish and Lenski 1998) which results in the wastage of beneficial mutations, the rate of beneficial mutations is observed to be smaller in large populations in microbial experiments (Perfeito et al. 2007). An understanding of the relationship between the population size and the rate of beneficial mutations would be an interesting avenue to explore. Other potential factors that can affect the correlation between the mutation rate and the population size include epistasis and recombination. Here, we have also ignored the cost of fidelity, and it remains to be seen how the results presented here are affected on including it (Kimura 1967;Kondrashov 1995;Dawson 1998). A more detailed understanding of the mutation rates, both empirically and theoretically, remains a goal for the future.