Effective population size and evolutionary dynamics in outbred laboratory populations of Drosophila

Census population size, sex-ratio and female reproductive success were monitored in 10 laboratory populations of Drosophila melanogaster selected for different ages of reproduction. With this demographic information, we estimated eigenvalue, variance and probability of allele loss effective population sizes. We conclude that estimates of effective size based on gene-frequency change at a few loci are biased downwards. We analysed the relative roles of selection and genetic drift in maintaining genetic variation in laboratory populations of Drosophila. We suggest that rare, favourable genetic variants in our laboratory populations have a high chance of being lost if their fitness effect is weak, e.g. 1% or less. However, if the fitness effect of this variation is 10% or greater, these rare variants are likely to increase to high frequency. The demographic information developed in this study suggests that some of our laboratory populations harbour more genetic variation than expected. One explanation for this finding is that part of the genetic variation in these outbred laboratory Drosophila populations may be maintained by some form of balancing selection. We suggest that, unlike bacteria, medium-term adaptation of laboratory populations of fruit flies is not primarily driven by new mutations, but rather by changes in the frequency of preexisting alleles.


Introduction
A comprehensive understanding of evolution in either laboratory or natural populations requires analysis of the impact of both natural selection and random genetic drift on evolutionary trajectories. Moreover, the interaction of genetic drift with natural selection is a topic of much interest in both evolutionary biology and conservation (Hartl and Clark 1989;Futuyma 1998;Allendorf and Luikart 2007). Experimental evolution studies in the laboratory can create conditions that lead to very strong natural selection under controlled environmental conditions, often in moderate sized populations, yielding systems particularly conducive to the study of how drift and selection interact in the course of adaptive evolution (e.g. Latter and Mulley 1995;Rose et al. 1996;Matos et al. 2002Matos et al. , 2004Simões et al. 2007;Garland and Rose 2009). However, there has been criticism of the use of laboratory populations, especially for the study of adaptive evolution (Promislow and Tatar 1998;Harshman and Hoffmann * For correspondence. E-mail: ldmuelle@uci.edu. 2000). These critiques have claimed that the methods used for maintaining these laboratory populations will permit the increase of deleterious mutations that ultimately compromise the value of these populations for the study of adaptation by natural selection. However, a rigorous assessment of these claims would require a detailed assessment of the force of random genetic drift in such populations, which can only be accomplished with some estimate of their effective population sizes.
There is often only a rough idea of the effective population size (N e ) in laboratory evolution studies, with N e usually being estimated by census population size. A variety of methods have been devised for estimating the effective population size of natural populations (Krimbas and Tsakas 1971;Nei and Tajima 1981;Pollak 1983;Waples 1989Waples , 2006Jorde and Ryman 1995;Wang et al. 2010). These methods typically rely on fluctuations in allele frequencies over time, based on the assumption that the loci screened are effectively neutral. Under this assumption, allele frequency variation should be determined by genetic drift alone, although even modest levels of some types of selection can lead to estimates of N e that are lower than the true value (Mueller et al. 1985). Indeed, these methods often arrive at estimated values of N e that are an order of magnitude less than the census population size (Krimbas and Tsakas 1971;Begon et al. 1980). Alternatively, it is possible to estimate N e without reference to genetic changes in populations, by making use of detailed demographic data (cf. Begon et al. 1980;Mueller et al. 1985). However, the collection of this type of demographic information is often difficult or impossible for many natural populations. Laboratory populations, by contrast, are typically amenable to the application of these techniques.
We have long used replicated laboratory populations of Drosophila melanogaster to study the experimental evolution of senescence and other functional characters. Some of these populations have now been under continuous laboratory selection for over 800 generations (vid. Rose et al. 2004). To gain insight into the expected magnitude of genetic drift in these populations, we undertook a multi-year study of population-size fluctuation and within-population variance in female reproductive success. Here, we use those data to evaluate alternative methods for estimating effective population size, as well as to make inferences concerning the dynamics of selected genes in these populations. This information can then be used to assess more directly the concerns raised about the utility of using laboratory populations to study adaptive evolution.
A major goal of the present study is to investigate the impact of genetic drift on the evolution of Mendelian laboratory populations. This goal is quite different from that of understanding the role of drift in natural populations, since we are concerned with relatively short periods of time and small populations. Further, because of the short periods of time and the small populations involved, new mutations are systematically less likely to play an important role in the adaptation of laboratory populations of Drosophila compared to the situation in bacterial laboratory evolution (cf. Lenski et al. 1991).
The genetic response of Drosophila laboratory populations to age-specific selection happens rapidly. For instance after only 15 generations of selection for high fitness in latelife, females showed a 33% increase in longevity and males a 15% increase (Rose 1984). These differences in longevity between selected and control populations have continued to increase over time until they are now nearly threefold. With breeding populations of about 1000 individuals, this evolution is too fast to be due to the de novo appearance of beneficial mutations. Likewise, the relatively similar response to selection in the five independent populations studied in this case further argues against the idea that new mutations have played an important role in the adaptation of these populations. Against this backdrop, we have focussed on the study of genetic mechanisms we believe to be the most important for the evolution of these populations, specifically the increase in frequency of extant alleles from the start of these experiments and the role of drift in slowing or impeding this evolution.

