Eurotatorian paraphyly: Revisiting phylogenetic relationships based on the complete mitochondrial genome sequence of Rotaria rotatoria (Bdelloidea: Rotifera: Syndermata)

Background The Syndermata (Rotifera+Acanthocephala) is one of the best model systems for studying the evolutionary origins and persistence of different life styles because it contains a series of lineage-specific life histories: Monogononta (cyclic parthenogenetic and free-living), Bdelloidea (entirely parthenogenetic and mostly benthic dweller), Seisonidea (exclusively bisexual and epizoic or ectoparasitic), and Acanthocephala (sexual and obligatory endoparasitic). Providing phylogenetic resolution to the question of Eurotatoria (Monogononta and Bdelloidea) monophyly versus paraphyly is a key factor for better understanding the evolution of different life styles, yet this matter is not clearly resolved. In this study, we revisited this issue based on comparative analysis of complete mitochondrial genome information for major groups of the Syndermata. Results We determined the first complete mitochondrial genome sequences (15,319 bp) of a bdelloid rotifer, Rotaria rotatoria. In order to examine the validity of Eurotatoria (Monogononta and Bdelloidea) monophyly/paraphyly, we performed phylogenetic analysis of amino acid sequences for eleven protein-coding genes sampled from a wide variety of bilaterian representatives. The resulting mitochondrial genome trees, inferred using different algorithms, consistently failed to recover Monogononta and Bdelloidea as monophyletic, but instead identified them as a paraphyletic assemblage. Bdelloidea (as represented by R. rotatoria) shares most common ancestry with Acanthocephala (as represented by L. thecatus) rather than with monogonont B. plicatilis, the other representative of Eurotatoria. Conclusion Comparisons of inferred amino acid sequence and gene arrangement patterns with those of other metazoan mtDNAs (including those of acanthocephalan L. thecatus and monogonont B. plicatilis) support the hypothesis that Bdelloidea shares most common ancestry with Acanthocephala rather than with Monogononta. From this finding, we suggest that the obligatory asexuality of bdelloideans may have secondarily derived from some other preexisting condition in earlier lineage of rotifers. Providing a more complete assessment of phylogenetic relationships and inferring patterns of evolution of different types of life styles among Syndermata awaits comparisons requiring mitochondrial genome sequencing of Seisonidea.


