Denis Roze
email: roze@sb-roscoff.fr

"Sex uncovered" special issue: Selection for sex in finite populations

Keywords: deleterious mutation, evolution of sex, genetic drift, modifier model, multilocus model, selection

Finite population size generates interference between selected loci, which has been shown to favour increased rates of recombination. In this paper, I present different analytical models exploring selection acting on a "sex modifier locus" (that affects the relative investment into asexual and sexual reproduction) in a finite population.

Two forms of selective forces act on the modifier: direct selection due to intrinsic costs associated with sexual reproduction, and indirect selection generated by one or two other loci affecting fitness. The results show that indirect selective forces differ from those acting on a recombination modifier even in the case of a haploid population: in particular, a single selected locus generates indirect selection for sex, while two loci are required in the case of a recombination modifier. This effect stems from the fact that modifier alleles increasing sex escape more easily from low-fitness genetic backgrounds than alleles coding for lower rates of sex. Extrapolating the results from three-locus models to a large number of loci at mutation-selection balance indicates that in the parameter range where indirect selection is strong enough to outweigh a substantial cost of sex, interactions between selected loci have a stronger effect than the sum of individual effects of each selected locus. Comparisons with multilocus simulation results show that such extrapolations may provide correct predictions for the evolutionarily stable rate of sex, unless the cost of sex is high.

INTRODUCTION

Almost 50 years ago, [START_REF] Crow | Evolution in sexual and asexual populations[END_REF] published a mathematical treatment of the idea previously proposed by [START_REF] Fisher | The genetical theory of natural selection[END_REF] and [START_REF] Muller | Some genetic aspects of sex[END_REF] that sex increases the rate of adaptation of species by allowing different beneficial mutations to be combined into the same genome. A controversy ensued shortly after, as a different model analyzed by [START_REF] Smith | Evolution in sexual and asexual populations[END_REF] led to a contradictory result: beneficial alleles at two different loci spread at the same rate in sexuals and asexuals. The debate was finally settled by [START_REF] Felsenstein | The evolutionary advantage of recombination[END_REF], who showed that a critical difference between the two models concerns population size: while Crow and Kimura assume that new mutations arise as single copies, Maynard Smith's model considers an infinite population, in which beneficial alleles are initially in linkage equilibrium. As demonstrated by [START_REF] Hill | The effect of linkage on limits to artificial selection[END_REF] a few years before, stochastic effects occurring in any finite population tend to generate negative associations between beneficial alleles (i.e., beneficial alleles tend to be found more often in different genetic backgrounds than combined within the same background). In sexual populations, these negative associations are broken by recombination, which therefore increases the rate of adaptation.

The same type of mechanism allows sexuals to escape Muller's ratchet, the gradual accumulation of deleterious alleles that occurs in non-recombining populations (Muller, 1964;[START_REF] Felsenstein | The evolutionary advantage of recombination[END_REF].

Since Crow and Kimura's paper, the question of the possible advantage of sexual reproduction has been one of the most hotly debated topics in evolutionary biology.