Populations
Ten populations of D. melanogaster were used to obtain information on population size, sex-ratio and female fecundity. Five of these populations, the B's, have been maintained on a two-week discrete generation cycle entirely in bottles or vials since February 1980 (Rose 1984;Rose et al. 2004). Each of these populations has been maintained for the past 33 years by starting each generation with eggs laid in 20 vials by about 80 adults each. Fourteen days after egg collection, the emerged adults in each vial are transferred to a fresh 8-dram (22 × 95 mm) vial and given about 1 h to lay eggs, with a target density of 60-100 eggs per vial. All B populations (B1. . .B5) have been kept genetically independent of one another, with no crossing. Genomic analysis indicates that these populations have indeed undergone independent evolution (A. D. Long, personal communication).
The second set of five populations, the O's, was under selection in our lab for late-life reproduction (Rose 1984;Rose et al. 2004). Fourteen days after egg collection, emerged adult O flies were transferred from up to 60 vials into 2-5 plexiglass cages (20.2 × 15.5 × 24.8 cm) per population, resulting in an initial number of young adults that has varied between about 2000 and 5000 per O population. Each cage is provided with an unyeasted food plate for food and moisture, at intervals of 2-3 days, for the next 7-9 weeks. Two to three days before egg collection, a yeasted food plate was placed in the cage, followed by a fresh unyeasted plate to collect eggs over a 24 h period. These eggs were transferred to 8-dram vials at a density of about 80 eggs per vial to initiate the next generation.
All the 10, B and O populations were maintained on banana-barley malt-corn syrup food, at about 25 • C and under constant light.

Census data
From June 1992 to June 1997, the total number of flies in each B population was sexed and counted immediately after egg collection. A B-census was taken approximately every other generation during this five year period. From July 1992 to August 1994, the total number of flies in each O population was sexed and counted immediately after egg collection. An O census was taken approximately every generation during this two year period. also taken from the O populations just prior to their normal days of reproduction after they had been exposed to yeasted food plates for 48 h, and placed in unyeasted vials for 24 h and the number of eggs in the vials subsequently counted. The normal O-culture regime involves egg laying in cages. Thus, the measurements made here in vials and on individual flies may exhibit genotype-by-environment effects that are different from female egg laying behaviour in cages. Historically, if the O flies had not laid a sufficient number of eggs in 24 h, the egg laying period would be extended. We expect that this resulted in fewer females laying no eggs and this would, to some extent, increase N e estimates relative to the values obtained from 24 h egg laying observation, although we have not attempted to quantify this. These samples were collected approximately every generation between January 1993 and April 1994.

Population size and fecundity
The sex ratio in the B and O populations is close to 1 ( total number of O flies per population also shows a gradual increase from the middle of 1993 to the end of the experiment. This suggests that we must treat each sample date as a random variable and sample evenly across the different census dates. Fecundity of the B and O females shows a similar pattern of variation between census dates (figure 3). The pattern is more pronounced in the O populations. For instance, we see that in the first and last samples of 1993 fecundity is uniformly low in the O populations, while fecundity is much higher in all O populations in the June (third) 1993 samples.

Estimates of effective population size
The standard theory for studying random genetic drift is the Wright-Fisher model (Ewens 1979, pp. 75). This model makes a number of useful predictions about the impact of population size (N) on genetic parameters like the variance in allele frequencies over time and the rate at which genetic variation (heterozygosity) is lost. However, many real populations do not conform to the assumptions of the Wright-Fisher model due to unequal numbers of each sex, or total numbers that vary over time, or variation in reproductive success, or selection at linked loci (Gillespie 2000). Population geneticists have addressed these problems by deriving a modified population size, called the effective population size (N e ), which behaves in a similar fashion to N in the Wright-Fisher model. In some cases, simple analytic formulas can be derived for N e . However, with the varying and complicated demography of the B and O populations studied here, that is not possible. An alternative approach, which we will use, is to simulate genetic change in a population that follows the demography of a B or O population and then use the dynamics of neutral genetic variation in these simulations to infer the effective population size.
In our experimental populations, we expect the effective population size to be impacted by sex ratio, varying population size and variation in reproductive success. However, we do not expect the effective population size to be affected by repeated sweeps of new, favourable mutants as described by Gillespie (1994aGillespie ( , b, 1997Gillespie ( , 2001 who noted that selected substitutions at one locus can induce random fluctuations at closely linked neutral loci that resemble drift even in populations of infinite size. These allele frequency dynamics, called pseudohitchhiking, require the appearance of favourable mutants that then drag along linked variation at nearby loci. This effect can lower the equilibrium levels of variation, giving the appearance that the population is smaller than traditional theory would suggest. In our laboratory populations, the number of individuals and generations is simply too small to expect that this phenomena will have been of any significance over the period of these experiments. As an example, consider a B population with N e = 646, 780 generations of laboratory evolution and a mutation rate at a locus to a favourable allele of 5 × 10 −7 (cf. Gillespie 2001). Then, over the entire 780 generations, less than one favourable mutant (0.50) would be expected and once the mutant appears, it's time to fixation would be many hundreds of generations. The chances of this type of hitchhiking would be even less in the O populations, due to the much smaller number of generations they have undergone.

Drift simulation:
The sex ratio, total population size and variance in female fecundity change over time in our data, making a simple estimate of effective population size difficult. Consequently, we have chosen to use the observations of demographic variation in combination with computer simulations to make estimates of effective population size. We also used computer simulations to examine the effect of selection on the dynamics of allele frequencies at linked neutral loci.
We carried out simulations of B and O demography to cover roughly the period of selection from 1980 until the end of 2010. For the O populations, this entails 160 generations of selection and for the B populations 780 generations. We initially started each population in our simulations with 500 males and 500 females. Genotypes of these individuals were determined by the starting allele frequencies which varied with each simulation. Simulations involved (i) a single locus; (ii) two linked loci and (iii) three linked loci. The single locus simulations incorporated the effects of drift only and were used exclusively to estimate N e in B and O populations. The two locus models studied drift and selection at both loci. The three locus simulations were similar to the two locus simulations except that the middle locus was assumed to be neutral and selection was based on the genotype at the two flanking loci. The two and three locus simulations were used to study the loss of rare genetic variants favoured by natural selection and the effects of selection on tightly linked neutral loci, respectively.
We let the numbers of males and females in the current generation, t, be N t,m and N t,f , respectively. For each of these individuals, we also kept track of their individual genotypes. We do not have specific data for the variation in male reproductive success among these populations during the period under study, but other experimental work with B males suggest that the distribution of male mating success is approximately Poisson in these populations (Joshi et al. 1999) with the majority of males mating only once. However, Long et al. (2007) have shown that in a B environment, many females mate more than once prior to day 14, the standard day of egg laying. This could present an opportunity for eggs to be fertilized by more than one male parent to the extent that the last male mate does not completely displace sperm from earlier matings. We know even less about the mating of O flies. Consequently, we modelled two different mating scenarios. In the 'sampling zygote' scenario we have assumed that each female in the B and O populations mates only once with a male chosen at random and then the next generation is created by sampling female zygotes.
Thus, a male mate for each female was chosen from the N t,m total with replacement (figure 4, sampling zygotes). While this gives rise to a binomial distribution in the number of matings, since the probability of each male being successful is small (1/N t,m ) this distribution will be approximately Poisson.
In the second scenario, called 'sampling gametes', we assume that at the time of egg laying, a typical female may have sperm from several males and thus their eggs may not be fertilized by only a single male. To accommodate this life-history, we have also done simulations by sampling gametes from males and females separately (figure 4, sampling gametes).
Our observations of female fecundity, however, show significant variation both between individual females and over the different census periods. Thus, the relative probability of each female reproducing was determined in two steps. First, a census date was chosen at random. We then sampled N t,f fecundity values from the observed fecundities at the chosen census period. These fecundities were normalized so that the ith female (f i ) had a probability, p t,fi , of contributing gametes with the requirement that 1 = N t,f i=1 p t,fi . When we sampled zygotes, the mate for each female had been chosen, and we would eventually sample eggs at random from the N t,f females with the probability of any particular female's fertilized egg being chosen equal to p t,fi . However, if a particular female was chosen twice, the zygotes would be fertilized by the same male. When we sampled gametes we took a second male gamete sample if the same female was chosen more than once.
The numbers of males and females in subsequent generations (N t+1,m and N t+1,f , respectively) were chosen at random from all the observed male-female pairs summarized in figure 2. The genotypes for each of these new individuals were determined by sampling a female parent from the current female population, utilizing the distribution p t,f . The genotype of the fertilized egg from this female was determined by random sampling of alleles from that female and her previously determined mate. In simulations with selection, these probabilities were further modified based on the fitness of the individual's genotype. Thus, a drift-only simulation following the life-cycles in figure 4 departs from a standard Wright-Fisher model in four important ways: (i) sex ratio departs from 1:1, (ii) the population size varies over time, (iii) the variance in female reproductive success is not Poisson and (iv) in the case of sampling zygotes, the females are monogamous. Most importantly, we use our empirical observations from the B and O populations to model sex ratio, population size variation, and female reproductive success.
In the case of the single-locus simulations, the chance of the maternal allele or paternal allele in a parent being chosen for the next generation is simply 1/2. For the two-locus and three-locus simulations, there were multiple gametes possible and their probabilities depend on the recombination fraction, r. With two loci, there were two nonrecombinant gametes that each has a probability of (1 − r)/2 of transmission. There were likewise two recombinant gametes each with probability r/2 of occurrence. For three loci, we assumed the gene order is: A, B, C. We further assumed that the recombination fraction between loci A and B is r 1 and between loci B and C, r 2 . We assumed that there is no interference, such that two nonrecombinant gametes occur with There are assumed to be two gametes that result from a recombination between loci A and B, with frequency r 1 (1 − r 2 ). Symmetrically, there are also two gametes that result from a recombination between loci B and C with probabilities r 2 (1 − r 1 ). Finally, we have two double recombinant gametes that occur with frequencies r 1 r 2 . After drawing zygotes for each of the N t+1,m males and N t+1,f females, the new allele (gamete) frequencies were calculated and the cycle was repeated.
Estimating N e : The main idea behind estimates of effective population size is the use of results from the simple Wright-Fisher model of drift, with analogous results derived for populations with more complicated demography (Ewens 1982). For instance, under the Wright-Fisher model, a population of size N, with two alleles at a single locus, with frequencies p and 1 − p, respectively, has an expected variance in the allele frequency in the next generation of p(1 − p)/2N. For a population that does not conform to the Wright-Fisher model, we can use this result to define the variance effective population size N (v) e , as In the single locus, drift-only simulations we produced 1000 single generation transitions from a variety of starting allele frequencies.
These 1000 values can be used to estimate Var(p) and estimate N (v) e using equation (1) (table 1). The estimates obtained vary depending on the initial allele frequency and sampling scheme but the estimates of variance effective population size based on sampling zygotes average around 348 for the O populations and 363 for the B populations (table 1). However, N (v) e estimates based on sampling gametes are 698 for the O population and 787 for the B population.
The leading eigenvalue of the Wright-Fisher transition matrix can also be used to predict the loss of heterozygosity over time. In a population with initial heterozygosity, H 0 , the expected heterozygosity at time t, H t , is 1 − 1 The estimated eigenvalue effective sizes also differ depending on the initial allele frequency ( The variance effective size estimates for the B populations obtained here are similar to those estimated by Teotónio et al. (2009), which were about 400. Teotónio et al. (2009) used changes in the allele frequencies over a defined period of time to estimate the effective population size, following the methods described by Nei and Tajima (1981). Teotónio et al. (2009) used only a single locus to obtain their estimates. In simulations of the Nei and Tajima method (results not shown), we have found that the results are sensitive to the number of loci used in the simulation and the starting allele frequency. Accordingly, it is unclear how much weight should be placed on the observation that the Teotónio et al. (2009) estimates are similar to our variance effective size estimates.

Probability of allele loss
From the single-locus drift results, we can also estimate the probability that an allele will be lost during the course of observed B and O evolution. We do not base our estimates of loss probabilities on any estimated effective population size. Rather these estimates of loss are based on the iteration of the life cycle outlined in figure 4 (sampling zygotes) for either 160 generations (O populations) or 780 generations (B populations). This simulated evolution was repeated 1000 times and the number of iterations that resulted in the loss of an allele was used to empirically estimate the probability of allele loss. These estimates ( figure 5) show that, as expected, the chance of rare alleles being lost is greater in the B populations than the O populations. Although the effective population size of the B populations is about the same as the O's, the B populations have undergone many more generations of drift than the O's and that difference in generation number explains this result. Diffusion approximations can be used to estimate the probability of allele loss after a fixed period of time given some starting allele frequency (Ewens 1979, eqn. 5.29). If an allele initially at frequency p, undergoes m generations of drift then the probability that the allele is lost (goes to frequency 0) is given by, where C = 1 4 [arcos (1 − 2p)] 2 and N (l) e is the effective population size. With the results in figure 5, we can estimate a 'probability-of-allele-loss' effective size as the value of N (l) e that best predicts the values of P 0 . Using eqn. (3), we found the least squares estimate of N (l) e for the O populations as 518 (±53, 95% confidence interval) and for the B populations 599 (±24) which are in good agreement with the range of eigenvalue effective size estimates (table 1).
Our simulations have sampled the numbers of male and female breeding parents independently each generation from all observed population sizes. However, as the results in figure 2 illustrate, there appears to be a slow increase in population size over the course of the study, especially in the O populations. In principle, this temporal population size increase could result in a different rate of allele loss relative to a random sequence of population sizes. To study this in the O populations, we estimated the parameters a 0 and a 1 from the linear equation, N t = a 0 + a 1 t + ε t , using the observations summarized in figure 2 separately for males and females. In our simulations, this regression was used between generations 65 and 76 to generate male and female breeding populations and we estimated the Var(ε t ) = σ 2 N from the residual sum of squares. Since the variation around the regression is both asymmetric and scales with the mean, we simulated random population sizes around the regression with a lognormal distribution. Thus, at generation t, we generated a random variable, The population size used in the simulation was round (exp (x t )), where the function round() rounds a number to the nearest integer. From these simulations, which sampled gametes at random, we estimated N (l) e as 960 ± 82, which is not significantly different from our previous estimate of the eigenvalue effective size (955).
We have noted previously (e.g. Borash et al. 2007) that genes that affect survival or fecundity after the first week of adult life are essentially neutral in the B populations. Borash et al. (2007) therefore looked for inbreeding depression for late-life fitness in the B populations. There was no inbreeding depression observed for either longevity or female age-specific fecundity. However, there was evidence of inbreeding depression in male virility. Borash et al. (2007) also showed that once deleterious recessive alleles affecting late-life fitness rise to a frequency of 0.4 or higher, there will be sufficient numbers of recessive homozygotes to allow for the detection of inbreeding depression. Our present simulations can be used to assess the likelihood of low frequency, neutral alleles in these populations rising to high frequencies of 0.4 or greater. For alleles with an initial frequency of 0.1, the probability of these alleles rising to a frequency of 0.4 or higher is 9.5%. Thus, if there were roughly 10 or more loci affecting late-life fitness, then each population would be expected to have one or more such loci with alleles at a frequency of 0.4 or more.
An important question for experimental evolution is whether a laboratory experimental system will result in the increase of all favourable genetic variants that exist in the experimental population. If favourable genetic variants are at a low frequency, then even with directional selection some may still be lost due to genetic drift. When selection began on the O populations, their ancestral B-type population, the 'IV' population of Rose (1984), had already been in the laboratory environments for five years and was most likely not too far from genetic equilibrium with respect to selection. Once O-type age-specific selection started, we assume that there was strong directional selection at some loci for traits that enhance net reproduction late in adult life. As we have mentioned previously, we believe that adaptation to this novel demographic selection regime occurred primarily by means of existing genetic variation, not new mutations. If we make standard simplifying assumptions that these selected phenotypes were affected by multiple loci, as suggested by the results of Burke et al. (2010) for closely related populations, and that the favoured phenotypes had a simple additive pattern of inheritance, we can then do some simple, but instructive, numerical experiments which we describe next.
The progress of directional natural selection will be limited by the number of loci with genetic variation for selected phenotypes as well as the strength of selection on those variants. We have studied the dynamics of selection at two loci in the O populations using a simple additive fitness model. At the A-locus, the three genotypes, aa, Aa and AA contribute h 1 , h 2 and h 3 , respectively, to fitness. At the B-locus, the three genotypes, bb, Bb and BB contribute k 1 , k 2 and k 3 , respectively, to fitness. Thus, the fitness of the AaBB genotype, for instance, would be h 2 + k 3 and so on. We express the strength of selection, β, as the ratio of the difference in fitness between the most and least fit genotypes to the fitness of the least fit genotype.
When zygotes are sampled rather than gametes, accounting for selection involves determining the genotype production of each possible mating pair of males and females. This has to be done since the contribution of each female to the zygote pool varies each generation. Thus, with two loci and two alleles at each locus, there are 10 different genotypes and thus 55 different male-female pairs. The frequency of each mating pair will depend on the genotypes of the females residing in the population and the randomly selected males chosen for them to mate with. In addition, the contribution that each mating pair makes to the pool of zygotes will depend on the randomly chosen fecundity assigned to each female.
When directional selection is relatively weak (β = 0.01), the probability that a rare allele (p 0 = 0.01 at each locus) will be lost is not different from that of a neutral allele (figure 6). When selection is roughly 10 times stronger (β = 0.1), there is a significant drop in the probability of a rare allele being lost, relative to the neutral case (figure 6). In the O populations, favourable genetic variation that was initially rare, e.g. there are only about 20 copies of the favourable allele in the population, is expected to have had a high The favoured alleles, A and B, had initial frequencies of 0.01. The figure shows the probability of these initially rare alleles being lost at only one of the two loci. The recombination fraction was set to 0.01 and 0.5. Since the results did not differ significantly with these two recombination fractions, the results have been pooled and the averages are shown here. The initial population was constructed from a random sample of gametes assumed to be Hardy-Weinberg and linkage equilibrium. Drift was modelled by sampling zygotes. probability of being lost over the course of selection if their fitness effect is weak, e.g. a 1% benefit or less. However, when the fitness effects are at least moderate, e.g. a 10% benefit, then this probability goes down dramatically (figure 6).

Genetic hitchhiking and soft sweeps
The dynamics of alleles at neutral loci may be affected by closely linked selected loci (Maynard Smith and Haigh 1974). We studied this by looking at three linked loci: two outer loci under selection and a central neutral locus. To guide our modelling efforts, we relied on results obtained from the recent genomic analysis of comparable laboratory populations from our laboratory which were selected for rapid development and early reproduction (ACO populations: Burke et al. 2010). Burke et al. (2010) found numerous examples of genes with alleles initially at low-tomoderate frequencies (10-20%) which then increased in frequency (to 70-80%) over the course of selection, but were not fixed. While the precise population genetic mechanism controlling these polymorphisms is not known, their observations are consistent with stable polymorphisms that have moved to a new equilibrium value due to the deliberately imposed change in the laboratory selection regime undergone by the ACO populations relative to their CO ancestors (vid. Chippindale et al. 1997;Rose et al. 2004).
We have modelled this situation by assuming there existed a two-locus polymorphism in the ancestral (B) condition and then these alleles came under selection due to the new O selection regime. We used a two-locus, multiplicative selection model to determine the dynamics at two selected loci. At the A-locus, the three genotypes, aa, Aa and AA contribute h 1 , h 2 and h 3 , respectively, to fitness. At the B-locus, the three genotypes, bb, Bb and BB contribute k 1 , k 2 and k 3 , respectively, to fitness. Thus, the fitness of the AaBB genotype, for instance, under the multiplicative model would be h 2 k 3 and so on. This model can generate stable polymorphisms where the equilibrium frequencies of the a and b alleles are (h 3 − h 2 )/(h 3 + h 1 − 2h 2 ) and (k 3 − k 2 )/(k 3 + k 1 − 2k 2 ), respectively. These allele frequency equilibria with linkage disequilibrium (D) equal to zero will be stable when (from Bodmer and Felsenstein 1967), where r is the recombination fraction. When r is less than this critical value, it appears that the allele frequency equilibria still hold, but D is no longer 0 at equilibrium (Ewens 1979). In our three locus simulations, we set h 1 = k 1 = 0.6, h 2 = k 2 = 1 and h 3 = k 3 = 0.9. With these fitness values, the equilibrium frequency of both the A and B alleles is 0.8. This value was chosen to be similar to the polymorphisms observed by Burke et al. (2010). We started the simulations with the A and B alleles at an initial frequency of 0.2 and r = 0.01. Gametes were assumed to be drawn from a population in linkage equilibrium, and random samples of gametes (which may not have been in linkage equilibrium) were generated to start each simulation. The neutral C-locus was equidistant from the A and B loci and had a recombination fraction of 0.005 with each of these two loci. This recombination fraction corresponds to around 125,000-250,000 base pairs in Drosophila. The C-locus also had two alleles and the simulations were started with several different initial frequencies of the C allele: 0.01, 0.02, 0.03, 0.04, 0.05, 0.1, 0.25 and 0.5. For this study, we only simulated drift through the sampling of gametes since sampling zygotes would have been involved following 231 mating pairs. Heterozygosity is expected to decline over time due to genetic drift, but this loss may be accelerated if there is strong selection on closely linked loci. We compared the distribution of heterozygosity after 160 generations of selection under two conditions: (i) a neutral locus with no linked loci and two alleles at frequency 0.5 and (ii) a neutral locus with flanking loci under strong selection using the multiplicative model and parameter settings described in the previous paragraph ( figure 7). In each case, most of the probability density for heterozygosity at the neutral locus is close to the starting value of 0.5. However, the left tail of the distribution is somewhat longer in the case with linked selected loci. In fact, the average heterozygosity after 160 generations is 0.46 for the single unlinked locus and 0.43 for the neutral locus linked to selected loci with r = 0.0005 (about 12,500-25,000 bp). This difference, although small, is statistically Figure 7. The distribution of heterozygosity after 160 generations of drift in an O population. The curve marked 'neutral' shows the results for a single locus with two alleles at an initial frequency of 0.5. The additional curves show the density function for a neutral locus linked to two selected loci on either side of the neutral locus. The recombination fraction between the neutral and selected loci is indicated by r. significant (P = 0.0003). When r = 0.005, the change in heterozygosity is not statistically significant. These results suggest that strong selection and very tight linkage at an existing polymorphic locus may cause small declines in heterozygosity, but would not be the cause of a substantial decline in heterozygosity.
Neutral alleles may, by chance, also develop linkage disequilibrium with a selected locus and then rise to high frequency as the selected alleles increase in frequency. To study this, we estimated the probability of alleles at a neutral-linked locus being lost in 160 generations and compared this to the chance of loss for a single-unlinked-neutral allele (figure 8). These results show no significant differences between the loss probabilities when neutral alleles are linked to a selected locus and when they are not. There is one exception, and that was observed when the starting neutral allele frequency was set at 0.01, the smallest starting allele frequency used in these simulations. In that case, the neutral locus linked to selected loci showed a slight, but statistically significant, elevated probability of allele-loss.
We have studied the one significant case described in the previous paragraph in more detail. After 10 generations of selection, we estimated the correlation between one of the selected loci and the neutral locus (figure 9). We saw that, when there is no selection at a linked locus, the distribution of the correlation coefficient is nearly symmetric around 0, indicating little linkage disequilibrium. However, when there is selection at the linked locus, the distribution spreads out and has a long negative tail. This suggests that there is occasionally a strong negative correlation between the favoured A allele and the neutral B allele. Such a negative correlation presumably leads to an accelerated decline in the frequency of the B allele as the A allele increases in frequency. When the initial frequency of the B allele is 0.01, this acceleration is sufficient to significantly increase the probability of loss of the B allele. We note that, even though this skewed distribution is seen when the initial frequency of the neutral B allele is 0.05, the occasional negative correlations that develop at this higher starting allele frequency are not sufficient to cause an increase in the fixation probability.  (B and b). Selection was determined by the two-locus genotype at the A and C loci that were separated by r = 0.01. If we let the four gametes, AB, Ab, aB and ab have frequencies x 1 , x 2 , x 3 and x 4 , and the allele A have frequency, p and the B allele frequency, q, then the correlation coefficient is estimated as (1−q) . The frequency of the B alleles is indicated as p 0 . In the simulation labeled 'neutral', the selection coefficients at the A locus were all set to 1, thereby eliminating the action of selection at the linked locus.

Maintenance of genetic variation
With reliable estimates of effective population size, the rate at which variation should be lost from these laboratory populations due to drift can be calculated. While, we do not have detailed information on the genetic variability of the B and O populations at present, we do have substantial information on somewhat similar populations in our lab called the ACO and CO populations (Burke et al. 2010). The methods for maintaining the ACO and CO populations are similar to the B and O populations in terms of the number of breeding adults, food resources, larval rearing conditions, etc. Consequently, it is reasonable to assume that the effective population size of the ACO and CO populations may be in the range of effective sizes we have estimated for the B and O populations.
The five ACO populations were each derived from a single CO population sharing a numerical subscript (e.g. ACO 1 was derived from CO 1 ), but these two sets of populations underwent different numbers of generations of selection prior to the genomewide study of Burke et al. (2010). Each ACO population underwent 605 generations of selection since sharing a common ancestor with its matched CO population, while the CO's underwent 252 generations of selection, prior to the study of Burke et al. (2010). If each of these populations had the same effective population size, N, we would expect the ratio of heterozygosity in a CO population to its matched ACO population to be, (4) It is known from genetic data on 688,000 polymorphic SNP's that the average heterozygosity in the ACO populations is 0.343 and in the CO populations 0.374 (Burke et al. 2010). The ratio of these values is 1.088 ± 0.0014 (95% confidence interval). Using equation (4) and N = 646 (the average B N (e) e ), the expected ratio would be 1.19. Thus, the ACO population seems to be losing less variation than expected. The effective population size of the individual ACO and CO populations would have to be close to 2000 to explain the observed difference in rate of heterozygosity loss by drift alone. However, during the first 150 generations of selection, the ACO populations were under intense selection for rapid development and extremely early fertility. If anything these condition may have resulted in a reduction of heterozygosity at a faster rate than expected under drift alone, making these observations all the more remarkable. An alternative explanation for these results is that some of the variation in these Drosophila populations is maintained by natural selection and thus is not lost over time, as suggested by the strong selection example in figure 6. While we do not have information concerning the possible phenotypes involved with such polymorphisms, many candidates exist. Certainly early vs late reproduction may be an important phenotypic trade-off also at the root of these genetic polymorphisms ). The relative balance of male and female fitness may also change in these demographic regimes causing different levels of sexual antagonism (Chippindale et al. 2001).

Discussion
Detailed demographic data have been used to estimate the variance, the eigenvalue and the probability-of-allele-loss effective population sizes for the B and O populations of Rose (1984). Two different sampling schemes were used. Assuming females mate once with a randomly chosen male, then the next generation is constructed from a random sample of zygotes from the mated females. We also studied the consequences of sampling gametes at random.
By sampling zygotes at random, the estimated variance effective sizes for the B and O populations were 363 and 348, respectively. The estimated eigenvalue effective sizes for the B and O populations were 646 and 648, respectively. Finally, the estimated probability of loss effective size was 599 for the B populations and 518 for the O populations. When drift is modelled by random sampling of gametes, the average variance effective size estimate for the O populations is 698 and the average eigenvalue effective size is 955. The corresponding numbers for the B populations are 787 and 1084, respectively.
We determined that in the B populations there has been sufficient time for neutral genetic variation to be lost, especially if initial allele frequencies were less than 0.2. Conversely, enough time has elapsed in the B populations for some rare neutral variants to increase to high allele frequencies. In fact, we suggest that precisely this type of dynamic may be responsible for the decreased virility of old B males, which show inbreeding depression at later ages (vid. Borash et al. 2007), unlike the O populations.
We can now assess the claims of Promislow and Tatar (1998) and Harshman and Hoffmann (2000) that the apparent difference in the longevities of the B and O populations may be due to deleterious late-acting alleles, made neutral in the B environment that rose to high frequencies by genetic drift. We have estimated that, even after 700 generations, the B populations may have had a few neutral alleles rise to high frequency (>40%) if the trait is affected by 10 loci. However, the observed differentiation of the B and O populations for major life-history traits like longevity and age-specific fecundity occurred in just a few dozen generations. These types of effects cannot possibly explained by drift induced rise to high frequency of deleterious alleles in the B populations. Further, recent studies of hybrid flies created by crossing all possible combinations of B lines failed to show the expected hybrid vigour for longevity or age-specific fecundity (Borash et al. 2007). Harshman and Hoffmann (2000) also point out with concern that long term selection can produce extreme and unusual changes in populations. To support their claim about unusual and extreme outcomes of selection Harshman and Hoffmann (2000) cite a paper reviewing >500 generations of selection for positive and negative geotaxis (Ricker and Hirsch 1988). However, the populations used for that study had a history of extreme population bottlenecks, including one generation of reproduction with one male and two females (Ricker and Hirsch 1988). Even under better times, these populations consisted of only 60 breeding pairs, and these were separated into six groups of 10 with the 10 pairs with the highest or lowest geotaxis scores being put into one bottle, and so on. Thus, the opportunity for inbreeding, if there was between-family variation in geotaxis, was increased by this protocol (Erlenmeyer-Kimling et al. 1962). Thus, it is not laboratory selection causing extreme outcomes, but severe inbreeding, drift and population bottlenecks in the problematic case which Harshman and Hoffmann (2000) raise.
Recent genomic studies of 10 other populations in our laboratory that were also selected for different ages of reproduction have revealed extensive allele frequency differentiation, without much loss of genetic variation (Burke et al. 2010). We have shown that for selection on balanced polymorphisms, where the newly favoured alleles have initial frequency much greater than (2N) −1 , the effects of selection on linked neutral variation in populations of the effective population size that we estimate here are to cause only very modest reductions in heterozygosity (cf. Teshima et al. 2006). Some type of soft sweep or shifts between stable balanced polymorphisms may be a primary factor influencing the evolutionary dynamics of our laboratory-selected populations, rather than fixation due to drift or hitchhiking effects.
Using the estimated genomewide heterozygosities in the Burke et al. (2010) populations, and the estimates of effective population size in this paper, we found greater average heterozygosity in the ACO populations than expected on the basis of genetic drift alone. One possible explanation for this observation is that some of the genetic variation in these Drosophila populations may be maintained by some form of balancing selection. On numerous occasions, our laboratory adapted populations have shown that they can respond rapidly to reverse selection (Teotónio and Rose 2000;Teotónio et al. 2002;Joshi et al. 2003). These observations support the idea that phenotypic adaptation in outbred laboratory populations of Drosophila can be accomplished by means of standing genetic variation at many segregating loci without extensive loss of genetic variation, rather than selective sweeps that purge genetic variation. This is not to discount the existence of effectively neutral genetic variation, given the small effective population sizes of these populations. Such genetic variation no doubt exists and will be gradually lost over many generations of laboratory evolution. However, our findings and those of Burke et al. (2010) imply the maintenance of more segregating, selectively significant, genetic variation than might have been expected based on conventional population-genetic thinking. A similar view of Drosophila evolution in natural environments has also been proposed recently by Karasov et al. (2010). For now, at a minimum, it is clear that a model of phenotypic evolution dominated by selective sweeps leading to fixation of newly arisen, large effect, favourable mutations is not an appropriate sole explanatory hypothesis for these cases of experimental evolution in fruitflies.