Background
The Rotifera (also called Rotatoria) is a group of aquatic micrometazoans, mostly less than a millimeter in size. It includes more than 2,000 described species that usually occur in freshwater and marine environments throughout the world, but some species are found in wet terrestrial habitats, such as moist soil, or mosses and lichens living on fallen trees and rocks [1][2][3]. This group is usually distinguished from other metazoans by presence of the corona ("wheel organ" = the crown of cilia) located in the cephalic region, which is used for locomotion and food gathering. Due to their great abundance and high reproductive potential, rotifers have been considered to play a significant role in the food webs of certain freshwater environments through their involvement in energy cycling and nutrient transfer [4][5][6]. Most authorities accept that the Rotifera consists of three classes, each having a unique reproductive strategy [2,7,8]: Bdelloidea (exclusively parthenogenetic), Monogononta (cyclic parthenogenetic with facultative sexual reproduction), and Seisonidea (exclusively bisexual). Of these, the class Bdelloidea is very unique in that its species (several hundred) are exclusively female. This is the largest metazoan asexual group where no sexual reproduction has ever been reported and represents an ancient origin of asexuality with great evolutionary success in the diversification of the species [9,10]. The class Monogononta, representing the largest group of rotifers, comprises more than 1,500 species. They are mostly found in freshwater, brackish and marine waters and are characterized by a unique life cycle, cyclical parthenogenesis, that alternates between two reproductive modes, i.e. amictic and mictic phases of reproduction according to the absence/presence of males, respectively [2,5]. The class Seisonidea consists of three epizoic species (belong to two genera, Seison and Paraseison [11]) on marine crustaceans (Nebalia species) and all of these are known to reproduce by amphimixis [12].
There is little doubt about the close relationship between Rotifera and Acanthocephala, the clade known as Syndermata because they share the feature of the syncytial epidermis [13,14]. This affinity has received broad support from many earlier studies utilizing different sources of phylogenetic information: morphology [15][16][17], SSU rDNA sequences [18][19][20], combined analysis of molecular and morphological characters [21,22], and combined analysis of SSU and LSU ribosomal DNA sequences [23]. In contrast to syndermatan monophyly, internal phylogeny within the clade (among Bdelloidea, Monogononta, Seisonidea, and Acanthocephala) has been the subject of relatively vigorous contention [24][25][26][27]. Most controversies are related to Eurotatorian monophyly versus paraphyly, although there are differences regarding the phylogenetic position of Seisonidea within the clade. The monophyly/ non-monophyly of eurotatorians (Monogononta and Bdelloidea) relative to Acanthocephala has received much attention because this issue plays a key role in understanding the evolution and ecological diversity found among major groups of Syndermata. Lorenzen [15] first recognized the lemnisci (paired projections of the neck epidermis into the body cavity) and proboscis (an invaginable anterior part of the body) as synapomorphic characters for unifying the Bdelloidea and Acanthocephala (inferred the latter as being the highly specialized sister group to the former; see [24] for different view), and this relationship has also been recovered in many subsequent molecular analyses (SSU+LSU+mtDNA cox1 [26]; SSU [19]; SSU+mtDNA 16S [20]; SSU+LSU+histone H3+mtDNA cox1 [28]). Based on morphological and molecular perspectives, Garey et al. [19] recognized Acanthocephala as a subtaxon of Rotifera and united these two groups under the superclass 'Lemniscea' (Lemniscea hypothesis-eurotatorian paraphyly). However, eurotatorian monophyly (Monogononta+Bdelloidea) and sister-group relationship of Seisonidea and Acanthocephala have also been suggested from some earlier studies based on morphological evidence (ultrastructural similarity [14,29]; cladistic analysis of morphological dataset [7,30,31]) or from molecular analysis (the partial sequence of nuclear heat-shock protein hsp82 [8]; SSU rDNA [32]; combined analysis of hsp82+SSU data [33]). Furthermore, employing different methods for phylogenetic analysis resulted in inconsistent tree topologies when SSU data were analyzed [34]. Although the most recently published work from ESTbased phylogenomic analysis supported Eurotatoria paraphyly [27], but the phylogenetic issue regarding the Eurotatoria monophyly/paraphyly still awaits other types of data for corroboration.
With a very few exceptions, metazoan mitochondrial genomes are circular DNA molecules (mostly less than 16 kb in size) that encode 37 genes: 13 protein-coding genes (atp8 is missing in many nematode and flatworm species so far reported), two ribosomal RNA genes, and 22 transfer RNA genes [35,36]. Due to its universality and remarkably stable feature in genome content across various metazoan phyla, comparisons of the mitochondrial genome information (e.g., nucleotide sequence, amino acid sequence and gene order rearrangement) have often proven useful for reconstructing the deep node phylogeny and for assessing the phylogenetic relationships among closely related species [36][37][38][39]. In recent years, there has been an unprecedented increase in mitochondrial genomic surveys in relation to phylogenetic comparisons of a variety of animal groups. Mitochondrial genome information has now become available for more than a thousand animal species (See NCBI metazoan mitochondrial genome resources). However, the distribution of completely characterized mitochondrial genome sequences has been strongly biased across the metazoan taxa: the subphylum Vertebrata and the phylum Arthropoda account for more than 80% of the metazoan mitochondrial genome data determined so far, whereas there are still a considerable number of metazoan phyla for which there is limited mitochondrial genome information, and some phyla have never been investigated [40,41]. Complete mitochondrial genome sequencing from poorly investigated groups is needed to supplement a gap in our current understanding of mitochondrial DNA evolution in the metazoa.
Mitochondrial genome information from the Syndermata was reported, for the first time, from the acanthocephalan species Leptorhynchoides thecatus [42] and recently thereafter from the monogonont rotifer Brachionus plicatilis [43]. Despite belonging to the same clade (Syndermata), the mitochondrial genomes reported for these two organisms did not share many characteristics in their organization such as chromosome structure, gene order, codon usage and the secondary structure of tRNA molecules. For example, the mitochondrial genome of B. plicatilis, a representative of the Monogononta, is encoded in two separate mitochondrial chromosomes, each having different gene content and copy number. In contrast, as in most metazonas, all genes of L. thecatus mtDNA are encoded in a single type of circular mitochondrial DNA molecule. The lack of common features in the mitochondrial genomes between these species necessitates additional characterization of mitochondrial genomes from the other major groups of the syndermata. Information from Bdelloidea is expected to supplement our understanding of syndermatan mitochondrial genome evolution and provide utility as a molecular marker in resolving internal phylogeny among the major groups of the Syndermata. To this end, we characterized the first complete mitochondrial genome sequence of the bdelloid species Rotaria rotatoria, and compared its mitochondrial genome information with other syndermatan species in order to investigate phylogenetic issues regarding the monophyly or paraphyly of Eurotatoria.

General features of the R. rotatoria mitochondrial genome
The mitochondrial genome of R. rotatoria is 15,319 bp in length (GenBank accession no. GQ304898) and contains 12 protein-coding genes (lacking atp8), two rRNA genes, and as is occasionally found in some other metazoans, all but one (trnC) of the typical 22 tRNAs genes. All genes are encoded in the same direction in a single circular mitochondrial DNA (Fig. 1), unlike that of the monogonont rotifer B. plicatilis, in which mitochondrial genes are encoded in two separated circular chromosomes [43]. The genome organization including the gene order, length of overlapping regions between the genes and length of intergenic spacer regions is shown in Table 1. The nucleotide composition of the entire R. rotatoria mtDNA sequence is 29.5% A, 43.7% T, 17.4% G, and 9.5% C ( Table 2). The overall A+T content (73.1%) approximates more or less that of the acanthocephalan L. thecatus (71.5% in A+T content), but is significantly higher than those found in the monogonont rotifer B. plicatilis (A+T contents of 63.9% and 62.9% for mtDNA-I and mtDNA-II, respectively). The average of A+T composition of protein-coding genes (73.7%) is very similar to that of the entire sequence, but there is noticeable compositional difference across the codon positions: the average A+T content of the third position (80.9%) is noticeably higher than those of the first (67.8%) and second codon (72.4%) positions, indicating that there is a strong bias toward T (47.2%) and A (33.7%) in third codon position.

Protein-coding genes
Twelve protein-coding genes were identified using NCBI ORF Finder and by comparing their sequence with those of homologous genes reported from the acanthocephalan L. thecatus [42] and monogonont B. plicatilis [43]. Even after an exhaustive search using BLAST, we failed to discover any atp8-like protein sequence. The lack of atp8 in the mtDNA genome is not very rare, and quite common in most nematode and platyhelminth species reported so far (cf. Trichinella spiralis, a nematode species where the atp8 exists [44]).
In many metazoan mtDNAs, the codon usage is strongly biased, especially in the third position of synonymous codons within amino acid families [35,[45][46][47]. In the protein-coding genes of R. rotatoria mtDNA, the higher A+T content is related to amino acid composition with higher abundance of A+T-rich codons. Amino acids encoded by T-rich and A-rich codons (more than two Ts and As in a triplet, respectively) account for 47.3% and 22.2% of the total amino acid composition, respectively, totaling 69.5% of the entire protein sequences (Table 3). Another factor that causes nucleotide bias is differential preference of thymine (T) or adenine (A) in the third position of codons. The average relative nucleotide frequency in the third position is 47.2% (T), 33.7% (A), 13.7% (G), and 5.4% (C), respectively. Indeed, both in twofold-and fourfold-degenerate codon families, the T-ending codons are greatly favored in most synonymous codons, whereas Cending codons are evidently avoided. The relative frequency of TTT (9.5%) coding phenylalanine is much higher than that of its degenerate codon TTC (1.5%). Another notable example of biased codon usage is that the codon ATT (8.0%) is much more frequently used than ATC (0.7%) in the isoleucine family. The propensity of a preference for T, with an apparent bias against C may reflect mutational bias at the third codon position as postulated in other metazoans [48]. Of 12 protein-coding genes, six (atp6, cob, nad4L, cox2, nad3, and nad2) use ATG as the initiation codon, whereas others initiate with ATT (cox1, nad6, nad1, and nad5), or ATA (nad4 and cox3), respectively ( Table 1). The majority of genes end with complete termination codon TAA (cox1, nad1, cob, nad4L, nad4, nad5, cox3, and nad3) or TAG (nad6 and atp6), but cox2 and nad2 terminate with the incomplete stop codon T, which is presumed to become a complete termination codon through post-transcriptional polyadenylation [49].

Transfer RNA and ribosomal RNA genes
Implementation of tRNAscan algorithm failed to find any tRNA-like secondary structure, but we manually identified them by eye. Of 22 tRNA genes normally found in most other metazoan mtDNAs, 21 tRNA-like nucleotide segments (ranging from 47 to 65 bp in size) except for trnC can be folded into a cloverleaf secondary structure with some mismatches or incomplete configuration of the DHU and/or TΨC arms (Additional file 1). In spite of the exhaustive tRNA search, we were not able to identify any candidates for trnC with high confidence. Out of 21 tRNAs found, only 10 (trnA, trnL1, trnM, trnF, trnP, trnS1, trnS2, trnT, trnY, and trnV) display a typical cloverleaf-like structure equipped with DHU and TΨC arms. Many other tRNAs appear to lack a TΨC arm (6 tRNAs; trnD, trnQ, trnE, trnG, trnH, and trnK) or a DHU arm (4 tRNAs; trnN, trnI, trnL2, and trnW), or both (1 tRNAs; trnR). The rarity of having a typical cloverleaf-like secondary structure was also found in the acanthocephalan species Leptorhynchoides thecatus [42] where almost none of the inferred tRNAs displayed a typical stem-and-loop configuration in the DHU and TΨC arms. This differs from most other metazoan mtDNAs, including the monogonont rotifer Brachionus plicatilis [43]. Indeed, the putative secondary structure of 22 tRNAs found in B. plicatilis contains a typical cloverleaf structure comprising a amino-acyl arm (a stem of seven nucleotide pairs; ntp), a DHU arm (a stem of 2-4 ntp with a 2-8 nt loop), an anticodon arm (a stem of 5 ntp with an anticodon) and a TΨC arm (a stem of 2-4 ntp with a 2-7 nt loop).
Based on sequence comparison with L. thecatus and B. plicatilis, two ribosomal RNA genes were identified: The rrnL (529 bp) is located between trnY and trnL1 as found in L. thecatus, but its size is remarkably smaller than any other metazoan large ribosomal subunit molecule (mostly larger than 1 kb) including those of other syndermatans reported thus far (e.g., 925 bp and 1,107 bp in L. thecatus and B. plicatilis, respectively). The rrnS (521 bp) is located between trnM and cox2. and trnV], and 236 bp-NCR5 [between trnV and trnS1, respectively). Despite an exhaustive search using the NCBI ORF Finder, we were not able to find ORF-like candidates of significant length from these spacer regions. The relatively larger genome size of the R. rotatoria mtDNA is attributed to the abundance of conspicuously long, unassigned spacer regions, compared to that of acanthocephalan L. thecatus (13,888 bp [42]). The region located between trnR and trnS2 (NCR3) contains 24 tandemly repeated units of a 18-nt sequences (5'-NRRWYBTYRNG-HRRYYYY-3'), each having potential to be folded into a stem-loop structure (a stem of 5-6 ntp with a 6-8 nt loop). Such repeats, known as variable numbers of tandem repeats (VNTRs), have been reported from a variety of animal mtDNAs and they have been utilized as molecular Circular representation of the mitochondrial genome of Rotaria rotatoria Figure 1 Circular representation of the mitochondrial genome of Rotaria rotatoria. All genes are encoded in the same direction and 21 tRNA genes are denoted by the oneletter abbreviation and two leucine and two serine tRNA genes are labeled, according to their anticodon sequence, as L1 (trnL-uag), L2 (trnL-uaa), S1 (trnS-ucu), and S2 (trnS-uga), respectively. The intergenic non-coding regions with notable size (>200 bp) are denoted as NCR1-NCR5, respectively.

Mitochondrial molecular phylogeny Phylogenetic position of Platyzoa within Lophotrochozoa
The Platyzoa refers to an assemblage of some acoelomate and pseudocoelomate animal taxa including Platyhelminthes, Rotifera, Acanthocephala, Gastrotricha, and Gnathostomulida [51]. Phylogenetic affinity of these groups with other lophotrochozoans has not yet been confirmed. Results of recent molecular analyses depicted Platyzoa as a sister to Trochozoa (morphology+SSU dataset [21]), or as a derived subclade of Lophotrochozoa (LSU+SSU dataset [23]. This disagreement still awaits unambiguous resolution. In this study, we performed phylogenetic analyses for the amino acid dataset of major representatives of Bilateria with two cnidarian sequences (Aurelia aurita and Acropora tenuis) as outgroups. Resulting trees from both ML and BI methods showed that the protostomes and deuterostomes were divided into two separate clades with very high support (Fig. 2). Within the protostome clade, monophyly of Ecdysozoa and monophyly of Lophotrochozoa were also found, but the dichotomy between these two groups was collapsed when the nematode sequences (Caenorhabditis elegans, Trichinella spiralis) were included in the analysis. In all cases, nematode sequences were nested within the Lophotrochozoan clade rather than grouped with ecdysozoan members (Additional file 2), and their position varied considerably according to different phylogenetic methods. Their unusual position appears to be artificial and due to the longbranch attraction (LBA) problem as documented in many previous papers [52][53][54]. Therefore, in order to reduce this factor in the tree reconstruction, we excluded nematode sequences from subsequent analyses. After the exclusion of nematodes, the concatenated amino acid sequence dataset consisting of 1,556 homologous positions was used for the subsequent analyses. Note that monophyly of Lophotrochozoa and monophyly of Ecdysozoa were recovered after removal of nematode sequences from the analyses, and support values for clades were noticeably improved in all analyses (see table description of Additional file 3 for the phylogenetic position of Chaetognatha within the metozoa). The monophyletic grouping of ecdysozoan members received 100% BP (bootstrap percentage) in ML, and 1.00 BPP (Bayesian posterior probability). The Lophotrochozoan clade was also well supported with 1.00 BPP in BI, but received a very weak bootstrap support value of 56% in ML analysis (Fig. 2). Although the Lophotrochozoan clade is not strongly supported, the current amino acid sequence data support the Ecdysozoa/Lophotrochozoa split. It is particularly noteworthy that both phylogenetic methods (BI and ML) recovered Platyhelminthes, Rotifera, and Acanthocephala (all platyzoan members included in this study) as monophyletic lineages in the Lophotrochozoan clade and nodal support for Platyzoa was very strong (100% in both BI and ML analysis). This result is in agreement with an earlier analysis based on a combined dataset of LSU and SSU sequences that recognized the Platyzoa as a derived monophyletic assemblage within the Lophotrochozoa [23]. As yet, the availability of complete mitochondrial genome information is limited to Platyhelminthes, Acanthocephala, and Rotifera, and further evidence from Gnathostomulida and Gastrotricha is required to confirm the monophyly of Platyzoa and its derived position within the Lophotrochozon clade (sensu [51]).

Phylogenetic implications for Syndermata: Eurotatoria paraphyly
Phylogenetic relationships among Syndermata, in particular with regard to Eurotatoria monophyly/paraphyly are still vigorously debated. The answer to this question is prerequisite for correct understanding of the evolution of different types of life styles found among major groups of Syndermata [27]. In this study, we revisited this issue based on phylogenetic analysis of amino acid sequence for eleven protein-coding genes (except atp6 and atp8) of the complete mitochondrial genome. The resulting trees both from Bayesian and maximum likelihood analyses recovered syndermatan members (Monogononta, Bdelloidea, and Acanthocephala) as a monophyletic group, which is in turn grouped with Platyhelminthes (Fig. 2). The phylogenetic affinity of the Syndermata with Platyhelminthes in this study is concordant with some earlier publications asserting their sister-relationship [23,55]. It is particularly significant to note that both BI and ML methods placed Bdelloidea and Acanthocephala as sister taxa and this relationship was very strongly supported by 1.00 BPP in BI analysis, but received relatively weaker support (65% BP) in ML analysis. Moreover, none of the phylogenetic methods recovered the two eurotatorian species, Rotaria rotatoria (Bdelloidea) and Brachionus plicatilis (Monogononta) as monophyletic, but instead always depicted them as a paraphyletic assemblage. Bdelloidea represented by R. rotatoria appears to share most recent common ancestry with Acanthocephala rather than with the monogonont B. plicatilis. The sister-group position of the Monogononta to the Bdelloidea+Acanthocephala, received strong support from BI (0.99 BPP) and moderate support from ML (70% BP), respectively. This result indicates a firm support for the 'Eurotatorian paraphyly' hypothesis. In all cases of the analysis, monophyly of Bdelloidea + Acanthocephala relative to Monogononta is always very robust. In order to assess whether our data rejects the alternative hypothesis (Eurotatoria monophyly), we performed statistical comparison of the likelihood scores for the best ML tree without constraint versus the best tree with the topological constraint of Eurotatoria monophyly using theTreefinder program. The result of topology test failed to detect a significant difference between these two competing hypotheses under various criteria (ELW, SH and AU; Table 4). Although Eurotatoria monophyly was not significantly worse interpretation of these data by tree topology comparison test, the optimal phylogenetic trees from different methods unanimously depicted B. plicatilis (Monogononta) sister to two other syndermatans (bdelloidean R. rotatoria and acanthocephalan L. thecatus), regardless of different options applied in the data analysis (including or excluding inferred gaps in the analysis). In addition, this relation (Eurotatoria paraphyly) is further corroborated by sharing greater similarity in gene arrangement between R. rotatoria and L. thecatus (next section).
A close relationship between Bdelloidea and Acanthocephala is also evident when gene arrangement patterns of R. rotatoria mtDNA are compared with other metazoan animal groups. Gene arrangement of R. rotatoria mtDNA is very different from those found in other metazoans, but shows the highest similarity to the acanthocephalan L. thecatus. The gene order comparison between R. rotatoria and L. thecatus identified the gene order of trnN-nad2-cox1-trnG, trnY-rrnL-trnL1-nad6, and nad4L-nad4-trnH-nad5 shared between the representatives of the Bdelloidea and Acanthocephala, respectively (Fig.  3). The gene order of nad4L-nad4-trnH-nad5 can also be found across many various metazoan mtDNAs, such as Bryozoa [56], Brachiura (Arthropoda [57]), Gastropoda (Mollusca [58]), Vertebrata (Teleostei [59] and Amphibia [60]), and Priapulida [54]. However, this conserved gene order is not shared by the B. plicatilis mtDNA and we suggest it is due to the idiosyncratic features of B. plicatilis mtDNA in which 36 genes are encoded separately into bipartite mitochondrial genomes [43]. In contrast, note that the gene order similarity found between R. rotatoria and L. thecatus (trnN-nad2-cox1-trnG and trnY-rrnL-trnL1-nad6) is very unique among the metazoans and this shared gene rearrangement lends another line of strong support for their common ancestry. The phylogenetic conclusion of Eurotatoria paraphyly suggested from the present study is also consistent with the most recent EST based phylogenomic analysis [27].
In general, the evolutionary potential for adaptive radiation and its subsequent diversification of extant biota are often associated with preferential adoption of different life styles during evolutionary history. In most cases, a robust phylogeny is a prerequisite to correctly understand the evolution of life styles among the groups. The Syndermata (Rotifera+Acanthocephala) is one of the best model systems for studying the evolutionary origin and persistence of different life styles because it shows a variety of lineage-specific reproduction modes. Of these, Bdelloidea is particularly among the most unique in that it contains exclusively several hundred species and all individuals are females, the largest metazoan group with no sexual reproduction that has ever been documented [9]. This group is also considered to represent an ancient origin of asexuality with a great evolutionary success in their diversifica-

-ln -ln ELW SH AU
Best ML tree 56700.07 Monophyly of Eurotatoria 56724. 49 24.42 0.199 0.197 0.196 tion of asexual species [9,10]. There are some additional cytological [61] and molecular evidence [62] suggesting that the evolution and maintenance of asexuality in Bdelloidea is of relatively ancient origin, dating back to more than tens of millions of years ago. Inferring the history of the obligatory asexual class Bdelloidea within the syndermatan phylogenetic framework is a central topic in evolutionary biology of asexuality because these organisms have long been regarded as perfect candidates to test the advantage of sexual reproduction [63]. The phylogenetic conclusion obtained from this comparative mitochondrial genome study provides essential clues in extrapolating the directionality of long-term asexual evolution of bdelloidean lineage. The sister group position of Monogononta to the Bdelloidea/Acanthocephala clade suggests that the obligatory asexuality of bdelloideans may have been derived secondarily from some other preexisting condition in earlier lineage of rotifers. However, because no mitochondrial genome information is available from Seisonidea, we are not able to postulate what is the most likely ancestral condition that gave rise to bdelloidean asexuality. There are only a few complete mitochondrial DNA sequences available from the Syndermata: B. plicatilis (Monogononta [43]), R. rotatoria (Bdelloidea; this study) and L. thecatus (Acanthocephala [42]). More information, especially from Seisonidea needs to be obtained to confirm the phylogenetic position of Seisonidea within Syndermata. This further phylolgenetic information would be a very substantial ingredient for precisely estimating the origin of asexuality in Bdelloidea and at the same time, contribute to a better understanding of mitochondrial genome evolution within the syndermatan phylogenetic framework.

Conclusion
In this study, we revisited phylogenetic relationships of Eurotatoria (monophyly versus paraphyly) based on comparative analysis of the complete mitochondrial genome information for major groups of the Syndermata. For this purpose, we determined the first complete mitochondrial genome sequence (15,319 bp) of a bdelloid rotifer, Rotaria rotatoria. Comparisons of inferred amino acid sequence and gene arrangement pattern with those of other syndermatan mtDNAs (from the acanthocephalan L. thecatus and monogonont B. plicatilis) supports the hypothesis that the Bdelloidea shares a most recent common ancestor with Acanthocephala, and not with Monogononta. From this finding, we suggest that the obligatory asexuality of bdelloideans may have secondarily derived from some preexisting condition in earlier lineages of rotifers. Definitively, testing this question of phylogenetic relationships and evolution of different types of life style among the major members of Syndermata awaits further investigation of mitochondrial genome sequencing from Seisonidea.

Sampling and molecular techniques
Rotifers were collected from a small pond at Inha University campus in South Korea (37°26'58.08"N, Linearized comparison of the mitochondrial gene arrangement of three syndermatan species Figure 3 Linearized comparison of the mitochondrial gene arrangement of three syndermatan species. Gene and genome size are not to scale. All genes are transcribed in the same direction (from left to right). The tRNAs are labeled by single-letter abbreviations. Gene clusters shared between L. thecatus (Acanthocephala) and R. rotatoria (Bdelloidea) are represented by shadowed areas.  Table 5). PCR reactions were conducted in a total volume of 50 μl reaction mixture containing 0.13 μg/μl of genomic DNA, 10× PCR buffer, 10 mM dNTP mixture, 10 pmole each primer, 25 mM MgCl 2 , and 2.5 units Taq polymerase (Roche Co.) using the following cycling conditions: 1 cycle (92°C for 2 min), 35 cycles (92°C for 1 min, 45°C to 55°C for 30 sec, 72°C for 1 min 30 sec), and 1 cycle (72°C for 10 min). The nucleotide sequence determined from these four fragments was then used to design R. rotatoria mtDNA-specific primers for long PCR amplification (see also Table 5 for details of the primer information). Four overlapping fragments (ranging in size from 1,559 bp to 5,481 bp) covering the entire mitochondrial genome of R. rotatoria were amplified using the Expand Long Template PCR System (Roche, USA) under the following conditions: 1 cycle of 2 min at 94°C (initial denaturation), 35 cycles of denaturation-primer annealingelongation (10s at 95°C, 1 min at 55°C, and 10 min at 68°C), and 1 cycle of the final extension (10 min at 68°C). A negative control (no template) was included for every PCR run to detect any potential contamination of the PCR products. The amplified PCR products were isolated on a 1.0% agarose gel containing crystal violet, excised and extracted according to the TOPO XL gel-purifying protocol (Invitrogen Co.). Purified DNAs were cloned into E. coli competent cells using TOPO XL PCR Cloning kit, as recommended by manufacturer. Sequenc-ing reactions of the target fragments were performed in both directions by 'primer walking' using a Big Dye Terminator Cycle-Sequencing Kit (Applied Biosystems), and overlapping fragments were assembled to complete the sequence of the entire genome.

Gene Annotation
Twelve protein-coding genes and two ribosomal RNA genes of R. rotatoria were identified by sequence comparison with the mtDNA sequences of the acanthocephalan L. thecatus and monogonont rotifer B. plicatilis and by using the NCBI ORF (open reading frame) Finder. We attempted to find tRNAs using tRNAscan-SE 1.21, but no tRNA-like structure was detected from this search. Therefore, tRNA genes were identified by searching for anticodon consensus motif sequences (TxxxR; xxx = anticodon) and by recognizing potential secondary structures by eye. In this process, we were aided by a web-based automatic annotation program for organellar genomes (DOGMA [64]).

Phylogenetic analysis
Eleven mitochondrial protein-coding genes (except atp6 and atp8, which are highly variable in their lengths among the groups) from 35 species (6 ecdysozoans, 22 lophotrochozoans, 5 deuterostomes, and two cnidarian outgroups) representing major clades of the Bilateria, including those of R. rotatoria mtDNA, were used in phylogenetic analysis (see Additional file 3 for details of taxon sampling). A multiple alignment of the amino acid sequences for each protein-coding gene was performed using ClustalX [65] with default options. Implementation of ClustalX for multiple sequence alignment does not always guarantee an unambiguous result due to the length and sequence variation among the species. This becomes more problematic when taxon sampling includes a wide range of animal taxa that evolve at very different rates. To improve reliability, we selected conserved blocks from aligned amino acid sequences for each of protein-coding genes using the Gblocks program [66] and a concatenated dataset was then prepared for the subsequent phylogenetic analyses. To reconstruct mitochondrial gene phylogeny, maximum likelihood (ML) and Bayesian inference (BI) were conducted for the concatenated amino acid sequence dataset. For ML analysis, the best-fit model of our amino acid sequence datasets was estimated using the Akaike Information Criterion (AIC) using ProtTest version 2.0 [67]. Maximum likelihood analysis was performed in Treefinder October version [68] using the MtArt matrix [69], which was selected as the best-fit model of amino acid substitution from ProtTest. Nodal support of the resulting ML trees was estimated by nonparametric bootstrap analysis with 500 random replications using Treefinder. We compared the likelihood scores of the competing hypotheses (the best tree versus alternative hypothesis) using various criteria (Expected-Likelihood Weight, ELW [70]; Shimodaira-Hasegawa [71]; Approximately Unbiased, AU [72] tests) implemented in Treefinder. The current version of the MrBayes program does not include the MtArt model of protein sequence, and thus we used the MtRev model with the likelihood parameter setting to "ngammacat = 4", "rates = invgamma" as an alternative best-fit model for Bayesian analysis. The analysis was run for 10 6 generations, sampled every 100 generations with four Markov Chain Monte Carlo (MCMC) chains using MrBayes 3.1.2 [73]. Bayesian posterior probability (BPP) values were estimated after the initial 200 saved trees (the first 2 × 10 5 generations) were discarded as burn-in. We also conducted maximum parsimony (MP) and neighbour-joining (NJ) analyses and nodal support in the resulting tree was estimated by nonparametric bootstrap analysis with 1,000 random replications in PAUP* 4.0b10 version [74].