Abstract
The accumulation of deleterious mutations due to the process known as Muller's ratchet can lead to the degeneration of nonrecombining populations. We present an analytical approximation for the rate at which this process is expected to occur in a haploid population. The approximation is based on a diffusion equation and is valid when N exp(−u/s) ⪢ 1, where N is the population size, u is the rate at which deleterious mutations occur, and s is the effect of each mutation on fitness. Simulation results are presented to show that the approximation estimates the rate of the process better than previous approximations for values of mutation rates and selection coefficients that are compatible with the biological data. Under certain conditions, the ratchet can turn at a biologically significant rate when the deterministic equilibrium number of individuals free of mutations is substantially >100. The relevance of this process for the degeneration of Y or neo-Y chromosomes is discussed.
NEW mutations arise continuously within populations, the vast majority probably being slightly deleterious (Crow 1993). An effectively infinite asexual population subject to recurrent deleterious mutations will achieve an equilibrium resulting from the continuous appearance of new mutations opposed by selection against them: deterministic mutation-selection balance. With independent and identical effects of each mutation, the equilibrium number of mutations in a haploid randomly mating population follows a Poisson distribution with mean u/s (Kimura and Maruyama 1966; Haigh 1978), where u is the per genome mutation rate and s is the selection coefficient against a deleterious mutation. At this deterministic equilibrium, the number of individuals free of mutations (the best or “least-loaded” class) in a large population of N breeding adults, is n0 = N exp(−u/s).
But random genetic drift plays a role in a finite population, and it may perturb this equilibrium, leading to the loss of the best class. In the absence of recombination, and with the reasonable assumption that back mutation is negligible for strongly selected mutations, the loss of this class is irreversible and mutations will continually accumulate in the population, leading to the decline of its mean fitness. This is the process known as Muller's ratchet (Muller 1964; Felsenstein 1974). Once the best class is lost, the new least-loaded class is now the one that has one mutation, but this is also subject to stochastic loss, so that a repetition of successive losses of the least-loaded class can be seen as successive clicks of the ratchet.
One of the important questions in this process concerns the rate or speed at which it operates, or: how much time does it take for the ratchet to click one notch? Because the degeneration of the Y chromosome (Charlesworth 1978, 1996; Rice 1994; Charlesworth and Charlesworth 1997, 1998) and the fate of asexual populations (Pamiloet al. 1987; Gabrielet al. 1993; Lynch et al. 1993, 1995) may involve the operation of Muller's ratchet, the quantification of its rate is of great biological importance.
Haigh (1978) suggested that the most important parameter for the ratchet mechanism would be n0, because it is the loss of the best class that drives the process; the smaller the n0, the faster the ratchet is likely to be. Although Haigh suggested an expression [Equation 9a in Haigh (1978)] for the time between clicks of the ratchet, this is basically a fit to his simulation results. Bell (1988) also suggested an expression based on a fit of the time to extinction of the best class as a function of n0 (roughly 10 n0). Several attempts toward the quantification of the process have subsequently been made, either using a quantitative genetics approach for estimating the rate of change of the average number of mutations (Pamiloet al. 1987; Gabrielet al. 1993; Lynchet al. 1993; Higgs and Woodcock 1995; Prügel-Bennett 1997) or using diffusion theory to calculate the mean time to loss of the least-loaded class (Stephanet al. 1993; Charlesworth and Charlesworth 1997).
In this article, it is shown that the size of the best class is not sufficient to predict the speed of the ratchet. For the same value of n0, the ratchet can turn at very different speeds. An approximation based on a diffusion equation, for the case n0 > 1, which almost certainly applies to the evolution of the Y chromosome (Charlesworth 1996), is presented, together with simulations of asexual haploid populations to check its validity. Gessler (1995) has derived an approximation that seems to work well for the case n0 < 1. Comparisons between the simulation results and the predictions from the analytical expression suggest that the formulas seem to predict the rate of the process better than the previous formulas for a region of parameter space that is of biological interest (Charlesworth and Charlesworth 1997).
APPROXIMATION BASED ON THE DIFFUSION EQUATION
We start with a haploid asexual population at equilibrium under mutation-selection balance, with x0 = exp(−u/s) being the frequency of individuals in the least-loaded class. The existence of this equilibrium requires that the number of individuals in the best class (n0 = Nx0) be >1. When n0 < 1, this equilibrium may not be approached in a finite population, because it requires the existence of individuals that have a very low probability of actually being present (Gessler 1995; Gessler and Xu 1999).
The way the frequency of the best class varies through time is dictated by mutation taking it below the equilibrium value, selection restoring it to that value, and by stochastic fluctuations due to drift. We wish to quantify how much time it takes for the frequency of the best class (from now on denoted by x), with initial value x0, to reach the value zero. One way to do this is to use a diffusion equation for the density function of the time until absorption occurs, subject to the condition that x = 0 is the only absorbing state (Ewens 1979, Equations 4.39, 4.40, p. 123). To solve this equation, one has to evaluate the deterministic change (drift coefficient) and stochastic variance (diffusion coefficient) in x.
Assuming a Wright-Fisher population, the diffusion coefficient is just the variance due to binomial sampling of N individuals from the previous generation (Stephanet al. 1993; Charlesworth and Charlesworth 1997),
The drift coefficient, representing the expected change in x due to mutation and selection (Stephanet al. 1993; Charlesworth and Charlesworth 1997), is
Let us now make the simplifying assumption that the changes in mean fitness are sufficiently small that they can be approximated by small perturbations from the equilibrium value (Stephanet al. 1993). This implies that the system is close to its equilibrium state most of the time, as is supported by our simulations (see Figure 9). We assume that the perturbations in w are mostly due to small fluctuations in the least-loaded class. This assumption was employed by Stephan et al. (1993) and Charlesworth and Charlesworth (1997), and is justified in practice by the observation that, for large N, the distribution among the classes that are present at any time remains close to the Poisson distribution given by the deterministic equilibrium formula [see Table 1 and Figure 4 of Charlesworth and Charlesworth (1997)]. We may express the mean fitness close to equilibrium as a Taylor expansion in x/x0:
Taking just the linear term in (2a), we obtain an approximation for the reduction in mean fitness below its equilibrium value as
When the frequency x is above its equilibrium value, the system responds with a reduction of mean fitness
In fact, the approach to the neighborhood of equilibrium takes some time (Haigh 1978), and just after a click the new best class is above its equilibrium value, so that
We now may write the drift coefficient as
Using these drift and diffusion coefficients, the time spent in the frequency interval [0, x0] (Ewens 1979, Equation 4.39) is
Using expressions 3a and 3b, and evaluating the integrals numerically for a given population size, mutation rate, and selection coefficient, we obtain the expected time to loss of the least-loaded class as T(N, u, s) = T0,x0 + Tx0,1.
SIMULATION METHODS
For a given population size (N), genomic mutation rate to deleterious mutations (u), and selection coefficient against each mutation (s), haploid asexual populations were simulated starting at mutation-selection equilibrium (Kimura and Maruyama 1966); i.e., the number of individuals in the class with m mutations is
Assuming that the sequence of events is mutation, reproduction, and selection, populations were then run for 100 generations. After this initial period, populations were run for >2000 generations and up to 100,000 generations for conditions under which the ratchet clicks slowly, so that the average time between clicks of the ratchet could be measured. Every generation, the number of mutations in every individual is counted and the number of individuals with the least number of mutations (least-loaded class) is recorded. If, at a given generation, the number of mutations in the least-loaded class increases, the ratchet has clicked. To form a new generation, individuals are sampled randomly from the previous generation, then subjected to the occurrence of mutations sampled from a Poisson distribution with mean u, and assigned probabilities of survival as (1 − s)k, where k is the number of mutations that an individual carries. A new generation of N individuals is constructed by comparing the probability of survival of each individual with a pseudorandom number drawn from a uniform distribution in the interval [0, 1]. Each run was repeated several times; generally five replicates were performed to obtain the results presented in the next section.
Average time between clicks of the ratchet, ±2 SE, for constant s as a function of the population size (as N increases, u increases so that n0 is maintained constant; s is set equal to 0.015). n0 is 202 for circles and 20.2 for triangles. The points joined with a line are the estimated times given by the approximation presented in the text.
Although this simulation procedure does not follow the fate of each mutation at a particular locus, which is extremely time consuming, it gives the same results as the multilocus stochastic simulations of Charlesworth and Charlesworth (1997) as far as the estimation of the time between clicks of the ratchet is concerned, for all parameter sets tested (results not shown).
SIMULATION RESULTS
If the loss of the least-loaded class is the determining factor in driving the ratchet (Haigh 1978), one would expect that the time between clicks of the ratchet would stay approximately constant over a range of parameter values that keep n0 constant. Figures 1, 2, and 3 show the simulation results for several parameter sets chosen such that n0 stays constant.
In Figure 1, s is kept constant at 0.015, and N changes with u to keep n0 constant (either 20.2 or 202). We observe that the time between clicks of the ratchet does not change significantly over an order of magnitude change in N. The increase in N seems to be compensated by the increase in u. In Figure 2, u is kept constant, at 0.1, and N changes together with s to keep n0 constant (with the same values as in Figure 1). Although for small values of n0 there seems to be no significant difference in the speed of the ratchet over an order of magnitude change in population size, for a higher value of n0 the difference is evident: as s becomes large (small N in the plot), the speed of the ratchet is greatly reduced. For example, for a population size of 705 individuals, with a selection coefficient of 0.08, we did not observe any click over 50,000 generations, but for a population of 19,000 individuals, with s = 0.022, the average time for a click is ~1560 generations. In Figure 3, N is kept constant and u changes with s (the mean equilibrium number of mutations, u/s, has the value 5), so that n0 is constant. We see that, for either small u and small s or for large u and large s, the speed of the ratchet is greatly reduced.
Average time between clicks of the ratchet, ±2 SE, for constant u as a function of the population size (as N increases, s decreases so that n0 is maintained constant; u is set equal to 0.1). n0 is 202 for circles and 20.2 for triangles.
These results show that, as noted previously by Stephan et al. (1993), the size of the best class is not sufficient to predict the speed of the ratchet, because the ratchet can turn at very different speeds for a constant n0. One observes that for the same n0, keeping N constant and varying u and s so that u/s is constant, there is a value of u and s for which the time of the ratchet has a minimum, i.e., for the same n0, increasing s can both slow down and speed up the ratchet (see the U-shaped curves in Figure 3). This shape is not predicted by any of the previous formulas. In contrast, such a minimum is predicted by the approximation presented here, although in the region of very small u and very small s the approximation gives lower times than the simulations. It also underestimates the time between clicks of the ratchet for small population size (or small n0).
The possible reasons why this minimum is observed, and why the present approximation underestimates the time for small population size, deserve some comment. Other things being equal, decreasing s should speed up the ratchet and decreasing u should slow it down. In the case of Figure 3, both s and u are changing to keep u/s constant, so that a minimum may occur, due to the fact that the dependence of the time on the mutation rate is different from that on the selection coefficient. In the region where s is very small, so that each mutation has very little effect on fitness, u is also very small. This means that the probability of a mutation occurring is very small, and the force of mutation that drives individuals from the best class to the next class is greatly reduced, leading to a slower ratchet. In the region where the mutation rate is large, the selection coefficient is also large, so that although mutations keep appearing at a high rate, selection is so efficient in restoring the best class that a great number of individuals come from the least-loaded class, which leads to a slower ratchet.
Average time between clicks of the ratchet, ±2 SE, for constant n0 as a function of the selection coefficient (as s increases, u increases so that u/s is constant). The population size is 30,000 for circles (n0 = 202) and 3000 for triangles (n0 = 20.2).
One observes from the comparison of the theoretical formula and the simulation results that, as long as s is not extremely small or large (Figures 1, 2 and 3), the predictions seem to approximate the simulations reasonably well, especially when N is big (or n0 is large). If N is small, it is more difficult for the system to maintain itself close to equilibrium, because drift is dominating. In this case, an approximation based on small perturbations becomes inadequate.
We now ask how the speed of the ratchet changes with population size, for a given mutation rate and a constant selection coefficient. The simulation results and the expected times calculated with the various approximations discussed above are shown in Figures 4, 5 and 6, for different values of N, u, and s, as a function of n0. The mutation rate and selection coefficient were chosen to lie in the parameter range that may be most relevant to the problem of the evolutionary degeneration of an incipient Y chromosome or neo-Y chromosome (Charlesworth 1978, 1996; Charlesworth and Charlesworth 1997). The values of u = 0.04 and s = 0.01 might be considered reasonable for the case of a Drosophila chromosome arm, if widely accepted estimates of mutational parameters are used (Charlesworth and Charlesworth 1998; Drakeet al. 1998). In each figure, the time between clicks of the ratchet is plotted against n0, for a fixed value of u and s, so that an increase in n0 is solely due to an increase in N. The population size varies in the interval [1000; 13,500] in Figure 4, [1000; 10,000] in Figure 5, and [1000; 30,000] in Figure 6.
Mean time between clicks of the ratchet (Trat), ±2 SE, as a function of n0 (this is varied by varying N). The mutation rate is 0.075 and the selection coefficient is 0.02. The results of the approximations of Charlesworth and Charlesworth (1997) are labeled Tc and of Stephan et al. (1993) are labeled Ts and Ts14 for their Equations 8 and 14, respectively. The result from the approximation presented here is labeled Tp.
The same as in Figure 4, but for a mutation rate of 0.05 and selection coefficient of 0.015.
The same as in Figure 4, but for a mutation rate of 0.04 and a selection coefficient of 0.01.
One observes that the approximation of Charlesworth and Charlesworth, contrary to what was previously thought (Charlesworth and Charlesworth 1997; Orr and Yuseob 1998), greatly underestimates the speed of the ratchet for large population sizes in the parameter range considered here. Although Stephan et al. could not establish exactly the range of validity of their two approximations, they suggested use of their Equation 8 for predicting the speed of the ratchet for the range of selection coefficients considered here. From comparison with the simulations presented, we see that their Equation 14 seems to describe the rate of change with N of the time between clicks of the ratchet better than their Equation 8, although it always underestimates the absolute time for this parameter range.
In simulations done to check the change in the ratchet's speed with different mutation rates (Figure 7), we can see that the range of parameters for which Stephan et al.'s Equation 14 gives a better quantification of the process than their Equation 8 is not only dependent on a large population size and strong selection, but is also modulated by the mutation rate. In simulations done to check the change of the ratchet's speed with different s (Figure 8) we see that, for a given N and u, as s increases, the estimated time given by their Equation 14 becomes an underestimate. In fact, Stephan et al. studied the process by dividing it in two separate phases: the establishment phase, during which the new least-loaded class reaches a value close to that of the deterministic equilibrium, and the extinction phase, during which the least-loaded class becomes extinct. Their Equation 14 is based on the assumption that the population size is large enough that the change in mean fitness of the population is sufficiently small to be approximated by small perturbations and that selection is strong enough that the process is mostly trapped in the extinction phase. This assumption is not supported by our simulations (see Figure 9). The time spent in the establishment phase is the one above the dashed line; the time spent in the extinction phase is the one below this line.
The time between clicks of the ratchet, ±2 SE, as a function of the mutation rate (N = 10,000, s = 0.015).
The time between clicks of the ratchet, ±2 SE, as a function of the selection coefficient (N = 10,000, u = 0.075).
The fluctuations in the size of the least-loaded class through time, for N = 10,000, u = 0.075, and s = 0.02, for a single simulation run. The dashed line indicates the equilibrium value of the least-loaded class.
When deriving our approximation, we assumed that the net loss of mean fitness due to a click of the ratchet would be K ≈ 0.6se−u, because the system did not recover its new equilibrium instantaneously. In Figure 10A, we plot the mean fitness of the population as a function of 1 − x/x0 for a set of simulation runs with parameters N = 10,000, u = 0.04, s = 0.01, after a click of the ratchet. As we assumed that the mean fitness of the population could be approximated as a linear function of 1 − x/x0, the slope of the linear regression line plotted in the figure corresponds to the value of K in the theoretical approximation. The value assumed in the derivation (K = 5.76 × 10−3) agrees well with the one from the regression. However, the agreement is not good for a population size of only 1000 (Figure 10B). This can be attributed at least to two factors: either the population size is so small that we cannot approximate the changes in mean fitness by small perturbations and/or the value of K is different from the one we are assuming. If we calculate the time by substituting the value of K from the linear fit in the theoretical expression, we find that the time obtained is still below the one measured in the simulations, so that an incorrect value for K is not the only source of the discrepancy.
Population mean fitness as a function of 1 − x/x0 for a set of simulation runs with parameters u = 0.04, s = 0.01, N = 10,000 (A) and N = 1000 (B). The line plotted is a linear regression line whose slope represents the value of K assumed in the theoretical approximation.
DISCUSSION
Speed of the ratchet: The equilibrium size of the least-loaded class, n0, is usually regarded as the chief parameter that determines the speed of the ratchet (Haigh 1978; Bell 1988; Gessler 1995; Charlesworth and Charlesworth 1998). The simulations presented here show that this parameter is not sufficient for this purpose, as noted previously by Stephan et al. (1993). For example, for an equilibrium size of the best class of ~200 we can get average times for one click of the ratchet varying from 900 to 8000 as a result of changes in selection coefficient (0.015–0.04) and mutation rate (0.075–0.2).
The parameter space considered here has been constrained to ensure n0 ⪢ 1, so that the Poisson distribution expected under mutation-selection balance can be approximately attained. Gessler (1995) has shown that this balance will not be met for conditions under which n0 < 1, because selection is too weak to counter mutation pressure. In this case, the distribution of the number of mutations is not Poisson but is close to a shifted negative binomial distribution, whose parametrization allows an estimation of the rate of the ratchet.
For ⪢ 1, the approximation for the advance of the ratchet based on a diffusion equation presented here seems to make a better prediction of the time between clicks of the ratchet than the previous approximations, for moderate selection coefficients in a range that is compatible with the biological data, provided the population size is not too small, so that the stochastic fluctuations are not too violent. As noted before by Stephan et al. (1993), for intermediate selection coefficients the establishment phase and the extinction phase are blurred and a separate analysis of these two phases does not well predict the outcome. An approximation based on the assumption that the mean fitness of the population is affected solely by fluctuations of the least-loaded class
The observation that, for a constant n0 and constant s, we did not observe significant differences in the speed of the ratchet, over an order of magnitude change in the population size (Figure 1), suggests that n0s is an important parameter, although the expression for the average time between clicks of the ratchet is not an explicit function of n0s. For the parameter range considered here we observe that an increase of 10-fold in n0s caused a decrease of ~10-fold in the speed of the ratchet.
Y and neo-Y chromosome degeneration: Because an incipient Y chromosome, or a neo-Y chromosome resulting from an autosome fusion or translocation, that fails to recombine with its homologue in the heterogametic sex is vulnerable to the ratchet, it is interesting to calculate the expected rate of its operation under the above approximation. The erosion of a proto-Y chromosome is very similar to the degeneration of a haploid asexual population if one replaces s by hs, where h is the dominance coefficient and s is the selection coefficient against homozygous mutations (Charlesworth and Charlesworth 1997). In the case of Drosophila, if one assumes an effective population size of males of 5 × 105, a deleterious mutation rate per Y chromosome of 0.04, and an average selection coefficient against a heterozygous mutation of 0.01, the present approximation gives a value of ~3 × 1025 generations for one click of the ratchet; if one sets 5 generations per year for Drosophila, this will correspond to 6 × 1024 years per click. However, if the mutation rate is slightly higher, say 0.06, this time would decay to 40,000 years per click. Our approximation suggests that an increase in the mean number of mutations from 4 (n0 = 9158) to 6 (n0 = 1239), just due to an increase in the mutation rate, increases the speed of the ratchet by 20 orders of magnitude. This comes from the fact that there is a more or less exponential increase of the time between clicks of the ratchet with an increase in n0 for a given selection coefficient. On the other hand, for a mutation rate of 0.04, if hs is ~0.005, n0 will be 168 and the estimated average time for a click is 697 generations. Because every click of the ratchet in a haploid population indirectly leads to the fixation of a deleterious mutation in the whole population (Charlesworth and Charlesworth 1997), this means that if the ratchet is clicking every 697 generations, in a period of 1 million years we would expect an incipient Drosophila Y chromosome to become fixed for ~7170 deleterious mutations (again assuming 5 generations per year).
The neo-Y chromosome system of Drosophila miranda constitutes an excellent clock to set the time scale over which the degeneration of a nonrecombining region is supposed to occur. The time estimated for the origin of the chromosomal rearrangement generating the neo-Y in D. miranda is ~1.25 million years ago (mya; S. Yi, personal communication). The neo-Y shows evidence for degeneration, and the neo-X is partially dosage compensated (Steinemannet al. 1993; Steinemann and Steinemann 1998, 1999). These observations seem to suggest that, if there is a general process responsible for the degeneration of the nonrecombining segment of the genome such as Muller's ratchet, it is expected to show its signature over a time scale of the order of 106 years.
In Table 1 we show the expected time for a click of the ratchet, under the approximation proposed here, for a population of half a million males for various values of u and hs. We observe that, for these values for the average effect of a mutation, the time for a click of the ratchet becomes biologically irrelevant when n0hs goes above 15. If the average heterozygous effect of a nonlethal deleterious mutation is of the order of 1% (Charlesworth and Hughes 1999), under the above approximation, for n0 < 1500 the ratchet may play a role in the degeneration of the neo-Y, but its rate is probably too small to explain the degeneration observed if n0 > 1000. For the ratchet to be the main process causing the degeneration of the neo-Y, n0 has to be probably <500. In the case of lethal mutations, which probably occur at a rate of ~0.0025 for an incipient Y chromosome (Fryet al. 1999) and with hs of ~2% (Crow 1993), the time for a click of the ratchet in a population of half a million chromosomes is biologically irrelevant.
Muller's ratchet is a priori more likely to be an important force in driving the degeneration of mammalian Y chromosomes, given that their effective population size is around one order of magnitude less than that for Drosophila (Charlesworth and Charlesworth 1997). There is some evidence that the evolution of the mammalian Y chromosome has been punctuated by at least four events that suppressed recombination between the X and the Y (Lahn and Page 1999), the first event having occurred ~300 mya. If the deleterious mutation rate and average effect of mutations for a mammalian proto-Y chromosome were the same as those estimated for the Y in Drosophila, then with a population of 5 × 104, proto-Y chromosomes would degenerate due to this process at an average rate of 1 click every 40 thousand generations, for a deleterious mutation rate of 0.04. Our ignorance of the value of these parameters does not allow any final conclusion about the ratchet being a leading process in the degeneration of the mammalian Y, although it seems more likely than for the Drosophila case.
Time for a click of the ratchet for a hypothetical neo-Y chromosome
Because of the assumptions we have made to derive these results, some caution has to be taken considering their implications. First we must note that we have assumed, for simplicity, that all deleterious mutations have the same effect. However, recent work has suggested that an equal effect of mutations assumption does not fit the data from mutation accumulation experiments, which are designed to measure the mutation rate to deleterious mutations and the selection coefficients against those mutations (Keightley 1996; Davieset al. 1999; Fryet al. 1999). The occurrence of many mutations with small effects, and a few with large effects, seems to be more consistent with the results. If this is the case, the ratchet is expected to turn at a much higher speed than for a single selection coefficient of the order of 1% (Gessler and Xu 1999), but each turn will cause a very small decline in mean fitness if many mutations have very low selection coefficients. It is likely that the degeneration of the Y chromosome and the evolution of dosage compensation are driven by selection to increase the activity of the X relative to the Y in males, in response to the decline in mean fitness of the Y (Charlesworth 1996; Charlesworth and Charlesworth 1997). The strength of such selection is determined by the rate of this decline and will be very weak if it is small.
We have also assumed independence of mutational effects and it has been shown that epistasis slows down the speed of this process (Charlesworthet al. 1993; Kondrashov 1994). However, if there is in fact a distribution of mutational effects with a more or less exponential shape, epistasis will not stop the ratchet (Butcher 1995).
Acknowledgments
We thank D. Gessler and E. Baake for helpful comments. I. Gordo also thanks G. McVean for helpful discussions concerning this work. I. Gordo is supported by the Gulbenkian Foundation and Program PRAXIS XXI of Portugal, and B. Charlesworth by the Royal Society.
Footnotes
-
Communicating editor: W. Stephan
- Received August 23, 1999.
- Accepted November 22, 1999.
- Copyright © 2000 by the Genetics Society of America