The fact that sex is advantageous in the long term does not in principle prevent an asexual mutant from invading a sexual population: indeed, many strong costs are associated with sexual reproduction, in terms of energy required to find or court a mate, increased risks of predation and disease transmission, or investment into the male function -since the males of most species do not contribute to the next generation in terms of resources (e.g., Maynard Smith, 1971;[START_REF] Lewis | The cost of sex[END_REF][START_REF] Lehtonen | The many costs of sex[END_REF]. Strong benefits of sex are thus needed to counterbalance these costs. Insights about the possible strength of indirect selection for sex (stemming from the reduction of negative genetic associations) can be obtained from recombination modifier models, such as the one originally proposed by [START_REF] Nei | Modification of linkage intensity by natural selection[END_REF]: these models represent a "modifier" gene that affects recombination rates between other genes, at which beneficial/deleterious alleles segregate. The direction and magnitude of indirect selection acting at the modifier locus can then be explored mathematically under different scenarios. Using this approach, different conditions that may favour sexual recombination have been described: selection acting independently at different loci within finite or spatially structured populations (due to the Hill-Robertson effect described above, e.g., [START_REF] Felsenstein | The evolutionary advantage of recombination. II. Individual selection for recombination[END_REF][START_REF] Otto | The advantages of segregation and the evolution of sex[END_REF]Barton, 1997, 2001;[START_REF] Iles | Recombination can evolve in large finite populations given selection on sufficient loci[END_REF][START_REF] Barton | Evolution of recombination due to random drift[END_REF][START_REF] Keightley | Interference among deleterious mutations favours sex and recombination in finite populations[END_REF][START_REF] Martin | Selection for recombination in structured populations[END_REF][START_REF] Barton | The Hill-Robertson effect and the evolution of recombination[END_REF][START_REF] Gordo | Sex and deleterious mutations[END_REF][START_REF] Hartfield | The role of advantageous mutations in enhancing the evolution of a recombination modifier[END_REF], negative epistasis between selected loci [START_REF] Feldman | Evolution of recombination in a constant environment[END_REF]Charlesworth, 1990Charlesworth, , 1993a;;[START_REF] Barton | A general model for the evolution of recombination[END_REF][START_REF] Otto | Deleterious mutations, variable epistatic interactions, and the evolution of recombination[END_REF], fluctuating selection, generated for example by biotic interactions [START_REF] Charlesworth | Recombination modification in a fluctuating environment[END_REF][START_REF] Barton | A general model for the evolution of recombination[END_REF][START_REF] Gandon | The evolution of sex and recombination in response to abiotic or coevolutionary fluctuations in epistasis[END_REF], and spatially varying selection within structured populations [START_REF] Lenormand | The evolution of recombination in a heterogeneous environment[END_REF]. In many of these models, selection for recombination can be strong in non-recombining populations, but decreases rapidly as the baseline recombination rate increases.

There are important differences between the evolution of sex and the evolution of recombination, however: in particular, reproductive modes are probably far more constrained than recombination rates within genomes. In many obligatory sexual species, reproducing asexually may require a combination of several mutations which may be individually deleterious: in mammals for example, parthenogenesis would require both the spontaneous development of an ovule without fertilization and the suppression of genomic imprinting (which prevents normal development in the absence of a paternally inherited genome), which may explain why it has never been observed.

Similarly, in cyclical parthenogenetic species such as aphids, cladocerans and rotifers, overwintering eggs are typically produced by sexual reproduction (e.g., [START_REF] Simon | Ecology and evolution of sex in aphids[END_REF], and it has been argued that sex may be maintained for this reason -although it is less clear why overwintering eggs could not be produced asexually, as this occurs for example in some Daphnia pulex populations [START_REF] Hebert | Clonal coexistence in Daphnia pulex (Leydig): another planktonic paradox[END_REF][START_REF] Innes | The origin and genetic basis of obligate parthenogenesis in Daphnia pulex[END_REF]. Based on this type of observation, several authors have proposed that sex may be maintained by selection at the species level [START_REF] Williams | Sex and Evolution[END_REF][START_REF] Nunney | The maintenance of sex by group selection[END_REF][START_REF] De Vienne | Lineage selection and the maintenance of sex[END_REF]: according to this hypothesis, only the lineages in which viable asexual mutants cannot occur (due to some constraint) persist over evolutionary time, while asexual lineages eventually go extinct due to mutation accumulation and/or failure to adapt. Although this type of process can certainly explain the maintenance of sex in some groups (such as mammals), it cannot explain evolution towards higher rates of sex (and in the extreme, obligatory sex) within a lineage: for this, sex must confer benefits on a sufficiently short timescale so that its direct costs can be outweighed. Assessing whether this is possible or not requires obtaining quantitative predictions on the possible strength of selection for sex under different scenarios.

In this paper, I explore the advantage of sex generated by selection acting in-dependently at multiple loci (no epistasis) within finite populations. To date, most analytical results on the effects of finite population size stem from haploid recombination modifier models representing selection acting at two loci (e.g., [START_REF] Otto | The evolution of recombination: removing the limits to natural selection[END_REF][START_REF] Barton | Evolution of recombination due to random drift[END_REF][START_REF] Martin | Selection for recombination in structured populations[END_REF][START_REF] Barton | The Hill-Robertson effect and the evolution of recombination[END_REF]. While indirect selective pressures acting on sex have been shown to be different from those acting on recombination in diploid organisms -since sex allows both recombination and segregation (affecting homozygosity/heterozygosity at single loci, e.g., [START_REF] Uyenoyama | On the origin of meiotic reproduction: a genetic modifier model[END_REF][START_REF] Otto | The advantages of segregation and the evolution of sex[END_REF]Agrawal, 2009a,b;[START_REF] Michod | Deleterious mutations and selection for sex in finite, diploid populations[END_REF], in haploids selection for recombination is often thought to be equivalent to selection for sex.

In the first part of the paper, however, I show that selective pressures acting on a sex modifier locus (affecting the rate of sexual vs. asexual reproduction) are different from those acting on a recombination modifier, even in the case of a haploid population: in particular, while two selected loci are required to generate indirect selection on a recombination modifier, a single selected locus can favour sex within a finite haploid population. In a second part, I use individual-based simulations to explore how the rate of sex evolves when deleterious mutations occur at a large number of loci (in the presence of a direct cost of sex), and compare the results with analytical predictions obtained by extrapolations from three-locus models. The results show that such extrapolations give correct predictions for the equilibrium rate of sex in different cases (in particular as long as the cost of sex is not too high) and capture the qualitative effects of the different parameters of the model (population size, number of loci, strength of selection against deleterious alleles).

METHODS

The sex modifier model. The effects of one or two selected loci on the evolution of a modifier gene affecting the relative investment into sexual vs. asexual reproduction can be studied using analytical methods (Table 1 shows the parameters used in the different models). The basic model represents a population of N haploid individuals with discrete generations. As in previous models [START_REF] Uyenoyama | On the origin of meiotic reproduction: a genetic modifier model[END_REF][START_REF] Otto | The advantages of segregation and the evolution of sex[END_REF]Agrawal, 2009a,b;[START_REF] Michod | Deleterious mutations and selection for sex in finite, diploid populations[END_REF][START_REF] Otto | Differential selection between the sexes and selection for sex[END_REF]) I assume that each individual i invests a proportion σ i of its resources into sexual reproduction (producing both male and female gametes) and a proportion 1 -σ i into clonal reproduction. A parameter c represents the overall effect of direct costs associated with sex: the probabilities that individual i is the maternal parent of a given individual of the next generation through asexual and sexual reproduction are proportional to 1 -σ i and σ i /c, respectively (c = 1 in the absence of a direct cost of sex). This cost may stem from the fact that a fraction of resources invested into sex are used to produce male gametes, which do not contribute to the next generation in terms of resources ("cost of males", c = 2 under equal investments into the sexual female and male functions) and/or other costs such as additional resources that may be needed to facilitate gamete encounter. Individuals produce asexual spores and gametes by mitosis; spores develop directly into new (haploid) individuals, while gametes fuse at random to form zygotes, which immediately undergo meiosis to produce new haploids. I assume that each parent produces a very large (effectively infinite) number of juveniles, N individuals being then sampled randomly among all juveniles to form the next adult generation.

Investment into sex is controlled by a single "modifier" locus with two alleles m and M : the rates of sex of m and M individuals are σ and σ + δσ, respectively.

Either one or two selected loci control the overall fecundity (total investment into reproduction) of individuals. In the one-locus case, the two alleles at the selected locus are denoted a and A, the fecundities of a and A individuals being proportional to 1 and 1 + s. With two selected loci, I denote b and B the two alleles at the second selected locus, and suppose for simplicity that allele B has the same selection coefficient s as allele A (although the model is easily extended to the more general case of different selection coefficients); the effects of the two loci on fecundity are supposed to be multiplicative (no epistasis). Finally, r ij denotes the recombination rate between loci i and j.

Quasi-linkage equilibrium (QLE) approximation. Approximate expressions for the expected change in allele frequency at the modifier locus in terms of the model parameters and expected genetic variances at the different loci can be obtained using a separation of timescales argument, which assumes that the different forces causing changes in allele frequencies (selection, drift, modifier effect) are weak relative to "effective" rates of recombination (which depend on the rate of sex in the population) breaking linkage disequilibria between loci. The analysis proceeds in two steps: the first step corresponds to the calculation of the expected change in frequency at the modifier locus over one generation, as a function of various moments of allele frequencies and linkage disequilibria in the parental generation; during the second step, a QLE approximation is used to express all moments involving linkage disequilibria in terms of moments involving allele frequencies only. For this, recursions for these different moments must be derived, and solved at equilibrium. A general method for deriving recursions on moments under a diploid life cycle is given in Appendix A of Roze and [START_REF] Michod | Deleterious mutations and selection for sex in finite, diploid populations[END_REF]. Appendix A of the present article briefly sketches the method for the case of a haploid life cycle, using the same general principles. Defining X i an indicator variable that equals 1 if allele i (M , A or B) is present in a given individual and 0 otherwise, the frequency of allele i in the population is given by p i = E [X i ], where E stands for the average over all individuals. Pairwise linkage disequilibria are defined as

D ij = E [ζ i ζ j ]
, where ζ i = X i -p i , while the three-locus linkage disequilibrium (e.g., [START_REF] Slatkin | On treating the chromosome as the unit of selection[END_REF] is given by

D MAB = E [ζ M ζ A ζ B ].
Throughout the paper, M t stands for the expected value of the moment M (such as D MAB D AB or p M D MA ) among individuals of generation t. Appendix A sketches how the life cycle can be decomposed into different steps in order to compute moments at generation t+1 in terms of moments at generation t (these recursions are implemented in a Mathematica notebook, available as Supplementary Material). When s, δσ and 1/N are small relative to recombination rates among loci multiplied by the rate of sex σ, one can show that moments involving linkage disequilibria quickly become small and change slowly over time [START_REF] Nagylaki | The evolution of multilocus systems under weak selection[END_REF][START_REF] Barton | Natural and sexual selection on many loci[END_REF][START_REF] Barton | Evolution of recombination due to random drift[END_REF][START_REF] Michod | Deleterious mutations and selection for sex in finite, diploid populations[END_REF]; in this case, these moments can be expressed in terms of the different parameters of the model and moments of allele frequencies, such as p M q M t or p M q M p A q A t (where q i = 1 -p i ). This QLE approximation simplifies considerably the equations, but is expected to break down when effective rates of recombination are small. However, we will see that accurate approximations can be obtained when recombination rates and the strength of selection s are of the same order of magnitude, in the situation where deleterious alleles are maintained at mutation-selection balance at selected loci: indeed in this case, changes in allele frequencies are only driven by the modifier effect and by drift, and the QLE only requires that δσ and 1/N are small relative to effective rates of recombination.

Multilocus simulations. As we will see in the next section, analytical results from three-locus models can be extrapolated to predict the overall strength of selection on a sex modifier when deleterious alleles segregate at a large number of loci; however, this supposes that the effects of higher-order interactions (between three or more selected loci) can be neglected, which is not guaranteed: even if the effect of each of these higherorder interactions is small, the number of such interactions increases very rapidly as the number of polymorphic loci increases. In order to test the predictions obtained from three-locus models, I use individual-based simulations representing deleterious mutations occurring at a large (quasi-infinite) number of possible loci. The simulation program (written in C++, and available from Dryad) is very similar to the programs used in previous papers (Roze [START_REF] Michod | Deleterious mutations and selection for sex in finite, diploid populations[END_REF][START_REF] Otto | Differential selection between the sexes and selection for sex[END_REF]. It represents a population of N haploid individuals whose genome consists in a linear chromosome with infinitely many loci. Every generation, the number of new deleterious mutations per chromosome is sampled from a Poisson distribution with parameter U , and the position of each mutation is sampled from a uniform distribution. The fecundity W i of an individual is given by (1 -s) n , where n is the number of mutations in its genome.

When forming the next generation, each parental individual i contributes as a mother (either through sexual or asexual reproduction) with a probability proportional to

W i [c (1 -σ i ) + σ i ]
(where σ i is the rate of sex of parent i). The new individual is then produced asexually with probability c (1 -σ i ) / [c (1 -σ i ) + σ i ] (in which case its genotype is a copy of the parental genotype), and sexually with the complementary probability. In this last case, a second individual is sampled to contribute as a father, the probability of sampling parent j being proportional to W j σ j . The two parental genomes are recombined to form a new haploid genome, the number of cross-overs being sampled from a Poisson distribution with parameter R (genome map length), and the position of each cross-over being sampled from a uniform distribution. During a given number of preliminary generations (usually 2000), the rate of sex is the same for all individuals (and is generally fixed to 1). Then, mutations are introduced (at rate 10 -4 per generation) at a sex modifier locus located at the midpoint of the chromosome.

Alleles at this locus take values between 0 and 1, corresponding to the rate of sex of the individual (which is thus a quantitative trait). When a mutation occurs at the modifier locus, with probability 0.5 the new value is sampled from a uniform distribution between 0 and 1, while with probability 0.5 it equals the value of the parent allele plus a random number between -0.1 and 0.1 (the new value being set to 0 or 1 if it is negative or greater than 1, respectively); this mutation scheme allows both large and small effect mutations to arise relatively frequently. The program runs for 2 million generations and records the average value of alleles at the modifier locus every 100 generations; the equilibrium rate of sex is obtained by averaging over the last 1.9 × 10 6 generations. In some simulations I only considered two alleles at the modifier locus, as will be explained in greater detail below.

RESULTS

All analytical results assume that selection is weak (s small) and are derived to the first order in δσ and 1/N (weak modifier, large population size). The effect of relaxing these hypotheses will then be explored by multilocus simulations. Through the following, σ eff denotes the effective rate of sex in the population (proportion of sexually produced progeny) which equals σ/ [c (1 -σ) + σ] to leading order, while ρ ij = σ eff r ij denotes the effective recombination rate between loci i and j (proportion of recombinants per generation). Finally, δσ eff = δσ/ [c (1 -σ) + σ] denotes the effect of allele M on the effective rate of sex. Two-locus model. Although the results shown in this section hold for both s > 0 and s < 0, I will assume for clarity that s > 0 (A is the favoured allele) when discussing the results. Using the methods described in Appendix A (implemented in the Supplementary Material), one obtains for the expected change in frequency of the modifier over one generation:

∆p M t = -δσ eff (c -1) p M q M sel t + s D MA t -s 2 p A D MA t (1)
where p M q M sel t is the expected genetic variance at the modifier locus among selected parents (weighting each parent by its relative fecundity, see Appendix A), given by:

p M q M sel t = p M q M t -2s p M D MA t -s 2 D MA 2 t .
(2)

As we will see below, p M D MA t is positive at QLE, and equation 2 represents the fact that the genetic variance at the modifier locus tends to be reduced by selection at the other locus (e.g., [START_REF] Smith | The hitch-hiking effect of a favourable gene[END_REF]Charlesworth, 1993b): p M q M thus tends to be lower when measured among selected parents than when measured before selection. The first term of equation 1 represents the effect of the direct cost of sex, reducing the frequency of the allele coding for more sex (this reduction being proportional to the genetic variance at the modifier locus among selected parents). The second term of equation 1 represents the effect of genetic hitchhiking, favouring the modifier allele that tends to be positively linked (on average) with the better allele (A)

at the selected locus. Finally, the last term of equation 1 represents the fact that for a given value of D MA , the strength of hitchhiking tends to decrease when the beneficial allele is more frequent in the population; this stems from the fact that hitchhiking is proportional to the ratio s/W (which determines the selection gradient at the selected locus), and that mean fecundity W is an increasing function of p A .

Further progress can be made by obtaining QLE expressions for the four moments D MA 2 , p M D MA , p A D MA and D MA that appear in the equations above.

The first of these moments is the variance in linkage disequilibrium, generated by finite population size. To leading order, it is given by:

D MA 2 t = p M q M p A q A t N ρ MA (2 -ρ MA ) (3)
(see also eq. B3a in [START_REF] Barton | Evolution of recombination due to random drift[END_REF]. The moment p M D MA is generated by finite population size and selection for allele A. To the first order in s, its expression at QLE is:

p M D MA t = s 1 -ρ MA ρ MA D MA 2 t (4)
(see also eq. B3b in [START_REF] Barton | Evolution of recombination due to random drift[END_REF]. This result can be interpreted as follows:

due to drift, allele M is sometimes associated with allele A (D MA > 0), sometimes with allele a (D MA < 0). Because A is beneficial, M tends to increase in frequency when D MA > 0, and decrease in frequency when D MA < 0, thus generating a positive covariance between allele frequency p M and linkage disequilibrium D MA . The same mechanism generates a covariance between p A and D MA in the presence of a direct cost of sex (c > 1), that has the sign of -δσ. At QLE, it is given by:

p A D MA t = -δσ eff (c -1) 1 -ρ MA ρ MA D MA 2 t . (5)
Finally, the expected value of the linkage disequilibrium at QLE is:

D MA t = 1 ρ MA 2 (1 -ρ MA ) δσ eff (c -1) p M D MA sel t -s p A D MA t + δσ eff [r MA + (c -1) ρ MA ] p M D MA sel t (6)
where p M D MA sel t is measured after selection, and at QLE is equal to s D MA

2 t /ρ MA .
The terms on the first line of equation 6 represent linkage disequilibrium generated by the Hill-Robertson effect between the modifier and the selected locus: finite population size tends to generate negative LD between deleterious alleles as different loci [START_REF] Hill | The effect of linkage on limits to artificial selection[END_REF]. Therefore, when allele M is deleterious due to the cost of sex (δσ > 0), it tends to be associated with the better allele A at the second locus. This effect disappears in the absence of a direct cost of sex (c = 1). By contrast, the term on the second line of equation 6 does not vanish when c = 1, indicating that linkage disequilibrium is generated between the two loci even in the absence of direct selection at the modifier locus. Like the first term, this second term has the sign of sδσ,

indicating that an allele increasing sex tends to be associated with the better allele (A) at the selected locus; furthermore, it is generated by the positive covariance between p M and D MA . This effect may be understood as follows: whenever D MA is negative, it tends to reduce the frequency of M , in which case an M individual that reproduces sexually has higher chances to mate with an m individual, and thus to recombine onto a better genome (carrying A, since m tends to be more associated with A). By contrast, a positive D MA increases the frequency of M , in which case M individuals have greater chances to mate with other M individuals, and thus to stay associated with the better allele. Averaging over these different situations, the fraction of alleles M that become disassociated from the beneficial allele (and associated with the deleterious allele)

when recombining with m individuals is thus lower than the fraction of alleles M that become disassociated from the deleterious allele (and associated with the beneficial one) when recombining with m (again because m tends to be less frequent when M is associated with A). Although both modifier alleles benefit from this effect, the allele that engages in sex more often (M if δσ > 0) benefits more, generating a positive association between this allele and the beneficial allele A. This mechanism should not work, however, when one of the modifier alleles never engages in sex (σ = 0).

Although equation 6 does not capture this (as it is derived under the assumption that δσ σ, which is necessary for the QLE approximation to hold), one can show that the recursion for D MA t becomes independent of δσ when σ = 0 and c = 1, confirming that no linkage disequilibrium is generated (on average) when one of the modifier alleles codes for fully asexual reproduction.

From the equations above, the expected change in frequency of the modifier can be written as:

∆p M t = -δσ eff (c -1) p M q M t -2s 2 2 -ρ MA ρ MA 2 D MA 2 t + δσ eff s 2 ρ MA σ eff D MA 2 t . (7)
The term on the first line indicates that interference between the two loci due to finite population size (through the variance in LD) tends to attenuate the effect of the direct cost of sex (Hill-Robertson effect). More interestingly, the second line shows that in the absence of cost, a modifier increasing sex tends to increase in frequency, due to the effect described just above. Note that this effect does not operate in the case of a recombination modifier, which is only expressed during the diploid phase, once individuals have "decided" to engage in sex (because recombination only has a genetic effect in double heterozygotes, the average benefits gained by allele M should thus be exactly the same as those gained by allele m).

As shown by Figure 1, two-locus simulations confirm that in the absence of a direct cost of sex (c = 1), a modifier increasing sex tends to increase in frequency during the sweep of a beneficial allele at a second locus. In the simulations the beneficial allele A is initially in frequency p A,0 = 0.05, while the sex modifier M is initially in frequency p M,0 = 0.5 and in linkage equilibrium with A. The program measures the final frequency of M at the end of the sweep, the whole process being repeated 50×10 6 times. A prediction for the total increase in frequency of M over the sweep can be obtained by integrating equation 7 over time (treated as a continuous variable), and assuming that p A remains close to a deterministic trajectory. In this case, using the fact that dp A /dt = sp A q A we have

∞ 0 p A q A dt = 1 p A,0 dp A /s = (1 -p A,0 ) /s, which
yields for the expected total increase in frequency of M :

∆p M total ≈ δσ s (1 -p A,0 ) p M,0 (1 -p M,0 ) N ρ MA 2 σ (2 -ρ MA ) (8)
assuming no direct cost of sex (c = 1). Figure 1B shows that this prediction fits well with the simulation results as long as the baseline rate of sex σ in the population is not too small -note that with a single selected locus, the expected change in frequency of the modifier stays very small. Furthermore, Figure 1A confirms that this benefit of sex vanishes when the rate of sex of allele m tends to zero; however, even a small baseline rate of sex generates an advantage for a mutant increasing sex (the total change in frequency of the modifier has a maximum at around σ ≈ 0.01, for the parameter values used in Figure 1).

Three-locus model. In the presence of a second selected locus (where B is the beneficial allele), the expected change in frequency of the modifier becomes:

∆p M t = -δσ eff (c -1) p M q M sel t + s ( D MA t + D MB t ) + s 2 ( D MAB t -p A D MA t -p B D MB t ) -s 3 ( D MA D AB t + D MB D AB t + p A D MAB t + p B D MAB t ) -s 4 D MAB D AB t (9)
where the genetic variance at the modifier locus among selected parents is given by:

p M q M sel t = p M q M t -2s ( p M D MA t + p M D MB t ) -s 2 D MA 2 t + D MB 2 t + 2 D MA D MB t + 2 p M D MAB t -2s 3 ( D MAB D MA t + D MAB D MB t ) -s 4 D MAB 2 t . (10)
QLE expressions for the different moments that appear in equations 9 and 10 are given in Appendix B. From these, it can be shown that all the different terms involving the selected loci generate selection for increased rates of sex. Although it is not possible for space reasons to provide an intuitive explanation for each of the moments involved, it is interesting to note that the extra benefits of sex that appear when introducing a second selected locus do not all stem from recombination between the selected loci (which again differs from the case of a recombination modifier). In the following I assume that allele M increases sex (δσ > 0), and consider the effect of the three-locus disequilibrium D MAB t . As we will see, this disequilibrium is positive on average, meaning that allele M tends to be more often associated with the extreme genotypes ab and AB than allele m. This positive D MAB favors M through two effects: the first corresponds to the fact that under multiplicative selection, the mean fitness of extreme genotypes (1 + s + s 2 /2) is higher than the fitness of intermediate genotypes Mutation-selection balance. The expressions shown in the previous section and in Appendix B were derived under the assumption that δσ, 1/N s ρ X (where again ρ X is the effective rate of recombination between loci in the set X): indeed, all results

were expressed to the first order in δσ and 1/N (assuming that these parameters are much smaller than all other parameters), while the QLE approximation requires that selection is much weaker than recombination (s ρ X ). However, this last condition becomes less stringent when selected loci are at mutation-selection balance (equilibrium between mutation generating deleterious alleles and selection removing these alleles): indeed, in that case all changes in allele frequencies are due to drift and to the modifier effect, and the QLE therefore only requires that δσ, 1/N ρ X . More accurate results can thus be obtained for the case where effective recombination rates are small (of order s). As can be seen from Appendix B, all moments (such as D AB 2 , D AB , p M D MA ...) generated by selection and drift become of the same order of magnitude when s and ρ X are of the same order; similarly, moments involving the modifier effect δσ (such as D MAB D AB , D MAB , D MA ...) all become of the same order. Extra terms must be added to the recursions for these different moments when ρ X coefficients are of order s, however the changes are relatively minor if we assume that deleterious alleles remain at low frequency, so that we can neglect all moments involving p i 2 , where p i is the frequency of the deleterious allele at locus i ("rare alleles approximation"). Because it is more convenient to work with moments involving the frequencies of deleterious alleles (rather than the frequencies of favoured alleles), I assume from now on that A and B are the deleterious alleles, and denote s the selection coefficient against these alleles (ab, Ab, aB and AB genotypes thus have fitnesses 1, 1 -s, 1 -s and (1 -s) 2 , where s > 0). Using the methods described in Appendix A and assuming that effective recombination rates may be of order s yields the following recursion for the squared linkage disequilibrium between selected loci D AB 2 , to leading order:

D AB 2 t+1 = pq AB t N + (1 -ρ AB ) 2 (1 -4s) D AB 2 t ( 11
)
where pq AB = p A q A p B q B (note that the factor 1 -4s was neglected in the previous analysis). At equilibrium, we thus have:

D AB 2 t = pq AB t N 1 -(1 -ρ AB ) 2 (1 -4s) (12)
which is approximately pq AB t / [2N (ρ AB + 2s)] when s and ρ AB are small. Appendix C provides approximations for other moments when ρ X coefficients are of order s. In particular, the linkage disequilibrium D AB generated by the Hill-Robertson effect becomes:

D AB t ≈ - 2s 2 pq AB t N (ρ AB + 2s) 2 (ρ AB + 3s) . ( 13
)
Assuming pq AB t ≈ (u/s) 2 , where u is the per-locus mutation rate towards deleterious alleles, one obtains:

D AB t ≈ - 2u 2 N (ρ AB + 2s) 2 (ρ AB + 3s) . ( 14
)
Figure 2A shows that this prediction fits well with two-locus simulation results: as long as s > 0.001 (for the parameters used in Figure 2A), the linkage disequilibrium between deleterious alleles at mutation-selection balance becomes weaker as s increases, which

is not predicted by the result obtained assuming s ρ AB , as it neglects the terms in s in the denominator of equation 14 -note that the discrepancy between equation 14 and the simulation results observed for s < 0.001 stems from the fact that equation 14 assumes 1/N s. Selection for a modifier increasing sex through its effect in breaking the linkage disequilibrium D AB is proportional to s 2 D AB (Appendices B, C) and should thus be maximised for intermediate values of s, as shown by Figure 2B.

From the expressions derived in Appendices B and C and assuming that the baseline rate of sex σ is of order , one obtains that the strength of selection for sex generated by each selected locus i is of order δσ s 2 p i / [N 3 ], where p i is the frequency of the deleterious allele at locus i, while selection for sex generated by the interaction between loci i and j is of order δσ s 4 p i p j / [N 5 ]. Because p i , p j are small at mutationselection balance (assuming u s), the effect of between-locus interactions may thus seem negligible relative to the individual effect of each locus. However, we will see in the next section that this is no longer true when deleterious alleles segregate at many loci, so that the number of between-locus interactions becomes very large.

Extrapolation to many loci. The results shown in the previous sections (and in Appendices B and C) can be extrapolated to predict the overall strength of indirect selection acting on a sex modifier when deleterious mutations occur over a whole genome (with a genomic deleterious mutation rate U , and assuming for simplicity that all deleterious alleles have the same selection coefficient s). Neglecting the effects of interactions between more than two selected loci, one obtains the following expression for the selection gradient at the modifier locus (see Appendix C):

s M ≡ ∆p M t p M q M t = δσ eff -(c -1) + U s N Θ 1 + (U s) 2 2N Θ 2 (15)
The first term of equation 15 represents the effect of the cost of sex, the second term is the sum of the individual effects of all selected loci, and the third term is the sum of all effects of interactions between pairs of loci. Presumably, interactions between three, four or more selected loci could also generate indirect selection for sex (through terms in (U s) 3 , (U s) 4 ...), but the computation of these terms would be tedious. The terms Θ 1 and Θ 2 that appear in equation 15 are complicated functions of the effective rate of sex σ eff , cost of sex c, strength of selection s and recombination rates, integrated over all possible genomic positions of deleterious alleles (which requires performing two numerical integrations over the genetic map). When effective recombination rates are large relative to s (which requires σ eff s), Θ 1 and Θ 2 become approximately independent of s, and can be derived from the expressions given in Appendix B: in this case, equation 15 indicates that the strength of indirect selection for sex depends on the U s product (independently of the individual values of U and s) which is also the variance in fitness in the population at mutation-selection balance. However, this is no longer true when effective recombination rates are of order s, as terms in s appear in the denominator of expressions for genetic moments at QLE (Appendix C). In general, the results obtained assuming small σ eff (Appendix C) are still approximately valid when σ eff s, as they often closely match those obtained from Appendix B in the regime where σ eff s. Finally, relatively simple expressions can be obtained for Θ 1 and Θ 2 assuming free recombination among all loci (equations C19 to C22); as we will see, these are often close to the results obtained by integrating numerically over the genetic map, as long as the mean number of deleterious alleles per genome is not too large, and map length R not too small (so that most loci segregating for deleterious alleles are loosely linked). recombination between all loci is given by (see Online Appendix S1):

N e ≈ N exp - 4U s (σ eff + 2s) 2 . ( 16
)
In a linear genome with restricted recombination, N e is given by a more complicated expression which depends on the position of the locus and on map length R. At the modifier locus (mid-point of the chromosome), one obtains:

N e ≈ N exp - 2U s Rα 2 2 log β 2s + e R -1 σ eff α sβ , ( 17
)
with α = 2s + σ eff , β = αe R -σ eff .
As can be seen on Figure 3A, replacing N by N e in equation 15 (either using 16 or 17) provides a better prediction of the strength of indirect selection for sex, as long as the effective rate of sex is not too small; we will see below that this generalises to other parameter values. Finally, Figure 3B shows the three components of selection on the sex modifier that appear within the brackets of equation 15 (after replacing N by N e ). As can be seen on the figure, the sum of the individual effects of all loci is small relative to the effects of pairwise interactions among loci, in the parameter range where indirect selection becomes important relative to the cost of sex.

Evolutionarily stable rate of sex. The analytical model can be used to predict the evolutionarily stable rate of sex, corresponding to the value of σ towards which the population should evolve in the long term (assuming that evolution proceeds by small steps, since the model supposes that δσ is small): this corresponds to the value at which the curves cross the x-axis in Figure 3A, that is, the value of σ for which s M = 0.

Such predictions can then be compared with the results of simulations representing deleterious alleles occurring at a large number of loci, and where the rate of sex is free to evolve (as described in the Methods). Figure 4 shows the equilibrium value of the effective rate of sex σ eff in the population as a function of the cost of sex, for different values of the deleterious mutation rate U . Simulation started from purely sexual populations, but additional simulations starting from purely asexual populations (for the same parameter values) led to the same results. Multilocus simulations confirm that the evolutionarily stable rate of sex first drops rapidly as the cost of sex increases from c = 1 (no cost), but then much more slowly as c reaches higher values. This is expected from the analytical results shown in the previous subsections, as the strength of indirect selection for sex increases very sharply below a given value of σ eff (see Figure 3): when the baseline rate of sex is low, indirect selection for sex can be substantial and compensate for strong direct costs. Overall, extrapolations from three-locus QLE models (equation 15, replacing N by N e as described in the previous subsection) provide correct predictions for the equilibrium rate of sex as long as the cost of sex is not too large, the fit being less good for higher values of the deleterious mutation rate U .

This likely stems from the fact that replacing N by N e does not adequately describe the effect of background selection on two and three-locus moments, in particular when the effective rate of sex is low (as was already apparent on Figure 3A): the discrepancy thus becomes more important for higher values of c (since the effective rate of sex σ eff is a decreasing function of c) and for higher values of U (since background selection is stronger when U is high). Note that the differences observed between the analytical predictions and the simulations is not due to the fact that large effect mutations occur at the modifier locus in the simulations: indeed, additional simulations were run in which mutations at the modifier locus had a maximum effect of ±0.03, and gave very similar results (not shown). Finally, Figure 4 also shows that integrating numerically the QLE results over the genetic map (in order to take genetic linkage into account)

gives results which are often similar to those obtained in the case of unlinked loci, unless U is high (compare dashed and solid lines).

As shown by Figure 5, smaller populations have higher rates of sex at equilibrium: this is expected from the analytical model, since all the different mechanisms favouring sex rely on finite population size, and become stronger in smaller populations. Different results were obtained by [START_REF] Iles | Recombination can evolve in large finite populations given selection on sufficient loci[END_REF], [START_REF] Keightley | Interference among deleterious mutations favours sex and recombination in finite populations[END_REF] and [START_REF] Gordo | Sex and deleterious mutations[END_REF], who found that the strength of selection for recombination or sex (measured either by the relative fixation probability or by the increase in frequency of an allele increasing sex/recombination) increases with population size.

However, these studies considered very low baseline rates of sex/recombination, for which the approximations used in the present paper break down. As shown by Fig-

ure 6, simulations indicate that the selection gradient at the modifier locus (measured from simulations representing only two alleles at the modifier locus, as described previously) is lower for N = 20000 than for N = 2000 when σ eff is higher than about 0.005, while the pattern is opposite for lower values of σ eff (in agreement with [START_REF] Iles | Recombination can evolve in large finite populations given selection on sufficient loci[END_REF][START_REF] Keightley | Interference among deleterious mutations favours sex and recombination in finite populations[END_REF][START_REF] Gordo | Sex and deleterious mutations[END_REF]. Understanding why selection for sex is stronger in larger populations when sex is rare would deserve more investigation; although [START_REF] Keightley | Interference among deleterious mutations favours sex and recombination in finite populations[END_REF] explained this result from the fact that larger populations maintain more polymorphic loci, this is not the case in the simulations shown here: the overall genetic variance at selected loci (measured by the sum of p j q j over all loci) is always higher when N = 2000 than when N = 20000, due to the fact that drift allow deleterious alleles to reach higher frequencies (results not shown). Figure 5A also shows that the strength of selection against deleterious alleles (s) has a non-monotonic effect on selection for sex (the equilibrium rate of sex being close to zero when mutations have very strong fitness effects). This effect is captured by the QLE model, and is due to the fact that the different genetic associations (such as D AB ) generating selection for sex are decreasing functions of s (see also

Figure 2); a similar, non-monotonic effect of s on the relative fixation probability of a recombination modifier was observed by [START_REF] Gordo | Sex and deleterious mutations[END_REF]. Overall, the analytical model provides more accurate predictions for larger values of population size, and higher values of s; Figure 7 shows that this remains true for different values of the deleterious mutation rate, up to U = 4 (as can be seen on Figure 7B, the effect of genetic linkage becomes more visible at higher values of U , since deleterious alleles are more abundant within genomes). Figure 8 shows that when s is sufficiently small, the equilibrium rate of sex is only weakly dependent on the values of U and s, as long as U s remains constant: indeed, the equilibrium rate of sex seems to reach a plateau as U increases, decreasing s so that U s remains constant. This effect is not well predicted by the QLE model, which tends to overestimate selection for sex for lower values of s.

According to the analytical model, the number of selected loci should have little effect on the results (for a fixed value of the genomic mutation rate U ), as long as it is sufficiently large so that sums over all loci can be approximated by integrals, and that the mutation rate per locus is small (u s, so that deleterious alleles remain at low frequency). A modified version of the program in which deleterious alleles occur at discrete locations along the chromosome was used to explore the effect of the number of loci: as shown by Figure S1, the results confirm that the equilibrium rate of sex does not depend on the number of loci, as long as it is higher than 10 3 (for the parameter values used in Figure S1). Finally, Figure S2 explores the effect of map length R (average number of cross-overs along the genome at meiosis): R has little effect on the results as long as it is not too small (roughly, R > 5 in Figure S2), while tighter linkage increases selection for sex by amplifying interference effects among selected loci and hitchhiking of modifier alleles increasing sex. As shown by Figure S2, integrating the QLE results over the genetic map overestimates selection for sex when R is small, which again is probably due to the fact that higher-order interactions among loci are not adequately represented by replacing N by N e .

DISCUSSION

Because every population is finite in size, drift-based explanations for the evolution of sex seem particularly compelling (e.g., Otto, 2009); however, few previous papers derived quantitative results on the possible strength of indirect selection for sex generated by this type of mechanism. As we have seen in introduction, the majority of analytical results on the subject are based on recombination modifier models, considering the spread of a mutant allele affecting recombination between two selected loci in haploid organisms. However, the questions of the evolution of sex and the evolution of recombination differ by several important aspects. A first (that was mentioned earlier) concerns the fact that it is probably easier to modify recombination rates than to switch to a different reproductive mode. Although this is certainly true in most obligate sexuals, this type of constraint may be less important in organisms that alternate between reproductive modes and become sexual in response to an environmental stimulus, as in cladocerans or rotifers [START_REF] Hebert | Clonal coexistence in Daphnia pulex (Leydig): another planktonic paradox[END_REF][START_REF] Innes | The origin and genetic basis of obligate parthenogenesis in Daphnia pulex[END_REF][START_REF] Becks | The evolution of sex is favoured during adaptation to new environments[END_REF]. A second difference is that strong intrinsic costs are often associated with sexual reproduction (cost of males, costs generated by the mating process...), while increasing or decreasing recombination should have less direct effect on fitness (although many or too few cross-overs per chromosome may result in aneuploidy). Finally, a third difference is that increasing sex may not yield the same indirect benefits as increasing recombination. While this has been shown previously in the case of diploid organisms, where a single selected locus may generate indirect selection for sex through segregation [START_REF] Uyenoyama | On the origin of meiotic reproduction: a genetic modifier model[END_REF][START_REF] Otto | The advantages of segregation and the evolution of sex[END_REF]Agrawal, 2009a,b;[START_REF] Michod | Deleterious mutations and selection for sex in finite, diploid populations[END_REF], we have seen in this paper that this is also true in haploids, as a single locus under selection tends to favour any mutation increasing sex in a finite population. This effect (which does not occur in the case of a recombination modifier) stems from the fact that alleles coding for more frequent sex escape more easily from low-fitness genetic backgrounds than alleles coding for lower rates of sex; the same type of benefit is also produced by interactions between selected loci, in addition to the benefits gained by breaking negative associations between those loci. This mechanism (which had not been described before) cannot provide an advantage to sexuals when competing against pure asexuals (since asexuals never "share" their genetic background), but is susceptible to select for increased rates of occurrence of sex in organisms using both reproductive modes. It is interesting to note that this mechanism presents some similarities with the concept of condition-dependent sex or recombination [START_REF] Agrawal | The evolution of plastic recombination[END_REF]Hadany andOtto, 2007, 2009), in which alleles that increase the rate of sex or recombination when present in low-fitness individuals benefit from a deterministic advantage due to the possibility of escaping from these low-fitness genotypes; however, an important difference is that the effect described in the present paper does not require condition-dependent rates of sex.

Can indirect selection generated by stochastic forces compensate for strong direct costs of sexual reproduction? This requires that sufficiently many selected loci are polymorphic at the same time, as the effect of a single selected locus (or a pair of interacting loci) is rather weak. As we have seen, predictions from analytical models can be extrapolated to situations where selection occurs at many loci (under certain conditions), but this neglects the effect of higher-order interactions (between three or more selected loci) which are difficult to explore analytically due to the large number of variables involved. Although the effect of each of these higher-order interactions should be relatively weak, the number of such possible interactions becomes very large rapidly as the number of selected loci increases. The multilocus simulation results presented in this paper show that extrapolations from three-locus models may nevertheless provide correct predictions for the evolutionarily stable rate of sex in some cases, once the overall effect of selection occurring in the background is accounted for by replacing population size by the effective population size N e , derived from background selection models. This suggests that selection on a sex modifier mainly stems from its effect on pairwise associations between selected loci, at least as long as the baseline rate of sex is not too low; however, effect of higher-order interactions should be explored more rigorously in order to obtain a better understanding of the different effects involved Overall, the results presented in this paper show that deleterious mutations tend to favour small rates of sex in the presence of strong direct costs of sex, but are unlikely to explain the evolution of obligate sex (again, due to the fact that indirect selection decreases very rapidly as the baseline rate of sex in the population increases).

Note however that population spatial structure should enhance indirect selection due to stochastic effects (while also slowing the spread of asexual mutants) and therefore allow higher rates of sex to be maintained [START_REF] Peck | The maintenance of sexual reproduction in a structured population[END_REF][START_REF] Martin | Selection for recombination in structured populations[END_REF][START_REF] Salathé | Mutation accumulation in space and the maintenance of sexual reproduction[END_REF]Hartfield et al., 2012). Generally, the magnitude of indirect se- in terms of M juv t ). Note that allele frequencies may change during each of these steps (due to differential fecundities, the direct cost of sex, and drift). The recursions describing the effects of these three events (selection, reproduction, drift) are in turn divided into two steps. A first step corresponds to the calculation of the expectation of moment M after the event, but using the values of allele frequencies before the event: in particular, M t corresponds to the expected value of M after selection, but where the allele frequencies that appear in M (note that allele frequencies also appear in the definitions of LDs) are replaced by their values before selection. For example, p A D MA t is the expected value of

p A E sel [(X M -p M ) (X A -p A )]
, where E sel stands for the average over individuals after selection (weighting each individual by its relative fecundity), and p M , p A are the frequencies of M and A before selection.

Similarly, M t and M t correspond to the expected values of M after reproduction and drift (respectively), but using the values of allele frequencies before reproduction (for M t ) and before drift (for M t ). Note that We have:

D M t = E sel [X M -p M ] t
p A D MA sel t = (p A + ∆ s p A ) E sel [(X M -p M -∆ s p M ) (X A -p A -∆ s p A )] t (A1)
where p M and p A are allele frequencies before selection. Expanding in terms of ∆ s p M , ∆ s p A (written as D M , D A ) yields:

p A D MA sel t = p A D MA t + D MA D A t -p A D M D A t -D M D A 2 t . (
D U = E W W ζ U (A3)
where U is any set of loci, ζ U = i∈U ζ i , W is the fecundity of an individual and W the average fecundity in the population. These can be written as

W = (1 + s p A + s ζ A ) (1 + s p B + s ζ B ) W = (1 + s p A ) (1 + s p B ) + s 2 D AB (A4)
and the ratio W/W can be expanded in a Taylor series in s, allowing one to compute D U in terms of allele frequencies and linkage disequilibria before selection. For example,

W/W = 1 + s (ζ A + ζ B ) to the first order in s, yielding p A D MA t = p A D MA t + s p A D MAA t + p A D MAB t + o (s)
. Linkage disequilibria with repeated indices (such as D MAA ) can be rewritten using the relation:

D Uii = p i q i D U + (1 -2p i ) D Ui (e.g., eq.
5 in [START_REF] Barton | Natural and sexual selection on many loci[END_REF] To compute the effect of reproduction on linkage disequilibria, I call σ ind the rate of sex of a given individual (equal to σ or σ + δσ depending on its genotype at the modifier locus) and σ the average rate of sex in the population. We then have: given by symmetrical expressions. The moments D MA D AB , D MB D AB , p A D MAB and p B D MAB are of order s δσ:

D U = E sel c (1 -σ ind ) c (1 -σ) + σ ζ U + ST=U r S,T E sel σ ind c (1 -σ) + σ ζ S E sel σ ind σ ζ T . (
D MA D AB t = (1 -ρ MA ) s (1 -ρ AB ) D MAB D AB t -δσ eff (r AB + c -1) D MAB D MA sel t 1 -(1 -ρ MA ) (1 -ρ AB ) (B6) p A D MAB t = 1 ρ MAB -δσ eff [r AB + (c -1) ρ AB ] p M q M p A D AB sel t + (1 -ρ MAB ) s D MAB D AB t -δσ eff (c -1) D MAB D MA sel t (
p M q M D AB sel t = -2s p M q M p A D AB t + p M q M p B D AB t + s p M q M D AB 2 t ρ AB (B8) D MA D MB sel t = s D MAB D MA t + D MAB D MB t + s D MAB 2 t 1 -(1 -ρ MA ) (1 -ρ MB ) (B9) p M D MAB sel t = s D MAB D MA t + D MAB D MB t + s D MAB 2 t ρ MAB (B10)
while the same moments before selection are obtained by multiplying the previous expressions by 1 -ρ AB , (1 -ρ MA ) (1 -ρ MB ) and 1 -ρ MAB (respectively). The threelocus linkage disequilibrium is generated by selection, drift and the modifier effect, and is of order s 2 :

D MAB t = 1 ρ MAB -δσ eff [r AB + (c -1) ρ AB ] p M q M D AB sel t + δσ eff [r MAB + (c -1) (2 -ρ MAB )] p M D MAB sel t + δσ eff [r AB + (c -1) (2 -ρ MA -ρ MB )] D MA D MB sel t -(1 -ρ MAB ) 2s D MA D AB t + D MB D AB t + p A D MAB t + p B D MAB t + 4s 2 D MAB D AB t . (B11)
Finally, the two-locus moments involving the modifier and a single selected locus are also affected by selection at the third locus. One obtains:

D MA 2 sel t = 1 ρ MA (2 -ρ MA ) pq MA t N + s 2 D MAB D MA t + s D MAB 2 t (B12) D MA 2
t (before selection) being obtained by multiplying the term in s by (1 -ρ MA ) 2 in the previous expression,

p A D MA t = 1 -ρ MA ρ MA s ( p A D MAB t + D MA D AB t ) + s 2 D MAB D AB t -δσ eff (c -1) D MA 2 sel t (B13) p M D MA sel t = 1 ρ MA s D MA 2 t + D MA D MB t + p M D MAB t + s 2 (2 D MAB D MA t + D MAB D MB t ) + s 3 D MAB 2 t (B14)
(which has to be multiplied by 1 -ρ MA to obtain the same moment before selection), and From equations C1 to C17, the selection gradient can be written under the form:

D MA t = 1 ρ MA δσ eff [(c -1) (2 -ρ MA ) + r MA ] p M D MA sel t + (1 -ρ MA ) s ( D MAB t -2 p A D MA t ) -2s 3 D MAB D AB t -s 2 2 p A D MAB t + 2 D MA D AB t + p B D MAB t + D MB D AB t . ( B15
s M ≈ δσ eff -(c -1) + s 2 N Θ 1 (p A + p B ) + s 4 N Θ 2 p A p B (C18)
where Θ 1 and Θ 2 are functions of recombination rates, the effective rate of sex σ eff , the cost of sex c and strength of selection s. Assuming free recombination between all loci, one obtains: Θ 1 = 2 (σ eff + 2s) 3 1 + 8 (c -1) (σ eff + 3s) σ eff + 4s (C19) and Θ 2 = 32 [(5σ eff + 16s) X 1 + 8 (c -1) (σ eff + 6s) X 2 ] 3 (σ eff + 2s) 3 (σ eff + 3s) (σ eff + 4s) 2 (σ eff + 6s) (σ eff + 8s) 3 (5σ eff + 12s) (5σ eff + 16s) (C20) with: (C22)

X 1 =1133568
Expressions for arbitrary recombination rates between the three loci can be obtained from equations C1 to C17. Extrapolating these results to the case where deleterious alleles occur over a whole genome can be done by summing equation C18 over all loci A and B, which yields equation 15 in the main text.

  (1 + s): this corresponds to the term in D MAB t on the second line of equation 9.The second effect (which is often much stronger) is the fact that the larger variance in fitness present in the sub-population carrying M increases the efficiency of selection: the fitter alleles A and B tend to increase in frequency within this sub-population (through the increase of the AB genotype), generating positive pairwise associations D MA , D MB (equation B15). Allele M thus tends to be associated with the fitter alleles A and B, and therefore increases in frequency (term in D MA t + D MB t on the first line of equation 9). Where does the positive D MAB come from? As shown by equation B11, different mechanisms are involved, of which I will only discuss two. A first source of D MAB is the negative linkage disequilibrium D AB generated by the Hill-Robertson effect (first line of eq. B11). This corresponds to the often described benefit of sex/recombination in a finite population: increasing recombination tends to break negative associations between selected loci, and thereby increases the efficiency of selection. A different source of positive D MAB (second line of equation B11), however, does not involve recombination between the selected loci. This second mechanism stems from a positive correlation between p M and D MAB , which can be explained as follows. Finite population size generates random fluctuations in D MAB ; when D MAB > 0, extreme genotypes (ab, AB) are more frequent on the M background and AB increases in frequency, generating positive associations D MA , D MB (by contrast, when D MAB < 0, AB increases in frequency on the m background, generating negative D MA , D MB ): drift and selection thus generate a positive correlation between D MAB and D MA (eq. B4), and also between D MAB and D MB . When D MA , D MB are positive, M tends to increase in frequency by hitchhiking, which in turn produces a positive correlation between D MAB and p M (eq. B10). This positive correlation generates a positive D MAB , in the same way as the positive correlation between p M and D MA generated a positive D MA in the two-locus model: during sex, alleles at the modifier locus tend to move more often from Ab, aB backgrounds to ab, AB backgrounds than the opposite, and allele M (which engages in sex more often) benefits more from this effect. Note that this last source of positive D MAB would operate even in the absence of recombination between the selected loci, and would not appear in the case of a recombination modifier. Indeed, using the same method to derive an expression for the expected change in frequency of a neutral recombination modifier (changing r AB to r AB + δr AB ) yields the same expressions as in Appendix B for c = 1 (no direct selection) and replacing δσ eff r AB by δr AB , except that the terms involving p M D MAB sel t and p M D MA sel t in equations B11 and B15 disappear.

Figure

  Figure 3A compares the prediction from equation 15 with multilocus simula-

(

  and more accurate predictions for low rates of sex). Multilocus simulations also confirm several qualitative predictions from the analytical model, in particular the fact that indirect selection increases very sharply as the baseline rate of sex in the population decreases (causing a relatively slow decrease of the equilibrium rate of sex as the cost of sex increases), and the fact that the strength of indirect selection is maximized for intermediate values of s, and is a decreasing function of population size.While the present model only considers haploid individuals, QLE results for diploids were derived in Roze and Michod (2010) using two and three-locus models (equivalent to Appendix B in the present paper). Analytical and multilocus simulation results indicate that the different parameters (N , c, U , s) have similar qualitative effects in haploids and diploids in the absence of dominance (additive deleterious alleles: h = 0.5). When deleterious alleles are partially recessive, however, sex is often disfavoured due to segregation load (generated by dominance within loci and associative overdominance between loci), and multilocus simulations indicate that populations may evolve towards complete asexuality when sex is costly (see Figures6 and 7in Roze and[START_REF] Michod | Deleterious mutations and selection for sex in finite, diploid populations[END_REF]. In general, conditions favouring sex thus seem more restricted in diploids than in haploids.

  A5)The first and second terms of equation A5 correspond to the effects of asexual and sexual reproduction on D U ; the second term involves a sum over all possible partitions of the set U into two subsets S and T, r S,T being the probability that among the loci in the set U present in a haploid juvenile (after meiosis), loci in the subset S are inherited from the maternal gamete, and loci in the subset T from the paternal gamete. To compute D U as a function of linkage disequilibria and allele frequencies after selection, one can replace σ ind and σ by σ + δσ p M + δσ ζ M and σ + δσ p M in equation A5, and compute all fractions to the first order in δσ.

Figure 1 .

 1 Figure 1. Average increase in frequency of allele M during the sweep of the beneficial

Figure 2 .Figure 4 .

 24 Figure 2. Average value of the linkage disequilibrium between two deleterious alleles

Figure 5 .Figure 7 .

 57 Figure 5. Effective rate of sex in the population at equilibrium as a function of the

  is the expected change in frequency of M during selection (also denoted ∆ s p M t ), while D M t is the expected change in frequency of M during reproduction (also denoted ∆ r p M t ). Although D M t (the expected change in frequency of M due to drift) equals zero, higher-order moments (such as D M 2 t ) are different from zero. Finally, the second step of the recursions corresponds to the update of allele frequencies, and is identical for the three life cycle events. Although this second step can be represented by a general expression (e.g., eq. A43 in Roze and[START_REF] Michod | Deleterious mutations and selection for sex in finite, diploid populations[END_REF], here I simply show the derivation of p A D MA sel t as an example (which is easily generalized to other moments).

  Similarly, the calculation of M t in terms of M t (selection), M t in terms of M sel t (reproduction) and M t in terms of M juv The general principles are the following. To compute moments M t in terms of moments before selection M t , one uses the fact that D i , D ij and D ijk after selection (but using allele frequencies before selection) are given by:

A2) p A D MA juv t and p A D MA t+1 are given by the same expression, replacing moments M t by M t and M t , respectively. A general expression that allows one to compute the effect of updating allele frequencies on any moment can be implemented in Mathematica (see Supplementary Material). t (drift) can be automated in Mathematica.

  . For example, p A D MAA t = p A D MA t -2 p A 2 D MA t , since any moment M t , M sel

t or M juv t involving D i (where i is a single locus) equals zero.

  B7) while D MB D AB t and p B D MAB t are given by symmetrical expressions. Moments p M q M D AB , D MA D MB and p M D MAB are generated by selection and drift, and are of order s 2 . After selection, they are given by:

  ) B D MB t , p M D MB t and D MB t are given by symmetrical expressions.D MA D AB t ≈ -s D MAB D AB t + δσ eff (c -1 + r AB ) D MAB D MA t ρ MA + ρ AB + 3s (C12) p A D MAB t ≈ -s D MAB D AB t + δσ eff [r AB pq M p A D AB t + (c -1) D MAB D MA t ] ρ MAB + 3s MAB + 2s -δσ eff r AB pq M D AB t + δσ eff [2 (c -1) + r MAB ] p M D MAB t + δσ eff [2 (c -1) + r AB ] D MA D MB t + 2s ( D MA D AB t + D MB D AB t + p A D MAB t + p B D MAB t ) (C14) p A D MA t ≈ -s ( D MA D AB t + p A D MAB t ) + δσ eff (c -1) D MA D MA t ≈ s (2 p A D MA t -D MAB t ) + δσ eff [2 (c -1) + r MA ] p M D MA t ρ MA + s . (C16)Under the same assumptions, the selection gradient at the modifier locus is given by:

	Moments D MB D MAB t ≈ ρ 2 2 1 ρ MA + 2s t , p (C13) t (C15) s M ≡ ∆p M t p M q M t ≈ -δσ eff (c -1) -s D MA t + D MB t pq M t . (C17)

  s 6 + 2412288 s 5 σ eff + 2142144 s 4 σ 2 eff + 1016320 s 3 σ 3 X 2 =921600 s 6 + 1920000 s 5 σ eff + 1660672 s 4 σ 2 eff + 763280 s 3 σ 3

	eff
	(C21)
	+ 271780 s 2 σ 4 eff + 38848 s σ 5 eff + 2319 σ 6 eff
	eff
	+ 196622 s 2 σ 4 eff + 26917 s σ 5 eff + 1530 σ 6 eff .

Finally, the effect of drift on moments involving linkage disequilibria is obtained using the following reasoning. D U t is the expected value of ζ U = i∈U (X i -p i ) (where p i 's are allele frequencies before drift, among juveniles) in an individual sampled randomly in the population after drift. Because drift corresponds to the random sampling of N juveniles, this is the same as the average value of ζ U among juveniles;

). This method is readily extended to more complex moments (e.g., [START_REF] Michod | Deleterious mutations and selection for sex in finite, diploid populations[END_REF].

Throughout the paper, I assume that population size is large and calculate all recursions to the first order in 1/N .

APPENDIX B: EXPRESSIONS FOR MOMENTS IN THE THREE-LOCUS

MODEL

The three-locus model involves different moments of allele frequencies and LD, which are of different order in s and δσ. In the following, r MAB is the probability that at least one recombination event occurs between the three loci, while ρ MAB = r MAB σ eff .

All moments are expressed at QLE to the first order in 1/N and δσ, and to leading order in s; derivations can be found in the Supplementary Material. The variance in the three-locus linkage disequilibrium D MAB is generated by finite population size:

where pq MAB t = p M q M p A q A p B q B t . Similarly, we have:

The moment D MAB D AB is generated by the modifier effect and drift:

with

tion and drift, and are of order s; when measured after selection, they are given by:

Moments before selection are given by the same expressions multiplied by (1 -ρ MA ) ×

(1 -ρ MAB ) and by 1 -ρ AB (respectively), while D MAB D MB and p M q M p B D AB are APPENDIX C: QLE UNDER WEAK RATE OF SEX

The following expressions for moments of allele frequencies and linkage disequilibria are derived for the case where alleles A and B are deleterious and maintained at mutation-balance (at equilibrium frequency p A = p B ≈ u/s, where u is the per-locus deleterious mutation rate), and for weak effective recombination rates ρ X (of order s); derivations can be found in the Supplementary Material. Moments generated by selection and drift are given by:

Moments generated by the modifier effect are given by: