Evolution of a transposon in Daphnia hybrid genomes
© Vergilino et al.; licensee BioMed Central Ltd. 2013
Received: 20 June 2012
Accepted: 6 November 2012
Published: 6 February 2013
Skip to main content
© Vergilino et al.; licensee BioMed Central Ltd. 2013
Received: 20 June 2012
Accepted: 6 November 2012
Published: 6 February 2013
Transposable elements play a major role in genome evolution. Their capacity to move and/or multiply in the genome of their host may have profound impacts on phenotypes, and may have dramatic consequences on genome structure. Hybrid and polyploid clones have arisen multiple times in the Daphnia pulex complex and are thought to reproduce by obligate parthenogenesis. Our study examines the evolution of a DNA transposable element named Pokey in the D. pulex complex.
Portions of Pokey elements inserted in the 28S rRNA genes from various Daphnia hybrids (diploids and polyploids) were sequenced and compared to sequences from a previous study to understand the evolutionary history of the elements. Pokey sequences show a complex phylogenetic pattern. We found evidence of recombination events in numerous Pokey alleles from diploid and polyploid hybrids and also from non-hybrid diploids. The recombination rate in Pokey elements is comparable to recombination rates previously estimated for 28S rRNA genes in the congener, Daphnia obtusa. Some recombinant Pokey alleles were encountered in Daphnia isolates from multiple locations and habitats.
Phylogenetic and recombination analyses showed that recombination is a major force that shapes Pokey evolution. Based on Pokey phylogenies, reticulation has played and still plays an important role in shaping the diversity of the D. pulex complex. Horizontal transfer of Pokey seems to be rare and hybrids often possess Pokey elements derived from recombination among alleles encountered in the putative parental species. The insertion of Pokey in hotspots of recombination may have important impacts on the diversity and fitness of this transposable element.
Transposable elements (TEs) are widespread in the living world. Few multicellular eukaryotes lack these genetic elements, which have the ability to move and multiply throughout the genome, but see for exceptions in endosymbiont genomes. There is substantial variation in the proportion of mobile genetic elements and especially TEs across genomes. Numerous studies have explored whether TE history reflects the evolutionary history of their hosts (that is, co-evolution[4–7]). Due to their high densities in eukaryotic genomes and their irreversible mode of insertion, some authors have proposed that SINEs (Short INterspersed Elements, non-autonomous Class I TEs) may be used for inferring the evolutionary history of their host. The phylogeny of a strictly vertically transmitted TE (from parents to offspring) should reflect that of its host(s). However, genetic and evolutionary processes, which may occur during speciation, may affect subsequent phylogenetic inferences based on TE sequences. For example, horizontal transfer (that is, transmission by vectors unrelated to host reproduction), incomplete lineage sorting (that is, the maintenance of ancient polymorphism in closely related species), lack of genetic variability or hybridization/introgression may account for incongruent phylogenetic patterns between TEs and their hosts[9–13]. Hybrid genomes should carry TE(s) from the genomes of their parental species. In contrast, if horizontal transfers occur between species that share ecological niches or are infected by the same parasites, we do not expect that TE phylogenies will match species phylogenies.
The Daphnia pulex complex has been intensively studied due to its dominance in freshwater habitats in North America, and its variation in reproductive mode and ploidy level. Daphnia usually reproduce by cyclic parthenogenesis, which is clonal reproduction interrupted by bouts of sexual reproduction. However, some lineages reproduce by obligate parthenogenesis[14–18]. The D. pulex complex includes seven species that have been distinguished on the basis of morphological, ecological, and genetic data[19–22]. In North America, two species of the complex, D. pulex and Daphnia pulicaria, are dominant in freshwater habitats, and are morphologically similar but ecologically distinct, although they hybridize in nature[18, 23–25]. Diploid hybrids always have D. pulex mitochondrial genomes (mtDNA) and have been found to reproduce by obligate parthenogenesis in nature based on patterns of allozyme variation[18, 26] and laboratory crosses, but see for experimentally produced hybrids capable of sexual reproduction. Analyses of allozyme variation indicate that introgression is rare or nonexistent in areas where the two species co-occur[16, 18, 25] but a recent study using mitochondrial and nuclear markers has shown that introgression between these species had a substantial impact on their evolutionary history.
Hybridization has also occurred among other species in the D. pulex complex[29–32]. For example, hybridization between North American species and circumarctic or European Daphnia species has been suggested based on allozyme, mtDNA or nuclear markers such as the ribosomal intergenic spacer (IGS)[21, 33, 34]. Some hybrids found in arctic and subarctic regions are polyploids (mainly triploids). Polyploids with D. pulex, D. pulicaria and Daphnia middendorffiana (a circumarctic species) mtDNA are thought to have arisen from hybridization between D. pulex and D. pulicaria, or with another species which no longer exists as a cyclic parthenogen[22, 29, 30, 32]. Another group of polyploids has been found in Daphnia tenebrosa, a circumarctic species and includes both triploid and tetraploid clones but their hybrid nature is unknown.
The D. pulex genome has recently been sequenced and numerous Class II DNA transposons have been identified from it. One Class II element, Pokey (Subclass I; sensu) of the piggyBac superfamily, has been extensively studied in non-hybrid populations. It inserts in a site-specific fashion into a unique TTAA site in the tandemly repeated 28S rRNA genes[39–41] and is also encountered in other parts of the genome. Two Pokey elements from D. pulicaria have been completely sequenced and show that the element has 16 bp imperfect terminal inverted repeats (TIRs) and encodes a putative transposase that is much longer than originally estimated and contains an intron (Y. Bigot, personal communication,). In addition, the non-coding region upstream of the transposase gene contains a repeated sequence that is similar to a sequence in the unique region of the host’s ribosomal IGS[34, 39, 40]. Based on patterns of polymorphism observed among natural populations, previous studies[41, 42] have suggested that Pokey may be active in sexual populations of D. pulex but not in obligate parthenogens. Although Class II DNA transposons are more prone to transfer horizontally between lineages than retrotransposons, horizontal transfer does not seem to play a major role in the evolution of Pokey as 28S rRNA genes and mtDNA phylogenies of species in the subgenus Daphnia are congruent with the phylogeny of Pokey elements leading to the conclusion that Pokey elements may co-evolve with their hosts.
The purpose of the present study is to describe Pokey sequences in hybrid isolates from the D. pulex complex. Specifically, we aim to determine if Pokey sequences inserted in rRNA genes reflect the putative parental species origins of the hybrid isolates. As Pokey elements insert in 28S rRNA genes, they are subject to genetic processes such as inter-allelic recombination, gene conversion and unequal crossing over that may lead to the concerted evolution of these genes. If Pokey elements are strictly vertically inherited and if recombination is negligible, we expect that a phylogeny of Pokey elements should follow the host phylogeny, and Pokey elements encountered in hybrid genomes should correspond to those encountered in their host parental species. If recombination is important in Pokey evolution and has occurred between elements showing recent activity and after hybridization, hybrid isolates will carry either sequences inherited from both parents or mosaic sequences between parental types. Alternatively, if recombination occurs between inactivated elements that diverged in the distant past, we could see Pokey mosaic sequences that do not necessarily reflect the progenitor species of the hybrid isolates in the D. pulex complex.
After cloning and sequencing partial Pokey elements from representatives of species in the D. pulex complex, we: (a) detected recombinant elements and located the recombination breakpoints in each, (b) estimated recombination rate parameters across sequences using a coalescent based method, (c) determined parental sequences from which the recombinant elements arose and (d) reconstructed the evolutionary history of the Pokey sequences. Moreover we used a PCR-RFLP protocol to search for the presence of different Pokey elements in additional hybrid and non-hybrid isolates of the D. pulex complex across freshwater habitats and geographical locations.
We analyzed 51 isolates representing six of the seven lineages of the D. pulex complex identified by Colbourne et al.. Only D. melanica was not included. Nine Pokey sequences were obtained from nine isolates studied by Penton and Crease. We established parthenogenetic lines (isolates) in the laboratory from 42 individual females sampled from nature between 2004 and 2010 (see Additional file1) and cultured them using standard techniques. For 50 of these isolates, genomic DNA was extracted from 10 to 30 individuals weighing approximately 100 mg (wet weight) using the DNeasy Tissue kit (QIAGEN Inc., Mississauga, ON, Canada) according to the supplier’s protocol. Species identity of these isolates was assessed by combining information on morphology, haplotype of the mitochondrial ND5 gene and genotype of the nuclear Lactate dehydrogenase (Ldh) gene as described in. Ploidy levels were previously assessed using microsatellite genotypes and flow cytometry.
RT-PCR was performed on RNA samples extracted from a sexual D. pulex isolate to determine if Pokey transposase transcripts were present and to provide further evidence for the presence of an intron in the transposase gene. RNA was extracted from the parthenogenetic offspring of one sexual D. pulex isolate PX2-ON-9 (provided by Dr. M. Cristescu at the University of Windsor). The tissue was stored in RNA Later (Qiagen) at −20°C and RNA was extracted using the RNAqueous-4PCR kit (Life Technologies Inc., Grand Island, NY, USA) and standard manufacturer’s protocols. The absence of DNA contamination was verified using standard PCR and Pokey transposase primers. Reverse transcription of the RNA used the SuperScript III One-Step RT-PCR with Platinum Taq Kit (Life Technologies) and standard manufacturer’s protocols. Primers used were Pok4065F (5’-TGATTCACCGAGGCCTCAGTTC) with Pok4488R (5’-GAATCGCTCGCGAGTCATGG) and Pok5026F (5′-TCGAACCTGCAGCCGGACGAATTTGCAG) with Pok5985R (5’-CACGTCGGTTAGAATATTCTGGCTCGTCGG). Numbers correspond to nucleotides in the 6.6 kb element from D. pulicaria (GenBank accession # AY115589.1,).
To examine the impact of hybridization on Pokey diversity, Pokey alleles inserted in 28S rRNA genes from various isolates of the D. pulex complex were cloned and sequenced. An approximately 1600 bp fragment of Pokey from 16 isolates (2 diploid hybrids, 5 polyploid hybrids, 1 introgressed polyploid (TE3-MB-4), 3 diploids, 3 polyploids with unknown hybrid state, 3 diploid non-hybrids) of the D. pulex complex (see Additional file1) were amplified using Pok5026F and the primer 28SR (5′-TCCATTCGTGCGCGTCACTAATTAGATGAC), which is located 46 bp downstream of the Pokey TTAA target site. PCR reactions were performed using the Phusion™ high-fidelity PCR kit (Finnzymes, Woburn, MA, USA) according to the supplier’s protocol. The thermocycling profile consisted of 1 cycle of 1 minute at 94°C, 35 cycles of 30 seconds at 94°C, 30 seconds at 55°C and 2 minutes at 72°C, with a final incubation of 5 minutes at 72°C. PCR amplicons were cloned using the StrataClone™ Blunt PCR Cloning Kit (Agilent Technologies Inc., Santa Clara, CA, USA) according to the supplier’s protocol. For each isolate, 10 to 12 plasmids chosen randomly from a minimum of 40 isolated colonies were purified using the E.Z.N.A. Plasmid Mini Kit II (Omega Bio-tek Inc., Norcross, GA, USA). Pokey sequences obtained from these plasmids were named using the labels of the isolate from which they were cloned plus the number that corresponds to the colony from which they were picked (for example PX2-QC-8_28 is the Pokey allele isolated from colony 28 produced during the cloning of PX2-QC-8 Pokey sequences). These Pokey alleles were compared to Pokey sequences from a previous study which focused on non-hybrid isolates to ensure that the elements analyzed in this study belong to the Pokey lineage. Partial Pokey sequences were manually aligned using the ClustalW module of BioEdit v.7.0.5; available athttp://www.mbio.ncsu.edu/BioEdit/bioedit.html.
Potential recombination events among Pokey alleles were explored in a dataset of 53 sequences in which ambiguous sites and indels were deleted to prevent identification of false recombination events. The phi-test, implemented in Splitstree v4.10, and GARD (Genetic Algorithm for Recombination Detection, available via the DataMonkey website athttp://www.datamonkey.org) were used to determine if recombination events may be detected in our Pokey sequence dataset as well as in partial transposase sequences included in them.
Hudson and Kaplan minimum recombination events, Rm were calculated with the DnaSP v5 software. Moreover we used kwarg (available at http://www.stats.ox.ac.uk/~lyngsoe/section26/), a heuristic alternative to the branch and bound method implemented in beagle, to calculate Rh, which gives an estimate of the lower bound of the minimum number of recombination events. Kwarg does not guarantee that the minimal recombination history will be found but returns a history with a low number of recombination events under the infinite sites assumption. We ran kwarg 2,048 times with the default scoring scheme and recorded the number of recombination events for each iteration.
As all the recombination analyses detected recombi nation events (see Results), we estimated Pokey‘s population recombination rate parameters along the sequence using a full likelihood coalescent approach based on a finite-sites model implemented in LAMARC v2.1.6. This approach estimates the parameter r LAM = c/μ where c is the recombination rate per site per generation and μ is the neutral mutation rate per site per generation. This model considers that (a) all recombination is homologous, (b) there is no association between recombination frequency and sequence divergence, (c) gene convergence and interference did not occur, (d) recombination events are selectively neutral and (e) recombination rate does not vary along the sequence. The F84 evolutionary model was used as it is a closer fit to the substitution model, Tamura 3-parameter with gamma shape and invariant sites, inferred with MEGA5 than the other model (GTR) offered by LAMARC. The F84 model is also computationally faster. To account for substitution rate differences between non-synonymous sites and all other sites, two categories of relative mutation rate were assigned. Our sampling strategy included 20 initial chains of 20,000 genealogies and two final chains of 1,000,000 genealogies with 2,500 samples discarded per chain. Adaptive heating was used to improve the search of parameter space with initial relative temperature of 1, 1.1 and 2. The entire analysis was replicated four times and the results were combined using the algorithm of Geyer. Coalescent-based analyses in LAMARC estimate both genetic diversity at the population level, θ LAM = 4N e μ (where μ is the mutation rate per site for a diploid population and Ne is the effective population size) and the recombination rate r LAM = c/μ (where c is the rate of recombination per site per generation). Therefore, the population recombination rate parameter can be calculated as R = θ LAM r LAM = (4N e μ)(c/μ) = 4Nec.
Nucleotide diversity among sequences (that is, π and θ) and mutation rate heterogeneity along sequences (the Gamma shape parameter, α) may have a dramatic impact on the performance of different methods to detect recombination breakpoints[56–58]. For example, Posada and Crandall have shown that different methods may have a different propensity to detect recombination when nucleotide diversity varies and mutation rate heterogeneity is fixed. Moreover, they have shown that mutation rate heterogeneity may lead some methods to detect false positives[56–58].
We estimated nucleotide diversity (π and θ) using DnaSP v5.10.01. Mutation rate heterogeneity (α) and a nucleotide substitution model were inferred using a maximum likelihood framework in MEGA5. To detect recombination breakpoints and identify which parental sequences formed the recombinant alleles, we applied a stepwise protocol of the maximum chi-square method[59, 60] using the stepwise package v0.1.1 (available athttp://stat.sfu.ca/statgen/research/stepwise.html). The maximum chi-square method is a sliding window approach that computes a chi-square statistic from a 2 x 2 table with counts of matches and mismatches between the two half windows separated by a proposed breakpoint. As high mutation rate heterogeneity in Pokey sequences (see Results) may lead to false negatives when window half-widths are too small, we used different values of window half-width (70, 80, 90 and 100 nucleotides) and 100,000 Monte Carlo replicates for the permutation distribution. The maximum chi-square analysis was reiterated including breakpoints detected by previous steps, as the stepwise protocol specifies, until no further breakpoint was detected. This maximum chi-square analysis identifies two sequences involved (one derived, one parental) for each recombination event detected. Recombination breakpoint detection by the maximum chi-square method is relatively more powerful and less prone to type I error (false positives) than other algorithms[56, 57] when sequences are subject to high mutation rate heterogeneity. Each event detected between different sequence pairs was checked by eye.
Dendrograms were produced for each fragment separated by recombination breakpoints highlighted by GARD analysis of the 53 sequences using the Neighbor-joining algorithm. We constructed a phylogenetic network using the NeighborNet algorithm implemented in Splitstree v4.10. To understand more accurately the process and evolutionary history of recombination in the Pokey fragments cloned here, we constructed an ancestral recombination graph using kwarg[51, 52, 62]. We chose the simulation with the lowest Rh value to build this graph.
Given the intra- and intergenomic variability of Pokey alleles in isolates from which elements were cloned, and the high number of recombination events in these sequences (see Results), we aimed to explore the presence/absence of these alleles in additional isolates belonging to the D. pulex complex. To do this, we developed a restriction enzyme analysis (that is, RFLP) based on Pokey sequences analyzed so far and the recombination breakpoints detected by the GARD analysis. We then screened Pokey alleles from 41 isolates (see Results) for RFLP haplotypes. PCR reactions (two replicates for each isolate) were performed in 25 μL reaction mixtures containing 1x Econotaq PCR buffer, 2 mM MgCl2, 50 μM dNTP, 0.1 μM of each primer (Pok5026F and 28SR), and 1 unit of EconoTaq polymerase (Lucigen Corporation, Middleton, WI, USA) or using the Phusion™high-fidelity PCR kit according to the manufacturer’s protocol (Finnzymes). PCR amplifications were performed using the thermocycling profile described above. Restriction enzymes that cut at specific sites in some Pokey alleles may be used to highlight Pokey sequences that have undergone recombination events. We used NebCutter v2.0 (http://tools.neb.com/NEBcutter2/) on each of the 53 sequences to choose enzymes that would highlight different pure and recombinant Pokey alleles. The restriction enzymes Dra I, Bsp HI and Bst EII (New England BioLabs Inc., Ipswich, MA, USA) were used in a two-step protocol. Dra I, Bsp HI and Bst EII enzymes cut at 380 bp, approximately 790 bp, and approximately 1080 bp. The eight possible conformations were differentiated on a 3% agarose gel. In the first step, Dra I and Bsp HI were used in the same mix containing NEBbuffer4 and digestion was conducted at 37°C for two hours. In the second step, NEBbuffer3 at a 1X final concentration, BSA at 1U/μl and Bst EII at 1U/μl were added to the solution and digestion was conducted at 60°C for one hour. Digested PCR products were separated on 3% agarose gels for two hours at 95 volts.
Cloning and PCR amplifications may produce chimeric sequences (for example, recombinants between alleles) and lead to an overestimation of recombination in nature. We assessed this hypothesis by re-amplifying Pokey sequences from isolates used in our phylogenetic survey using the PhusionTM high-fidelity PCR kit and by increasing elongation steps of the thermocycling profile to three minutes instead of two minutes. Increasing the elongation steps should reduce the probability of producing chimeric sequences. Results were compared to PCR-RFLP results with two-minute elongation steps to determine if the recombination pattern is reduced. Moreover, if recombination events identified in the data are not artifacts, PCR-RFLP results should be concordant with expectations based on the cloned elements that we sequenced.
Nucleotide diversity and mean similarity between sequence pairs for each part of the Pokey sequences
Recombination events in which groups of sequences (see text for more explanation) are involved according to maximum χ 2 and GARD dendrogram analyses
Events involving sequences according to the maximum χ2test (other sequences involveda)
Pairs of probable parental sequences estimated using GARD dendrograms based on two different recombination breakpointsb
I (h, g, c, PX2-QC-8_29), IV (b, EPC2-SP-2_1), V (h), VII (g)
h : c or i : c
I (i), VII (g, PX2-SF-8_29, PX2-SF-8_1, f)
h : c or i : c and c : h
II (i), IV (a)
d : a or EPC2-SP-2_1 : a or EPC2-DE-3_1 : a
d : c or EPC2-SP-2_1 : c or EPC2-DE-3_1 : c
d : i or EPC2-SP-2_1 : i or EPC2-DE-3_1 : i
II (i), IV (a, h)
b : i, b : c or b : a
b : i, b : c or b : a
I (a), IV (i), V (h), VI (g)
d : i
I (a), VII (TE3-MB-4_9)
i : f or h : f
g : f
I (a), VI (c), VII (a, TE3-MB-4_9)
unknown : f
g : h or g : i and i : g or h : g
IV (EPC2-SP-2_1), V (a, c)
i : Unknown or a : Unknown
I (TE3-MB-4_9), II (b, EPC2-SP-2_1), III (PX3-QC-1_20), IV (c)
h : Unknown
Groups of sequences with the same recombination site may represent a single recombination event between ancestral sequences. For example, the recombination breakpoint of events A (a : c), B (a : h), C (a : PX2-QC-8_29), D (j : TE3-MB-4_9) and E (i : a) were estimated to be bp 592 to bp 607 for all pairs. Thus, we considered them to be the same recombination event, I. This approach results in seven different recombination events named I (bp 592 to bp 607), II (bp 748 to bp 800), III (bp 815 to bp 856), IV (bp 825 to bp 856), V (bp 917 to bp 920), VI (bp 952 to bp 957) and VII (bp 1,038 to bp 1,073) (see Additional file4). Although it is still possible, it seems unlikely that all recombination events between sequence pairs were independent and not due to past events in the common ancestors of groups of sequences.
The Hudson-Kaplan minimum number of recombination events (Rm) in our partial Pokey sequences dataset is 30 according to the algorithm implemented in DnaSP v5.10.01. According to the kwarg algorithm, Rh, the estimated lower bound of the minimum number of recombination events, follows a normal distribution with mean 96.39 and a standard deviation of 3.10.
Using LAMARC to co-estimate mutation rate θ LAM = 5.94 × 10-2 (with 95% support intervals of 3.63 × 10-2 and 8.16 × 10-2) and the overall recombination rate r LAM = 5.09 × 10-1 (with 95% support intervals of 3.57 × 10-1 and 7.68 × 10-1) allowed us to estimate a population recombination rate of 3.02 × 10-2 recombination events per site per generation for Pokey elements.
Non-synonymous polymorphism at the nucleotide level leads to a high level of polymorphism at the aa level and, coupled with recombination events, produces a diversity of partial aa sequences of the Pokey transposase (see Additional file2).
Additional file5 shows an ancestral recombination graph with 89 recombination events (Rh) identified by the Song and Hein algorithm. The evolutionary history of Pokey allele TE3-MB-4_12, which contains 38 recombination events, is highlighted in red. This ancestral recombination graph shows putative past recombination events that were not detected by the maximum chi-square method. For example, the evolution of EPX2-DE-1_1, encountered in a European D. pulex isolate, seems to involve 25 past recombination events.
The v1.1 genome sequence of D. pulex was scanned for our partial Pokey sequence using TBLASTN to explore which alleles are present. Sequences from contigs with E-value = 0.0 and without indels more than 10 nucleotides long were aligned with our 53 sequences to generate a dataset of 58 sequences. A NeighborNet phylogenetic network (see Additional file6) shows that Pokey elements from the genome sequence group with sequences encountered in diploid and polyploid isolates from the pulicaria group (PC2-SK-5 and PC3-QC-3), in diploid D. tenebrosa (TE2-MB-1) and in European D. pulicaria (EPC2-CZ-1).
Taxonomic, habitat and geographical distributions of RFLP haplotypes
Percentage and taxonomic occurrence
33.33% of isolates, mainly with D. pulex mitochondria
QC, ON, MB (CAN), MI (USA)
38.09% of isolates, D. pulicaria, polyploids with D. pulicaria mitochondria, diploid/triploids with D. pulex mitochondria and one D. tenebrosa isolate
26.19% of isolates, D. pulex, diploid hybrids with D. pulex mitochondria and triploid hybrids D. pulicaria mitochondria, European D. pulicaria and D. tenebrosa isolâtes
30.95% of isolates, D. pulex, diploid hybrids with D. pulex mitochondria and triploid hybrids D. pulicaria mitochondria, D. middendorffiana, European D. pulicaria and D. tenebrosa isolates
One introgressed isolate (TE3-MB-4)
rock bluff pools
Churchill MB, (CAN)
9.76% of isolates, D. pulex, D. pulicaria or D. middendorffiana and D. pulex-D. pulicaria hybrids
ponds or rock bluff pools
Kujjuarappik and Metis, QC (CAN) Churchill, MB (CAN)
One D. pulex-D. pulicaria hybrid (PX2-MI-11)
26.8% of isolates, D. pulex and D. pulicaria
Our results show that recombination has a dramatic impact on the evolution of Pokey alleles in the rDNA of Daphnia. The recombinant sequences do not seem to be due to cloning artifacts as some isolates may have just one recombinant Pokey allele in their genome (for example, triploid isolates PC3-SK-5 and PC3-MB-4 or diploid isolates EPC2-CZ-1, EPC2-SP-2 and PX2-QC-1, Additional file7). Increasing the extension time of the PCR cycles to avoid artefactual recombination events during amplification did not change the RFLP results (data not shown). Moreover, Daphnia isolates show RFLP patterns that are concordant with those expected from the cloned Pokey sequences (see Additional file7).
Pokey sequences show substantial nucleotide diversity (π = 0.048 and θ = 0.037, Table 1), undergo recombination at a rate of 3.02 × 10-2 events per site per generation and rarely spread by horizontal transfer. As predicted in the case of recombination events between Pokey elements that co-evolve with their host and have undergone recent activity, hybrid isolates from the D. pulex complex carry either sequences inherited from putative parental species or mosaic sequences produced by recombination between parental types. Most Pokey sequences from pulex-pulicaria hybrids (for example, PX2-MB-1, PX3-QC-1_1, PC3-QC-1_28 or PX2-QC-8_28) either cluster with Pokey alleles from D. pulex (group h) or D. pulicaria (PC2-QC-4 and PC-SK-5_2 in group g) in the phylogenetic network (Figure 2) or are recombinants between sequences from these groups (that is, PX3-QC-1_20). Pokey sequences from polyploid isolates with mtDNA from either D. pulicaria or D. middendorffiana from Manitoba, Canada seem to be recombinants between sequences found in D. tenebrosa (group d) and D. pulicaria, D. pulex or D. arenata (that is, groups g, h or i, respectively, Table 2). This suggests that D. tenebrosa may be one of the parental species of these polyploid isolates, which are thought to originate from hybridization between D. pulex and D. pulicaria or with another species that no longer exists as a cyclic parthenogen[22, 29, 30, 32]. Phylogenetic analysis of Ldh A sequences is also concordant with a hybrid origin of polyploid isolates with D. middendorffiana mtDNA (Figure 3 in).
The introgressed triploid isolate, TE3-MB-4, which has D. tenebrosa mtDNA and most likely a D. pulex nuclear genome[22, 31]), carries alleles from group ‘a + TE3-MB-4_9’ with a minor trace from D. tenebrosa and major traces from D. pulex or D. pulicaria alleles. The presence of alleles with traces of D. pulex or D. pulicaria in D. tenebrosa or European D. pulicaria isolates (TE3-MB-1, TE3-MB-3, TE2-MB-1, TE2-MB-3, EPC2-MB-2 and EPC2-MB-3) may be a remnant of hybridization events between these four species. Similarly, Pokey alleles from European D. pulicaria isolates (EPC2-SP-2_1 and EPC2-DE-3_1) also seem to be recombinants between alleles from North American D. pulex or D. pulicaria and D. tenebrosa (recombination events G, K and L in Table 2 and Figure 2). Although Vergilino et al. did not detect evidence of hybridization between species in the pulicaria group (North American D. pulex and D. pulicaria) and the tenebrosa group (D. tenebrosa and European D. pulicaria) using microsatellite loci, Weider et al. did find evidence for it using allozyme data. Moreover, Ambrose and Crease have shown that IGS segments from some isolates of European D. pulex may have originated from North American D. pulicaria. As North American D. pulex and D. pulicaria have been found in Europe, introgression of Pokey alleles from North American species into European species is possible and may generate recombinants that persist for long periods of time and are able to spread over broad geographic areas.
Although most Pokey sequences from hybrid isolates are mosaics, which is consistent with recent recombination between alleles from different parent species, some recombinant alleles may be more ancient and represent retained ancestral polymorphism (that is, incomplete lineage sorting). The clustering of recombinant alleles in group b, which are encountered in D. tenebrosa isolates, a hybrid D. pulex (PX2-QC-8) and a polyploid with D. pulicaria mtDNA (PC3-QC-3) suggest that these recombinant alleles are ancient.
Pokey alleles from the isolate whose genome was sequenced cluster with alleles from group f (Figure 3 and Additional file6) whereas Pokey alleles from other D. pulex isolates cluster in group h. This may be the result of incomplete lineage sorting or a signature of past hybridization events that have occurred in the population from which the genome isolate was sampled. Alternatively, the genomic Pokey alleles may have a different evolutionary history than their paralogs in rDNA. Moreover, these genomic alleles show a limited signature of recombination according to the GARD analysis. As the recombination rate is variable across the Daphnia genome[65, 66], we may expect this result if Pokey elements are inserted in genomic regions showing low recombination rates.
The number of recombination breakpoints detected in our sampled sequences by different methods was highly variable. Most recombination breakpoint detection algorithms underestimate the number of recombination events and their accuracy varies depending on whether the method uses summary statistics or the maximal information in the sample. The maximum chi-square method performed in a stepwise manner and for different window sizes detected a minimum of seven unique recombination events between Pokey alleles from different species of the D. pulex complex (see Additional file3, Table 2). This value is likely to be an underestimate as a high recombination rate may lead to overlapping recombination events. As the maximum chi-square is a substitution-based method that compares substitution distributions from two windows flanking a tested site, recombination breakpoint detection will be dependent on the size of the windows used when mutation rate heterogeneity is high and/or nucleotide diversity is low. Small windows should increase the probability of detecting recombination breakpoints in high diversity regions but may be unable to detect recombination events in low nucleotide diversity regions. Increasing the size of the windows to increase recombination detection accuracy in regions with low nucleotide diversity may lead to less accuracy in the regions with multiple recombination events. The estimate of Hudson-Kaplan minimum recombination events (Rm) is 30 and the haplotype lower bound (Rh) estimate is 96.39 ± 3.10. Rm is a parsimony-based method (that is, the four gamete test) and may underestimate the actual lower bound of recombination events, particularly when the recombination rate is high[67, 68]. In contrast to the former estimators, Rh, estimated using a maximum-likelihood framework that uses the maximal information in the sample, is probably more reliable than Rm and shows the highest lower bound of recombination events. The accuracy of Rm and Rh to estimate the lower bound of recombination events has not been tested on sequences subjected to low and high heterogeneous mutation rates. If numerous recombination events occurred during the evolutionary history of Pokey elements, we may need sophisticated phylogenetic tools and analyses, such as the Ancestral Recombination Graph, to deepen our understanding of their evolutionary history.
The recombination rate of Pokey sequences estimated in our survey is high (3.02 × 10-2 events/site/generation) but comparable to the value (2.0 to 6.0 × 10-2 events/generation) calculated by McTaggart et al. based on changes in the relative frequency of rDNA length variants in four apomictically-propagated Daphnia obtusa lines over 90 generations. To our knowledge, this is the first time that the rate of recombination has been explicitly estimated for a DNA transposon. This estimate has to be viewed with caution as some of the assumptions of the algorithm used, such as the absence of gene conversion and interference, and the selective neutrality of recombination events, may be violated. Even so, Rh and Rm estimates lead to the conclusion that the recombination rate is high between Pokey elements. Such a high rate of recombination between Pokey elements may explain why no significant differences were found between the amount of variation in 28S rRNA genes with or without Pokey insertions in cyclic and obligately asexual D. pulex isolates as recombination may take place between the Pokey elements as well as the 28S rRNA genes themselves.
Numerous studies have focused on the relationship between recombination rates in host genomes and TE distributions[70–73] or dynamics[37, 41, 74]. The distribution and abundance of TEs may be influenced by selective forces[75–77], the stochastic process of mutation and/or the availability of insertion sites. One theoretical model based on the action of selection posits that the impact of recombination events between TEs inserted in non-homologous loci is deleterious (that is, ectopic exchange model) and, therefore, insertions will accumulate in regions of the genome with low rates of recombination. However, empirical data from sequenced genomes show that the relationship between TE density and recombination rate depends on the organism and the TE class or family[4, 39, 71, 79–82].
Recombination between TE copies of the same family has been reported before[83–87], but the significance of this for TEs in general is rarely discussed. Our survey focused on variation in an approximately 1,600 bp sequence. This dataset does not include the TIRs known to be essential for transposition or the IGS sequences upstream of the transposase gene that may influence expression of the transposase. If the TIRs or transposase gene have diverged in different species and recombination events produce mosaic sequences, the transposase may not be able to recognize TIRs with which it has not co-evolved. Thus, we might expect that many recombination events between divergent TE alleles will inactivate their mosaic product. Alternatively, Schaack et al. recently suggested that TEs might benefit from recombination in a similar manner to the host, by generating new variants that are able to evade sequence-specific host suppression machinery. Recombination between different insertions could bring together independent, beneficial nucleotide substitutions that may then be favored by intragenomic selection and result in an increase in copy number of the new, recombinant insertion.
Fragments of the Pokey transposase gene have different evolutionary histories. Moreover, non-synonymous polymorphism coupled with recombination has produced a substantial level of Pokey variability at the aa level (Table 1, Additional file2). This variability is primarily seen in exon 2, which is thought to contain a nuclear localization signal and a sequence similar to a Plant Homeo Domain or cysteine-rich zinc finger that is possibly involved in chromatin or protein-protein interactions[40, 88, 89]. For example, the fourth cysteine residue of the putative zinc finger is replaced either by a serine or a tyrosine residue in 3 of the 40 different aa sequences we analyzed (see Additional file2). This variation might translate to flexibility in target site interactions or changes to dimerization dynamics during the formation of the transposition complex, which could have implications for mobility and transposition efficiency.
DNA transposons in Caenorhabditis elegans are located preferentially in recombination hotspot regions whereas retrotransposons are not, and Duret et al. suggested a role for recombination in the transposition process. Glass et al. did find indirect evidence that Pokey insertions in rDNA could spread through unequal crossing-over, so perhaps there is some credence to this suggestion. In addition to Pokey in Daphnia, the retrotransposons R1 and R2 in Drosophila and numerous other arthropods are also found inserted into rDNA. Using the ectopic exchange model as a framework, Zhang et al. demonstrated that eukaryotic hosts may tolerate a high load of retrotransposable elements in their rDNA because they generally have many more rDNA copies than the minimum required for rRNA synthesis. The location of Pokey elements in rDNA may be advantageous for both element and host, as frequent recombination between different copies located in these loci may increase the efficiency of intragenomic selection, and insertions may be only mildly harmful due to the multi-copy nature of rDNA. Pokey haplotype 2 is a recombinant allele (Figures 1, Table 2) that was found in numerous hybrid and non-hybrid Daphnia isolates having different origins and occupying different locations and ecosystems (Table 3). This might be such a successful haplotype due to the amalgamation, through recombination, of separately beneficial mutations into a single lineage that was then able to spread to different populations and species within the D. pulex complex.
Most of the hybrid isolates from the D. pulex complex analyzed in this study carry Pokey sequences inherited from putative parental species or mosaic sequences produced by recombination between parental types. In addition, recombination may play an important role in generating life history variation in this TE. Future studies should test the activity of recombinant Pokey alleles to assess if recombination events have had a significant effect on their transposition capacity. Additionally, testing for recombinants among the other families of TEs found within D. pulex would be useful to determine whether rates of recombination among TEs depend in part on the region of the genome they inhabit or on other properties such as mode of transposition. Moreover, the possible beneficial effects of recombination on TE-level evolution in general should be investigated more thoroughly by estimating rates of recombination for different TEs in various genomes.
bovine serum albumin
Genetic Algorithm for Recombination Detection
polymerase chain reaction
restriction fragment length polymorphism
short interspersed elements
terminal inverted repeats
RV acknowledges a scholarship from Centre d’études Nordiques and a travel fellowship from UQAR. PDP is supported by a research grant from the Canada Research Chair program. This work was supported by Discovery Grants from the Natural Sciences and Engineering Research Council (NSERC) of Canada to FD and to TJC. We thank Yves Bigot for pointing out the existence of the intron in the Pokey transposase gene.
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.