Androdioecious populations
All populations reported here were derived from a population with standing genetic diversity previously evolved for 140 generations under our standardised laboratory environment [55],[70]. A schematic of the lab environment is shown in Figure 1B. This lab-adapted population is designated EEVA6140, where EEV stands for the Wormbase lab acronym, A for androdioecy, and 6140 for 140 generations of experimental evolution of replicate population #6. We drop EEV henceforth.
In order to maintain similar levels of genetic diversity among the ancestor populations used here, 500 A6140 hermaphrodites were individually mated with an excess of A6140 males. All F1 lineages were then expanded in numbers for two generations. A total of 444 of these were recovered, mixed in equal proportions and cryogenically frozen at -80°C at high densities [102]. This androdioecious population is named A00.
Trioecious and dioecious populations
The fog-2(q71) allele located in chrV:25°CM was introgressed in A6140 in two stages, the first with the goal of obtaining homozygous lineages for this allele in a A6140 genetic background, and the second with the goal of introducing A6140 genetic diversity.
In the first stage of the introgression, 50 fog-2(q71) homozygous females were recovered from a 100 generation lab-adapted male-female population (D2100; described in [55],[71], and individually crossed to an excess of A6140 males. F1 hermaphrodites were selfed to generate F2s, while the F1 male progeny were kept at 4°C in order to be used in subsequent crosses. F2 progeny were collected two days later as immature L4 staged larvae and scored for fog-2(q71) homozygosity by the accumulation of unfertilised oocytes during the next two days [80]. From each of the fifty F1 lineages, eight F2 females were outcrossed with an excess of F1 males from a different lineage, and their F3 female progeny mated with one sibling male to generate the F4s. F3 parents were PCR genotyped to confirm homozygosity at the fog-2 locus (see below). Only a small region surrounding the fog-2 locus is expected to contain the N2 wild isolate genetic background from where the fog-2(q71) allele is ultimately derived [80] (results not shown). This is because the fog-2 locus is located in the telomere and the origin of the allele for the introgressions is from a lab-adapted population with standing genetic diversity.
In the second stage of the introgression, 12 females from each of 50 F1 lineages were outcrossed with an excess of A6140 males. A total of 600 F1 hermaphrodites were selfed and F2 progeny scored for fog-2(q71) homozygosity. A total of 600 F2 females were next mated with males randomly coming from the 50 F1 lineages. The F3s were PCR genotyped and the wild type allele was estimated to be at 5.2 × 10-2. F3 females were mated with sibling males to generate the F4s. Since we did not control for the homozygosity of fog-2 in males, the wild type allele in the F4s is expected to segregate at 2.6 × 10-2. A total of 444 F4 lineages were recovered, separately grown to high densities and mixed in equal proportions to constitute the ancestral trioecious population, named T00. Samples were frozen at -80°C.
To derive the ancestor dioecious population (D00), 350 adult T00 females were individually mated with single T00 males, after scoring for the accumulation of unfertilised eggs in them. Twenty-four hours later, males were genotyped at the fog-2 locus and 318 lineages were recovered where both parents were fog-2(q71/q71). These lineages were grown to high densities, mixed at equal proportions and samples frozen at -80°C.
Monoecious population
To derive the ancestor monoecious population (M00), the xol-1(tm3055) allele located in chrX:-0.45°CM was introgressed in A6140 in two stages.
In the first stage, parental hermaphrodites from strain FX03055 [103] were mated with an excess of A6140 males. F1 progeny were separately selfed, and with PCR genotyping, F2 hermaphrodites were confirmed to be xol-1(tm3055) homozygous after laying of the F3 embryos (see below). Lineages were then kept for another generation. This four-generation cycle was repeated another two times, with the parental hermaphrodites coming from a homozygous xol-1(tm3055) F3 lineage obtained in the first cycle. Ten crosses were done in the first cycle and twenty in the second and third cycles. At the F1 generation, 70, 286 and 378 hermaphroditic lineages were obtained in the first, second, and third cycle, respectively. At the F2 generation, 3,800, 1,228 and 1,219 were obtained in the first, second, and third cycle, respectively. For the F3 generation we obtained 752, 350 and 1,144 lineages in the first, second, and third cycle, respectively. At the end of each cycle, F3 individuals were genotyped at 57 bi-allelic SNPs spanning the X chromosome (see below and Additional file 3 for SNP information). Two lineages containing less than 0.5 Mbp of the N2 wild isolate genetic background surrounding the xol-1 locus were recovered (results not shown).
In the second stage of the introgression, 300 hermaphrodites from each of the two lineages recovered after the first stage were individually mated with an excess of A6140 males. Their F1 progeny were selfed to generate 4,800 F2 hermaphrodites, each then being selfed in order to obtain the F3s and, subsequently, the F4s. Each of the 4,800 F2s was collected after reproduction for PCR genotyping. A total of 444 F2 homozygous xol-1(tm3055) lineages were recovered. They were grown to high densities, mixed in equal proportions and samples frozen at -80°C.
Green fluorescent protein testers
As previously described [70], a fully penetrant and genome-integrated GFP marker with a myo-3 promoter was introgressed in the A6140 population to obtain a tester population with standing genetic diversity (A6140GFP). Heterozygous and homozygous genotypes similarly express GFP in all muscle cells. A6140GFP is here used in the (population-wide) fitness assays (see below; Figure 4, Additional file 1: Figures S3-S5).
As previously described [70], A6140GFP hermaphrodites were selfed for 12 generations, in order to obtain highly inbred strains. One of these strains (A6140GFPL1) is here used in the hermaphrodite/female fertility and male fitness assays (see below; Figures 2 and 6).
Experimental evolution under different salt regimes
Frozen samples with >104 individuals of the A00, T00, D00 and M00 populations were thawed and after one passage, 104 individuals were seeded for each of rR1-# replicate populations; where r designates the salt regime (Sudden, Gradual or Control), R the reproduction system (Androdioecy, Trioecy, Dioecy, Monoecy), and # the replicate population number. Dioecious populations were cultured only in the sudden regime. A total of 43 populations underwent experimental evolution (see Additional file 1: Table S1 for full designations).
Following our standard laboratory environment (Figure 1B) [55], populations were kept in ten 9°Cm Petri plates with 28 mL of solid NGM-lite media (Europe Bioproducts, Cambridge England) covered by an overnight grown lawn of HT115 Escherichia coli. Bacteria provided ad libitum food for worm development from the L1 larval stage until adult reproduction. NGM-lite media contains NaCl at 25 mM (0.14% w/v). At 24 ± 2 hours of the life cycle, each population was seeded with 1,000 first larval staged (L1) individuals in each of the 10 Petri plates. After development to maturity for 72 ± 2 hours at constant 20°C and 80% relative humidity (RH), all 10 plates were mixed and harvested worms exposed to 1 M KOH: 5% NaOCl `bleach’ solution for 5 minutes ±15 seconds, to which only embryos survive [102]. After repeated washes with the M9 isotonic solution, embryos were maintained in a shaker in 3 to 5 mL M9 in 15 mL conical tubes, at 20°C and 120 rpm. After 24 ± 2 hours, adult debris was removed after centrifugation at 200 rpm and the density of live L1s estimated under a dissection scope. The appropriate M9 volume with live L1s was then placed in fresh NGM-lite plates to complete one life cycle.
The control regime was the same in which the lab-adapted population had been cultured for 140 generations. In particular, NGM-lite media in the plates were not supplemented with NaCl. Replicates were cultured for 68 generations. The sudden regime was characterised by the same standard conditions, except that the NGM-lite plates were supplemented to 305 mM NaCl from the start of evolution. NaCl was dissolved in ddH2O and the appropriate volume added to NGM-lite before autoclaving (1.78% w/v, minimum 99% purity, Roth P029.3). ST, SM and SA populations were cultured for 38 generations, SD populations were cultured for 50 generations. For the gradual regime NGM-lite plates were supplemented with increasing concentrations of NaCl from 33 mM at generation 1 to 305 mM NaCl at generation 35 and onwards. Replicates were cultured for 68 generations. With the exception of SD1-4, all other populations concurrently underwent experimental evolution. Samples from each population were periodically stored at -80°C at high densities (>103) for posterior characterisation (see Additional file 1: Table S1).
Male frequency assay in A6140
After thawing frozen samples of A6140 in 25 mM NaCl NGM-lite plates, worms were passaged twice as during experimental evolution. The parental P generation was then seeded into either 25 mM plates or 305 mM NaCl plates, in triplicate, at the L1 larval stage. At 96 ± 4 hours of the life cycle, 300 adult individuals were sexed under a dissection scope at 15× magnification, following random trajectories and covering the whole surface of the plate. Samples from each plate were separately cultured to generate the F1 generation, which was similarly maintained at either 25 mM NaCl or 305 mM NaCl and sexed when individuals reached adulthood. The assay was repeated eight times (blocks). As the data is proportional, for analysis we did generalised linear mixed effects models (GLMM) with binomial (logit-link) errors to estimate fixed differences between the P and F1 generation at each environment or to estimate fixed differences between environments at the F1 generation [104]. Only the later analysis is reported. Block was modelled as the random independent variable. Significance of effects was assessed with z-ratio tests. The function glmer in the package lme4 within R was used for computations [105],[106]. In one additional assay, we further detailed the effects of NaCl on male proportions (see Additional file 1: Figure S1). These assays were done over two generations in triplicate in NGM-lite plates containing 150 mM, 250 mM, 275 mM or 305 mM NaCl.
Male frequency assay in experimentally evolved populations
Ancestral (A00, T00) and evolved (rA1-3, GT1-3, ST1/3/5) populations were measured for male proportions once experimental evolution was completed (see Additional file 1: Table S1). The assay was done in 12 blocks, to include samples for generations 5, 15, 25 and 35 from the sudden regime, generations 35, 50, 56 and 68 from the gradual regime, and finally, generations 15, 35 and 68 from the control regime. Frozen samples of >105 individuals in each population were thawed and passaged for two generations in 25 mM NaCl NGM-lite plates. In the second generation, at 96 ± 4 hours of the life cycle, 300 individuals were sexed per each of three plates per population sample. We did not statistically model male frequency experimental evolution.
fog-genotype frequency assay
Frozen samples were thawed and passaged once in 25 mM NaCl NGM-lite plates. In the second generation, 48 L3 or L4 larval staged individuals per population and experimental regime were PCR genotyped at the fog-2 locus (see below). Samples from generations 5, 15, 25 and 35 were chosen for ST1/3/5 populations, generations 35, 50, 56 and 68 for GT1-3 populations and generation 23 for CT1-3 populations (see Additional file 1: Table S1). In each PCR, positive and negative control samples were included with gDNA of individuals with known genotypes that were constructed prior to the assays. After quality control of the data, we had an average of 42.7 ± 4SD observations per population and at each time point. Assuming binomial sampling distributions, we can thus detect wild type allele frequencies above 2% with a statistical power of 0.8. Note that since we did not detect the fog-2(wt) allele in CT1-3 at generation 23, it is impossible for transitions to occur by genetic drift in subsequent experimental evolution [70].
The `transition fitness’ of the fog-2(wt) allele (s), the invading allele in the trioecious populations, can be calculated as its expected deterministic frequency change over the resident fog-2(q71) allele [70],[82]. Note that transition fitness is different from `invasive fitness’ since the invader allele is at relatively high frequencies during most of experimental evolution, compare with [82],[83]. We employed ANOVA to estimate s as the fixed (continuous) effects of generation on the natural log ratio changes of fog-2(wt)/ fog-2(q71) allele frequencies, assuming random mating and selfing and no sex ratio segregation distortion (see Additional file 1: Table S2). When we failed to detect the fog-2(wt) allele (for example, generation 5 of the sudden treatment), we considered it at a frequency of 1/97, while when all samples were positive, we considered it at a frequency of 96/97. This allowed us to calculate the log ratios for an initial sample size of 48 individuals while being conservative. The lm function in the stats package of R was used for calculations [105].
Competitive fitness assay
We assayed fitness of all evolved populations at generation 35 (rT1-#, rM1-#, and rA1-#), and the three ancestral populations (T00, M00, A00; Additional file 1: Table S1). For this, we performed head-to-head competitions between experimental populations and the A6140GFP tester population. For each competition, more than 103 individuals from each competitor were thawed from -80°C in parallel and expanded in numbers for two generations in 25 mM NaCl NGM-lite plates. On the third generation, L1-staged GFP individuals were mixed with wild type experimental individuals at expected relative frequencies of 50%; see [70] for a density calibration curve of this protocol. Worms then developed, matured and reproduced in NGM-lite plates supplemented to 305 mM NaCl. For T00, M00 and A00 we also followed the same protocol and assayed them at 25 mM NaCl (results presented in Additional file 1: Figure S5). GFP scoring was done in the following generation, after bleaching and 24 ± 2 hours of maintenance in M9, by photographing 5 μl of M9 containing live L1 larvae on a glass slide, at a resolution of 1.5 pixel/μm under a fluorescent dissection scope equipped with a digital camera.
The competition assays were done in 12 blocks, each having a different thawing day. In each block, the wild type experimental populations, GFP tester and the A6140 populations were thawed in parallel, such that a total of 12 independent competitions between the GFP tester and the A6140 population were also obtained (see Additional file 1: Figure S3). Each competition was replicated between four and six times, with 3 × 103 individuals (three plates) being employed per replicate. All replicates started from the same mix of GFP and wild type experimental individuals. Setup and scoring was randomised across treatments and experimenters within block. The mean number of GFP/wild type L1s scored per population sample after one generation of competition was of 330/902 with the 2.5% and 97.5% quantiles being of 92/243 and 1,031/2,618 individuals, respectively. Frequency estimates of wild type and GFP alleles were also obtained for the mixes at setup (see Additional file 1: Figure S4).
Relative changes in the wild type allele frequencies over the GFP allele in each competition, measured at the L1 stage and in two consecutive generations, are expected to reflect variation in all fitness components: L1 to adult developmental time and survivorship, mating, fertilization, fecundity and embryo to L1 hatching success. Without loss of generality, a haploid fitness coefficient can be estimated as the change in ratio in one generation of the wild type non-GFP allele over the tester GFP allele [70],[83],[85]: w = ln (pwt.t1/pGFP.t1) = ln (pwt.t0/pGFP.t0); where p
wt
and p
GFP
are the wild type and GFP allele frequencies, respectively, at setup t0, and after the competition, t1. Note that this fitness coefficient is defined in a similar way as the transition fitness coefficient (see above) [82],[83].
Since the scoring procedure employed in the assay was based on the presence/absence of GFP expression, we were unable to score for progeny heterozygosity. To correct for this, an algorithm was written in R taking into consideration the male frequencies from each of the two competing populations (which were separately measured during the assay following the above described protocol) and the starting and final GFP counts in order to retrieve the expected number of alleles after the competition. This algorithm assumed random mating and selfing, and no sex segregation distortion [52] (see Additional file 1: Table S2). The GFP allele does not have a fitness cost over the wild type allele of the A6140 population when in high salt (see Additional file 1: Figure S3).
Due to the high heterogeneity among blocks and because at a given block not all experimental populations were assayed (see Additional file 1: Figure S3), the fitness estimates of the experimental populations obtained at a given block were transformed by subtracting the average w of the A6140 population for the corresponding block (w
t
). Quality control of w
t
data involved removing from analysis all estimates from the ST2 population because there was a significant deviation from the expected 50:50 L1 mix of wild type: GFP at assay set up that likely determined the observed negative fitness (see Additional file 1: Figure S4; see [70] for a demonstration of a change in the sign of fitness coefficients because of frequency dependence, using similarly designed assays).
Analysis of w
t
involved two steps, the first with the goal of obtaining the ordinary least-square estimates of the three ancestral populations (T00, A00, M00) and the second with the goal of determining the extent of fitness responses after 35 generations of evolution in all derived populations (G35; rT1-#, rM1-# and rA1-#; Additional file 1: Table S1), which we define here as `adaptation’ or `adaptive rates’ to the high salt environment. For the first step, w
t
was modelled as a function of the fixed reproduction system with ANOVA (see Additional file 1: Figure S5). We employed the function lm within the stats package in R for calculations. The least-square estimates obtained per ancestral population were then subtracted from the transformed data of the corresponding G35 experimental populations for subsequent analysis (Δw
t
; Figure 4). For each reproduction system, Δw
t
data were modelled as a function of the fixed salt regime (sudden, gradual or control) and random replicate population, with linear mixed effects models (LMM) using maximum likelihood methods [104]. We did not model the ratio of wild type to GFP alleles at setup since preliminary models did not find it to be significant (not shown). Significance of fitness responses at each regime were assessed with z-ratio tests. We employed the lmer function within the lme4 package in R for calculations [106]. Post hoc comparisons among salt regimes within each reproduction system were done with Tukey tests and assuming Student t distributions where the degrees of freedom were asymptotically determined, using the lsmeans package in R [107].
Selfed versus outcrossed hermaphrodite fertility assay
We measured the fertility of A6140 hermaphrodites when selfing and when outcrossing occurs in high or low salt (Figure 2). Stocks were thawed and in parallel passaged once in 25 mM NaCl NGM-lite plates. In the second generation, 50 L4 larval stage (reproductively immature) hermaphrodites were transferred to 9°Cm 25 mM NaCl plates and mated with an excess of males coming from the A6140GFP inbred tester strain. GFP outcrossed progeny grew on these plates and as L4s they were individually transferred to 6°Cm NGM-lite plates, with either 25 mM or 305 mM NaCl, and placed together with two A6140GFPL1 tester males or allowed to self without the presence of males. Parents were removed from the plates 24 ± 2 hours later, and three to four days later the number of adult viable progeny and GFP status was scored. The assay was done in two blocks. Setup and scoring was randomised across treatments and experimenters.
In order to compare selfed and outcrossed treatments, we eliminated all observations where less than five progeny were scored (eliminated 66 observations), and when less than 10% of male progeny in the outcrossed treatment were observed (eliminated 14 observations). We tested for the fixed effects of breeding mode and salt environment on fertility with GLMMs by taking block as the random independent variable. We considered that the model error followed a Poisson distribution by using the log-link function, since the data were scored as counts. Least-square estimates were obtained by maximum likelihood. z-ratio tests were done to test for significance of fixed effects and Tukey post hoc z-tests were done to test for differences in breeding mode in each salt environment. For calculations we employed the glmer function within the lme4 package in R [104],[106]. For plotting, ordinary least-square estimates are presented.
Hermaphrodite and female fertility assay
A6140 hermaphrodites and D00 females were measured for fertility when outcrossed to D00 males after development since the L1 larval stage in high salt conditions (Figure 3). A6140 and D00 stocks were thawed and passaged twice in parallel in 25 mM NaCl NGM-lite plates. On the third generation, L1 individuals were seeded in two 305 mM NaCl NGM-lite plates. After 72 ± 2 hours of development, A6140 hermaphrodites and D00 females were each individually placed with two D00 males in 6°Cm 305 mM NaCl NGM-lite mating plates with 5 uL of E. coli. All adults were removed 24 ± 2 hours later and the plates left to incubate in standard conditions. After four days all progeny was scored for number.
Preliminary analysis showed that data including zero fertility followed a negative binomial distribution (for A6140, central tendency `mean’ ± CI =1.5 ± 0.3; for D00, mean =3.74 ± 0.83; both similar dispersion parameters; calculated with the fitdistr function in the MASS package in R [104]). For this reason, we decided to eliminate all zero fertility observations (for A6140, we eliminated 24 observations; for D00, 12) and then tested for differences among populations with a generalised linear model (GLM) employing a Poisson distributed error. A z-ratio test was done to test for significance of effects. For plotting, ordinary least-square estimates are presented.
Developmental time to maturity and population fertility rate assay
By measuring the population fertility rate at several time points during the standard life cycle (Figure 1B), we were able to estimate developmental time to maturity in A6140 under conditions of low and high salt (Figure 2). For D00 and SD1-4, we only estimated fertility rates at one time point (Figures 3 and 7). All populations were concurrently assayed.
A6140, D00 and SD1-4 stocks were thawed and passaged twice in parallel in 25 mM NaCl NGM-lite plates. On the third generation, 1,000 L1s were seeded either in three 25 mM NaCl NGM-lite plates or in three 305 mM NaCl NGM-lite plates. For A6140, we hand-picked 24 hermaphrodites after 48 ± 2 hours or 64 ± 2 hours from the L1 seed, from each of the three plates in both salt conditions, and placed them in 12-well cell culture plates with 3 mL 25 mM or 305 mM NaCl NGM-lite and with 5 uL of E. coli. Two hours later, wells were scored for the presence of fertilised embryos (identified by their non-transparent appearance). For D00 and SD1-4, similar numbers of individuals were picked at 64 ± 2 hours after L1 seed and growth in 305 mM NaCl NGM-lite plates. Setup and scoring was randomised across treatments.
For analysis of A6140 data, we employed GLM on the count data of fertile to infertile individuals per assay plate and tested for the fixed effects of time of life cycle and NaCl condition (Figure 2). A6140 and D00 data were similarly compared at 305 mM NaCl with GLM (Figure 3). Experimental evolution responses in 305 mM NaCl in SD1-4 were compared to the values of D00 with a GLMM for the fixed effects of generation and the random effects of replicate population (Figure 7). For all models binomial errors were considered by using the logit link function (as the data is proportional data) and z-ratio tests were done to test for significance of effects. We employed glm or the glmer functions within the stats or the lme4 packages in R, respectively, for calculations. For plotting, ordinary least-square estimates are presented.
Male fitness assay
Generation 50 SD1-4 and ancestral D00 stocks were thawed and grown in parallel under standard 25 mM NaCl conditions for two generations. At 72 ± 2 hours of the third generation, nine experimental males were transferred to 6°Cm 305 mM NaCl NGM-lite plates spotted with 5 uL E. coli, together with nine A6140GFPL1 tester males and twenty-two day 4 adult fog-2(q71/q71) females (from the inbred line D0L27, described in [70]). An average of 16.9 ± 4.2SD females per mating plate were transferred to a new plate 24 ± 2 hours later and killed with 30 μl of the bleach solution. These plates were incubated during the next three to five days at 20°C and 80% RH, and the viable adult progeny scored for GFP expression at 30× magnification under a stereoscope. Five replicate mating plates were done per population and all populations assayed over four blocks. Setup and scoring was randomised across treatments.
Quality control of the adult progeny count involved discarding mating plates where less than six females were transferred and also those where less than twenty progeny were scored. Similarly to transition fitness and (population-wide) fitness (see above) per mating plate, a male haploid fitness coefficient was defined as: wm = ln(pt1/(1-pt1))- ln(pt0/(1-pt0)); with p
t0
being the fixed 0.5 wild type allele frequencies at the set up generation and p
t1
the observed wild type allele frequencies after the competition [85]. We assumed random mating and no sex ratio segregation distortion. For analysis we conducted LMM by taking the number of females from which progeny was scored as a continuous covariate, generation as a fixed independent variable and block as the random independent variable (Figure 7). z-ratio tests were done to determine the significance of evolutionary responses. The lmer function within the lme4 package in R was employed for calculations.
Simulations of fog-genotype frequency dynamics
We numerically simulated the expected change in the frequency of the fog-2(wt) allele during trioecious experimental evolution in the sudden and gradual regimes as a function of the relative outcross-fitness parameter (α). Despite following the previous C. elegans modelling of androdioecy [52],[90], we updated it to the expected reproduction under trioecy, as detailed in Additional file 1: Table S2.
Simulations started with N = 104 individuals, where the proportion of the different sexes and breeding modes were given by the observed allele frequencies at the fog-2 locus in the rT1/3/5 populations (see Additional file 1: Table S1). fog-2 genotype frequencies in males were assumed to be the same as those of hermaphrodites/females. αm matings occurred, with α being outcross-fitness and m the observed male frequency, indiscriminately between males and hermaphrodites or between males and females. Mating pairs were randomly selected by collecting either αmN females and hermaphrodites or the total number of females and hermaphrodites (the number that was lower), together with the same number of males. Males were sampled with replacement. During mating, we followed genotype identity for subsequent reference. If αmN was less than the total number of hermaphrodites and females, then unmated individuals were allowed to constitute pseudo mating pairs with probability 1 if they were hermaphrodites and 0 otherwise. This means that the β parameter of [52], which defined the success rate of the non-outcrossed gamete partition, was set to one under selfing. Accounting for these pseudo mating pairs was necessary to ensure Mendelian segregation at fog-2 under selfing. At the next stage, all of the mating pairs and pseudo mating pairs were sampled with replacement until a new population with N individuals was obtained. The fog-2 genotypes were reconstituted by random sampling of alleles from each of the parental genotypes. Half of the individuals that resulted from outcrossing events were defined as males. The remaining individuals, irrespective of being generated by outcrossing or selfing, were defined as hermaphrodites or females depending on their fog-2 genotype.
Estimates of the outcross-fitness parameter (α) congruent with the observed fog-2 genotype frequency dynamics were obtained under a likelihood framework, whereby fog-2(wt) allele counts over time provided the probability of the α-values. Simulations were done independently for each time period sampled during the evolution. For the sudden populations there were three periods: from generation 5 (G5) to G15, from G15 to G25, and from G25 to G35. For the gradual populations there were three periods: from G35 to G50, G50 to G56, and G56 to G68. Twenty replicate simulations were performed per starting conditions and for each combination of α between 0 and 2 in a grid of 51 points. The probability of the observed fog-2(wt) allele count data given the simulated mean frequencies was obtained from a binomial distribution. The natural logarithms of these probabilities were summed across periods, and the compound value with the higher likelihood was taken as the ML estimate. CIs were defined by -2lnLk around the ML (Figure 6).
Simulated fitness of T00
With the same simulation algorithm, fitness (w
t
) of T00 was estimated for one life cycle of competition with the GFP tester population, assuming fixed setup of 50 wild type to 50 GFP alleles (see Additional file 1: Figure S8). For several α values we obtained the sex ratios after one generation of competition and estimated allele frequencies assuming random mating and selfing and no sex ratio segregation distortion. For the results presented in Additional file 1: Figure S7, an inbreeding depression parameter (δ) was introduced as a sampling weight before reproduction of GFP hermaphrodites relative to the wild type trioecious females, according to the C. elegans androdioecious model (see main text and [52]).
Estimates of outcross-fitness in A6140
Outcross-fitness (α) was calculated for the lab-adapted A6140 population from the male frequency changes observed after two generations of culture in 305 mM NaCl (from Figure 2A). We followed the `male maintenance function, of [52]: mt1 = αm t0/(2βm t0 + β(1-αmt0)), with mt0 being the male frequency after one generation in high salt (P) and mt1 the male frequency after two generations (F1).
DNA extraction
Genomic DNA was prepared from single immature L3 or L4 larval staged individuals using the prepGEM Insect kit (ZyGEM), following standard protocols [108].
PCR genotyping
For the introgression of xol-1(tm3055) in A6140 (see above), we used PCR genotyping with the forward primers XOL2L: 5′-GGGATTGATAGGAGCGAAA-3′, XOL4L: 5′-ATGATTGATGATTTACCGAAGC-3′ and the reverse primer XOL2R: 5′-GGGATTGAGCACCAAACTT-3′. xol-1 (wt/wt) individuals yielded a single 471 bp band, xol-1 (tm3055/tm3055) one 375 bp band, and the heterozygotes both. PCR amplifications were carried out under the following conditions: 50 to 150 ng of genomic DNA were used in a 15 μl reaction that contained 3 μl Green GoTaq© Flexi Buffer 10× (Promega, Madison USA), 1 mM MgCl2, 0.3 mM of primer XOL2L, 0.5 mM of primer XOL4L, 0.5 mM of primer XOL2R, 0.2 mM of each dNTP, and 0.625 units of GOTaq DNA polymerase (Promega, Madison USA). Typically, reactions went through 35°Cycles, after an initial denaturation of 3 minutes at 95°C under the following conditions: 30×seconds at 95°C, 30×seconds at 59°C, 45 seconds at 72’C, followed by an extension step at 72’C for 3 minutes.
To score fog-2 genotype frequency dynamics (see above), we used PCR genotyping with the forward primers F1RFLP: 5′-CTGTCCAGATACGCCTCTCGTCT-3′, FogR4short: 5′-CTGATTGAGCAATATGTCGAATT-3′ and the reverse primer FogR3: 5′-ACGCCTGTGTGAAATTGGGCAAAAGATTAGACTGATTGAGCAATATCGATAATC-3′. fog-2(wt/wt) individuals yielded one 295 bp DNA band visualised in agarose gels, while fog-2(q71/ q71) gave one 264 bp band. Heterozygotes exhibited a double band pattern. PCR amplifications were carried out under the following conditions: 50 to 150 ng of genomic DNA were used in a 15 μl reaction which contained 3 μl Green GoTaq© Flexi Buffer 10×, 1.5 mM MgCl2, 0.5 mM of primer F1RFLP, 0.75 mM of primer FogR4short, 0.25 mM of primer FogR3, 0.2 mM of each dNTP, and 0.625 units of GOTaq DNA polymerase. After 3 minutes at 95°C, 35°Cycles were done each with 30×seconds at 95°C, 30×seconds at 58°C, 30×seconds at 72’C. A final step was done at 72’C for 3 minutes.
Wormbase.org genome release WS220 was used to design the PCR oligos.
SNP genotyping
To determine standing levels of genetic diversity, SNPs were genotyped in 39 populations at generation 22 of experimental evolution and in the ancestral lab-adapted population (see Additional file 1: Table S1). Sixty-nine bi-allelic SNPs along chromosome IV were chosen and genotyped with the iPlex SequenomTM MALDITOF platform following previously described protocols [55],[71]. WS200 was used for the oligonucleotide design employed for PCR amplification followed by allele-specific extension. Oligonucleotide information is available from the authors upon request.
For data analysis, we excluded SNPs with more than 80% missing data across all samples followed by removal of individuals with more than 60% of missing data. In a final step, SNPs with more than 50% of missing data were also removed. The resulting data include genotypes for 58 SNPs. SNP information and sample sizes are presented in Additional file 2.
For the introgression of the xol-1(tm3055) allele in A6140 (see above), 57 bi-allelic SNPs along chromosome X were also genotyped with the iPlex SequenomTM MALDITOF platform (see Additional file 3).
SNP diversity analysis
For each population sample, we calculated homozygosity as one minus the average observed number of heterozygous genotypes across SNPs (1-Ho). We also calculated pairwise SNP linkage disequilibrium (r2) as the composite genotype disequilibria, assuming that the genotype probabilities are the products of the gametic probabilities [109]. SNPs with minor allele frequencies of <0.05 were removed prior to analysis to prevent low sample size bias [110]. Because of this correction, we failed to calculate r2 for SM4 and SM7 (see Additional file 1: Table S3).
To test the dependency of fitness responses at generation 35 (Δw
t
) on homozygosity (1-Ho) or linkage disequilibrium (r2) at generation 22, we did ANCOVAs (Table 1). 1-Ho or r2 were taken as the continuous independent covariate and reproduction system (trioecy, monoecy, androdioecy) as the categorical fixed independent variable. The ST2 genotype data were not included in the analysis. Modelling reproduction system and salt regime as categorical fixed independent variables gave similar results (results not shown). Expanding the scale of 1-Ho or r2, for example by angular transformation, also produces similar results (results not shown). The lm function within the stats package in R was used for calculations.
Availability of supporting data
The data supporting the results of this article are available in the Dryad.org repository, doi:10.5061/dryad.4013. The data are composed of the male and fog-2 genotype frequency data, fitness data, male fitness, fertility data and SNP genotype data. The R scripts used for numerical modelling can be obtained from the authors upon request.