- Open Access
Comparative mitogenomics supports synonymy of the genera Ligula and Digramma (Cestoda: Diphyllobothriidae)
Parasites & Vectorsvolume 11, Article number: 324 (2018)
After observing differences in the number of reproductive complexes per proglottid within the genus Ligula, the genus Digramma was erected. However, the validity of Digramma has been previously questioned due to a low variability in the cox1, nad1 and ITS rDNA sequences between the two genera. We undertook a study to greatly increase the amount of sequence data available for resolution of this question by sequencing and characterizing the complete mitogenomes of Digramma interrupta and Ligula intestinalis.
The circular mtDNA molecules of Digramma interrupta and Ligula intestinalis are 13,685 bp and 13,621 bp in size, respectively, both comprising 12 PCGs, 22 tRNA genes, two rRNA genes, and two mNCRs. Both mitogenomes exhibit the same gene order and share 92.7% nucleotide identity, compared with 85.8–86.5% to the most closely related genus Dibothriocephalus. Each gene from D. interrupta and L. intestinalis is almost of the same size, and the sequence identity ranges from 87.5% (trnD) to 100% (trnH, trnQ and trnV). NCR2 sequences of D. interrupta and L. intestinalis are 249 bp and 183 bp in length, respectively, which contributes to the main difference in length between their complete mitogenomes. A sliding window analysis of the 12 PCGs and two rRNAs indicated nucleotide diversity to be higher in nad5, nad6, nad2, nad4 and cox3, whereas the most conserved genes were rrnL and rrnS. Lower sequence identity was also found in nad2, nad4, nad5, nad6 and cox3 genes between the two diphyllobothriids. Within the Diphyllobothriidae, phylogenetic analysis indicated Ligula and Digramma to be most closely related to one another, forming a sister group with Dibothriocephalus.
Owing to higher nucleotide diversity, the genes nad2, nad4, nad5, nad6 and cox3 should be considered optimal candidates to use as molecular markers for population genetics and species identification between the two closely related species. The phylogenetic results in combination with the comparative analysis of the two mitogenomes, consistently support the congeneric status of L. intestinalis and D. interrupta.
Based upon its paraphyly, differences in the position of the genital pore, the presence of an external seminal vesicle and the absence of a uterine sac in the Diphyllobothriidea Kuchta, Scholz, Brabec & Bray, 2008, the order Pseudophyllidea van Beneden in Carus, 1863 was suppressed and the Diphyllobothriidea and Bothriocephalidea Kuchta, Scholz, Brabec & Bray, 2008 were proposed [1,2,3].
The order Diphyllobothriidea includes 70 species considered valid, classified into 18 genera across three families [2, 4]. Adult diphyllobothriideans are found only in tetrapods, never having been recorded in fish , and the plerocercoids of groups such as Spirometra, Diphyllobothrium Cobbold, 1858 (syn. Diplogonoporus Lönnberg, 1892) and Dibothriocephalus Lühe, 1899 (a recently resurrected genus including some species from Diphyllobothrium, i.e. Dib. dendriticus (Nitzsch, 1824), Dib. nihonkaiensis (Yamane, Kamo, Bylund & Wikgren, 1986) and Dib. latus (Linnaeus, 1758) belonging to the family Diphyllobothriidae are often the principal agents of food-borne cestodosis . Ligula spp. also belong to the Diphyllobothriidae; these species use copepods as first intermediate hosts, freshwater fish as second intermediate hosts and birds as definitive hosts . Ligula intestinalis (Linnaeus, 1758) is a tapeworm of veterinary importance worldwide that reduces the fecundity of the cyprinid fishes by parasitic castration [6, 7] and induces mass mortalities of the carp Chanodichthys erythropterus .
Cholodkovsky  erected the genus Digramma Cholodkovsky, 1915 after observing differences in the number of reproductive complexes contained within each proglottid when studying the genus Ligula Linnaeus, 1758. In China, Ligula spp. are distributed in the Qinghai-Tibet Plateau with Schizothoracinae fishes serving as the primary second intermediate host. Digramma spp. are found across the rest of China where the goldfish Carassius auratus acts as their common second intermediate host . However, the validity of Digramma has been questioned due to a low level of difference between the species of two genera in the cox1, nad1 and ITS rDNA sequences [11,12,13]. Thus Digramma is considered to be synonymous with Ligula . Owing to the fact that only one gene and a limited number of isolates were included in that study, more sequence data and a greater range of taxa from different genera are required for a more comprehensive phylogenetic analysis of the Diphyllobothriidea .
This study, therefore, aimed to sequence and characterize the complete mitogenomes of Digramma interrupta and Ligula intestinalis and to perform phylogenetic analysis to investigate whether or not these two diphyllobothriids are congeneric using mitogenomic data. Differences within the mitochondrial genes were also compared to determine which genes would be suitable for the design of molecular markers as a means to differentiate D. interrupta from L. intestinalis.
Specimen collection and DNA extraction
Plerocercoids of D. interrupta and L. intestinalis were isolated from the body cavity of Carassius auratus collected from Liangzi Lake in Hubei Province and Gymnocypris selincuoensis from Siling Lake, Tibet, China, respectively. The tapeworms were preserved in 80% ethanol and stored at 4 °C. Total genomic DNA was extracted from the posterior region of a single tapeworm using a TIANamp Micro DNA Kit (Tiangen Biotech, Beijing, China), according to the manufacturer’s instructions. DNA was stored at -20 °C for subsequent molecular analysis. The morphological identification of specimens was confirmed using sequence data from the ITS2 rDNA region [13, 14] and the cox1 gene .
Amplification and DNA sequencing
PCR was carried out as described previously [15, 16], with minor modifications. Five degenerate primer sets (Additional file 1: Table S1) were designed to primarily amplify partial sequences of the nad5, cytb, nad2, cox1 and rrnS genes. The sequenced fragments were subsequently used to design primers specific for the amplification and sequencing of the whole mitogenome. PCR reactions were performed in a 20 μl reaction mixture, containing 7.4 μl dd PCR grade H2O, 10 μl 2× PCR buffer (2 mM Mg2+, 8 μl dNTP plus, Takara, Dalian, China), 0.6 μl of each primer (12.5 μM), 0.4 μl rTaq polymerase (250 U, Takara, Dalian, China), and 1 μl genomic DNA template. Amplification was conducted under the following conditions: initial denaturation at 98 °C for 2 min, followed by 40 cycles at 98 °C for 10 s, 48–60 °C for 15 s, 68 °C for 1 min/kb, and a final extension at 68 °C for 10 min. PCR products were sequenced on an ABI 3730 automatic sequencer using the Sanger method at Sangon Company (Shanghai, China) using the primer walking strategy.
Sequence annotation and analyses
The mitogenome was annotated broadly following the procedure described previously [15, 16]. The amplified fragments were initially checked by BLASTN , before being assembled manually in a stepwise manner. The annotation was recorded in a Word document with the help of the Geneious program , using the mitogenome of Dibothriocephalus latus (syn. Diphyllobothrium latum) (NC_008945) as the reference sequence. PCGs were found by searching for ORFs (employing genetic code 9, echinoderm mitochondrial; flatworm mitochondrial) and the nucleotide alignments against the selected reference genome in Geneious. rrnL and rrnS were annotated in a similar way, via comparison with homologs using Geneious. ARWEN  and MITOS  web servers were employed to identify and fold all tRNAs. Similarly, the NCBI submission file (*.sqn) and tables of statistics for mitogenomes were generated using a home-made GUI-based program, MitoTool . MitoTool was also used to calculate codon usage and relative synonymous codon usage (RSCU) for the 12 PCGs of D. interrupta and L. intestinalis and the results were sorted and imported into ggplot2  to draw the RSCU figure. A Tandem Repeats Finder  was used to predict tandem repeats (TR) in the major non-coding regions (mNCRs), and the secondary structures of NCR1 and TR were folded by Mfold software . Non-synonymous (dN) and synonymous (dS) mutation rates among the 12 PCGs of D. interrupta and L. intestinalis were computed using KaKs_Calculator  utilising a modified Yang-Nielsen algorithm. DnaSP v5  was adopted to conduct sliding window analyses. A sliding window of 500 bp and step size of 25 bp was implemented to estimate nucleotide divergence Pi (π) between the alignments of the mitogenomes of D. interrupta and L. intestinalis.
The mitogenomes of 35 cestodes, covering five orders and ten families (Additional file 2: Table S2) were obtained from GenBank and were used, along with the two new mitogenomes generated in this study, to create the phylogenetic reconstruction. Two trematodes, Dicrocoelium chinensis (NC_025279) and Dicrocoelium dendriticum (NC_025280), were used as outgroups. All 36 genes (12 PCGs, 2 rRNAs and 22 tRNAs) were used for phylogenetic inference and were extracted from GenBank files using MitoTool. PCGs were aligned in batches using MAFFT and integrated into our own in-house GUI-based program, BioSuite , adopting codon-alignment mode. All RNA genes (rRNA and tRNA) were aligned using a structural alignment algorithm Q-INS-i incorporated into MAFFT-with-extensions software . Gaps and ambiguous sites were deleted using GBlocks  integrated by BioSuite with default settings. BioSuite was subsequently used to concatenate the sequences into a single alignment and generate phylip and nexus format files.
GTR+I+G was chosen as the optimal model of nucleotide evolution for all datasets based on the Akaike information criterion by ModelGenerator . Two analytical methods were performed: maximum likelihood (ML) and Bayesian inference (BI). The ML analysis was performed in RAxML GUI  using a ML+rapid bootstrap (BP) algorithm with 1000 replicates. BI analysis was performed in MrBayes 3.2.6  with default settings, 6 × 106 metropolis-coupled MCMC generations, and 1000 sample frequency. Stationarity was considered to have been reached when the average standard deviation of split frequencies was below 0.01, ESS (estimated sample size) was above 200, and PSRF (potential scale reduction factor) approached 1. Bayesian posterior probability (BPP) values were calculated in a consensus tree, after discarding the first 25% of samples as burn-in. Finally, the resultant trees were visualised and annotated by iTOL  with the help of several dataset files generated by MitoTool.
Genome organization and base composition
The length of the circular mtDNA molecules of D. interrupta (GenBank accession number: MF671697) and L. intestinalis (GenBank accession number: MF671696) was 13,685 bp and 13,621 bp, respectively. Both mitogenomes were composed of 12 PCGs, 22 tRNA genes, two rRNA genes and two mNCRs (major non-coding regions), all of which were transcribed from the same strand (Fig. 1). As commonly reported for flatworms , the two mitogenomes lacked the atp8 gene. The gene order of the two mitogenomes was identical, conforming to the synapomorphic gene arrangement of the order Diphyllobothriidea . The A+T content of the mitogenomes of D. interrupta and L. intestinalis were 67.9% and 67%, respectively, which is in accordance with that of other cestodes (Additional file 2: Table S2).
The mitogenomes of D. interrupta and L. intestinalis shared 92.7% nucleotide identity (Table 1), compared with 85.8% and 86.2% to Dib. latus, 86.1% and 86.5% to Dib. nihonkaiensis, species of the most closely related genus, respectively. Nucleotide identity was 92.3% between Dib. latus and Dib. nihonkaiensis, 99.2% between Diphyllobothrium grandis (Blanchard, 1894) and D. balaenopterae (Lönnberg, 1892), 99.3% between Spirometra decipiens (Diesing, 1850) and S. erinaceieuropaei (Rudolphi, 1819), respectively. Amongst all 36 genes, the majority were equal in size between D. interrupta and L. intestinalis, with the exception of rrnS, trnY and trnS2 which had only one base difference. Sequence identity ranged from 87.5% (trnD) to 100% (trnH, trnQ and trnV) (Table 1). Seven overlapping regions and 20 intergenic sequences (Gap1–20) were found in both genomes, identical in size and position, with the exception of GAP9, GAP11 and GAP17, which differed in size (Table 1).
Protein-coding genes and codon usage
Concatenated PCGs of the mitogenome of D. interrupta and L. intestinalis were both 10,062 bp in size, with an A+T content of 67.5% and 66.6%, respectively (Table 2). The high A+T content was mainly concentrated on the third codon position (72.4% for D. interrupta and 70.1% for L. intestinalis). The start and termination codons within the mitogenome of D. interrupta and L. intestinalis were identical to one another. GTG was identified as the initial codon for cox3, and ATG for the rest of the 12 PCGs (Table 1). For each PCG, however, all selected Diphyllobothriidea species shared the same start codon (Additional file 3: Table S3), indicating that it may be a synapomorphy within this order. Amongst the termination codons, five out of 12 were identified as TAA, six as TAG, while cox3 used a truncated T stop codon.
Codon usage, RSCU, and codon family proportion (corresponding to the amino acid usage) among D. interrupta and L. intestinalis were investigated (Additional file 4: Figure S1). The most abundant codon families were Phe, Leu2, and Ile within the two mitogenomes, which show a preference for the A+T-rich synonymous codons (Additional file 4: Figure S1). This corresponds to the high A+T bias of the two diphyllobothriid mitogenomes.
We measured the selective pressure acted upon amino acid replacement mutations by the ratio of non-synonymous (dN) to synonymous (dS) substitutions for all 12 PCGs of D. interrupta vs L. intestinalis. Although the values (dN/dS) of atp6 (0.113), nad5 (0.111) and nad2 (0.110) genes were higher than cox1 (0.007) and cox2 (0.008) genes (Fig. 2a), all PCGs were under strong negative (purifying) selection (dN/dS < 0.12).
Transfer and ribosomal RNA genes
The sizes of rrnL in both the D. interrupta and L. intestinalis mitogenome was 967 bp, with 68.4% and 67.9% A+T content, respectively. Similarly, the lengths of rrnS were 742 and 743 bp, with an A+T content of 68.2% and 67.5%, respectively (Table 2). All 22 commonly found tRNAs were present in the mitochondrial genome of D. interrupta and L. intestinalis, adding up to a 1,410 bp total concatenated length in both mitogenomes (Table 2). All tRNAs could be folded into the conventional cloverleaf structure, with the exception of trnS1(AGN) and trnR, which lacked DHU arms. The absence of DHU-arms in trnS1(AGN) and trnR has also been reported in the Caryophyllidea  and the Anoplocephalidae .
The two major non-coding regions (mNCRs), NCR1 and NCR2, were located between trnY and trnL1 and between nad5 and trnG, respectively. The mNCRs were situated in the same location as all diphyllobothriideans surveyed to date (see Additional file 3 in our recent paper ). The NCR1 sequences of the two mitogenomes of D. interrupta and L. intestinalis were 223 and 224 bp in length with a heightened A+T bias of 78.9% and 75%, whereas the NCR2 sequences were 249 and 183 bp in size with 73.9% and 74.9% A+T content, respectively (Table 2). The NCR2 of D. interrupta contained six TRs (tandem repeats). Repeat units 1–5 were identical in nucleotide composition and size (34 bp). Repeat unit 6 was truncated with 29 bp (Fig. 3). TRs were also found in the NCR2 of L. intestinalis, and the consensus repeat (35 bp) was almost identical to that of D. interrupta, with an insertion of a single nucleotide A at the 17th position. Only four repeat units, however, could be found in the mitogenome of L. intestinalis, which contributed to the main difference in length of the complete mitogenome between the two genera. The last repeat unit was also truncated with 22 bp. Both NCR1 and the consensus repeat sequence in NCR2 of the two mitogenomes were capable of forming stem-loop structures (Fig. 3, predicted by Mfold web server).
Sliding window analyses and nucleotide diversity
A sliding window analysis of the 12 PCGs and two rRNAs of D. interrupta and L. intestinalis indicated the nucleotide diversity Pi (π) to be higher in nad5 (0.099), nad6 (0.088), nad2 (0.088), nad4 (0.087) and cox3 (0.087), whereas the most conserved genes were rrnL (0.036) and rrnS (0.04) (Fig. 2b).
Both methods (BI and ML) produced phylograms with concordant branch topologies, thus only the latter was shown (Fig. 4). The phylogenetic tree indicated the ordinal topology to be Caryophyllidea + (Diphyllobothriidea + (Bothriocephalidea + (Proteocephalidea + Cyclophyllidea))). Within the family Diphyllobothriidae, Ligula and Digramma clustered with maximum nodal support (BP = 100 and BPP = 1), which formed a sister group with the genus Dibothriocephalus. This clade clustered together with Diphyllobothrium, then forming a sister group with Spirometra.
The ordinal topology of Caryophyllidea + (Diphyllobothriidea + (Bothriocephalidea + (Proteocephalidea + Cyclophyllidea))) was consistent with previously identified interordinal relationships of tapeworms, when reconstructed from the dataset of nucleotide or amino acid sequences from partial mitogenomes, large and small subunit rRNA genes, and a combination of the former two . Additionally, this relationship identified between tapeworm groups was congruent with latest mitochondrial phylogenomics . However, a sister-group relationship between the orders Diphyllobothriidea and Bothriocephalidea has been suggested based on mitochondrial phylogenomics [38, 39]. This inconsistency may be due to the different methods of phylogenetics employed.
Within the family Diphyllobothriidae, the phylogenetic relationship of the three genera Spirometra, Diphyllobothrium and Dibothriocephalus was congruent with that of recent studies on the phylogenetics of Eucestoda based on mitogenomes, with the topology of Spirometra + (Dibothriocephalus + Diphyllobothrium) [35, 38,39,40,41]. In the present study, however, Ligula and Digramma were closely related to one another with maximum nodal support, then forming a sister group with the genus Dibothriocephalus. Inferred by the ITS2 rDNA sequence  and the 18S rDNA gene , the complexes of Ligula and Digramma have also been shown to be closely related to Dibothriocephalus. Further studies using the concatenated nucleotide sequences of 18S rDNA + 28S rDNA + rrnL + cox1, have again demonstrated the genus Dibothriocephalus to be the sister group of Ligula . These phylogenetic results suggest that Dibothriocephalus is the most closely related genus to Ligula and Digramma.
However, mitogenome sequence identity between D. interrupta and L. intestinalis is 92.7%, which is much higher than between either of these species and the represented members of Dibothriocephalus (85.8–86.5%). Furthermore, high mitogenome sequence identity was also found between the congeners in Dibothriocephalus (92.3%), Diphyllobothrium (99.2%) and Spirometra (99.3%). These results suggest that sequence differences between D. interrupta and L. intestinalis are of a degree expected between members of the same genus.
In one study, sequence identity of the cox1 and nad1 genes between D. interrupta and L. intestinalis has been shown to be 100% and 92.6% ; however, identity was deemed at 93.5% and 93.0% in the present study, respectively (Table 1). This inconsistency may be due to the use of partial sequence of cox1 and nad1 genes or resulting from the use of formalin-preserved specimens . The gene cox1 is considered to be a useful barcode for metazoans , and widely employed for cestode studies [44,45,46,47]. The two mitochondrial genes cox1 and cytb have also been used to study the population genetic structure of L. intestinalis on a local and global scale . However, a lower sequence identity was found in the nad2, nad4, nad5, nad6 and cox3 genes between D. interrupta and L. intestinalis (90–92%), in comparison to the moderate variation seen between the cox1, cytb and nad1 genes (Fig. 5). Additionally, the relatively looser selection pressure of nad5 (0.111) and nad2 (0.110) may accelerate the accumulation of non-synonymous substitutions, which would increase variation of the two genes . These results suggest that the nad2, nad4, nad5, nad6 and cox3 genes should be considered as optimal candidates for genetic markers to be used for population genetics and species identification studies between the two closely related species, D. interrupta and L. intestinalis.
The complete mitogenomes of Ligula intestinalis and Digramma interrupta were sequenced and characterized. The mitogenomes of these two species show a higher identity to each other than to any species in closely related genera. The two mitogenomes consistently support D. interrupta to be a congeneric species with L. intestinalis. High sequence variation in the nad2, nad4, nad5, nad6 and cox3 genes between the two diphyllobothriids suggest that these five genes should be considered as optimal candidates for genetic markers when studying population genetics or looking to differentiate the two closely related species, D. interrupta and L. intestinalis.
major Non-coding regions
Brabec J, Kuchta R, Scholz T. Paraphyly of the Pseudophyllidea (Platyhelminthes: Cestoda): Circumscription of monophyletic clades based on phylogenetic analysis of ribosomal RNA. Int J Parasitol. 2006;36:1535–41.
Kuchta R, Scholz T, Brabec J, Bray RA. Suppression of the tapeworm order Pseudophyllidea (Platyhelminthes: Eucestoda) and the proposal of two new orders, Bothriocephalidea and Diphyllobothriidea. Int J Parasitol. 2008;38:49–55.
Waeschenbach A, Webster BL, Bray RA, Littlewood DTJ. Added resolution among ordinal level relationships of tapeworms (Platyhelminthes: Cestoda) with complete small and large subunit nuclear ribosomal RNA genes. Mol Phylogenet Evol. 2007;45:311–25.
Kuchta R, Scholz T. Diphyllobothriidea. In: Caira JN, Jensen J, editors. Planetary Biodiversity Inventory (2008–2017): Tapeworms from vertebrate bowels of the earth. The University of Kansas Natural History Museum Special Publication No. 25; 2017.
Waeschenbach A, Brabec J, Scholz T, Littlewood DTJ, Kuchta R. The catholic taste of broad tapeworms-multiple routes to human infection. Int J Parasitol. 2017;47:831–43.
Dubinina MN. Tapeworms (Cestoda, Ligulidae) of the Fauna of the USSR. New Delhi: Amerind Publishing; 1980.
Cowx IG, Rollins D, Tumwebaze R. Effect of Ligula intestinalis on the reproductive capacity of Rastrineobola argentea in Lake Victoria. J Fish Biol. 2008;73:2249–60.
Sohn WM, Na BK, Jung SG, Kim KH. Mass death of predatory carp, Chanodichthys erythropterus, induced by plerocercoids of Ligula intestinalis (Cestoda: Diphyllobothriidae). Korean J Parasitol. 2016;54:363–8.
Cholodkovsky N. Notes helminthologiques. Ann Mus Zool Acad Sci Russ, Petrohrad. 1915;20:164–6.
Liao XH, Liang ZX. Distribution of ligulid tapeworms in China. J Parasitol. 1987;73:36–48.
Li JJ, Liao XH. The taxonomic status of Digramma (Pseudophyllidea: Ligulidae) inferred from DNA sequences. J Parasitol. 2003;89:792–9.
Luo HY, Nie P, Yao WJ, Wang GT, Gao Q. Is the genus Digramma synonymous to the genus Ligula (Cestoda: Pseudophyllidea)? Evidence from ITS and 5' end 28S rDNA sequences. Parasitol Res. 2003;89:419–21.
Logan FJ, Horák A, Štefka J, Aydogdu A, Scholz T. The phylogeny of diphyllobothriid tapeworms (Cestoda: Pseudophyllidea) based on ITS-2 rDNA sequences. Parasitol Res. 2004;94:10–5.
Bouzid W, Štefka J, Hypša V, Lek S, Scholz T, Legal L, et al. Geography and host specificity: two forces behind the genetic structure of the freshwater fish parasite Ligula intestinalis (Cestoda: Diphyllobothriidae). Int J Parasitol. 2008;38:1465–79.
Zhang D, Zou H, Wu SG, Li M, Jakovlic I, Zhang J, et al. Sequencing, characterization and phylogenomics of the complete mitochondrial genome of Dactylogyrus lamellatus (Monogenea: Dactylogyridae). J Helminthol. 2017; https://doi.org/10.1017/S0022149X17000578.
Zhang D, Zou H, Wu SG, Li M, Jakovlić I, Zhang J, et al. Sequencing of the complete mitochondrial genome of a fish-parasitic flatworm Paratetraonchoides inermis (Platyhelminthes: Monogenea): tRNA gene arrangement reshuffling and implications for phylogeny. Parasit Vectors. 2017;10:462.
Altschul SF, Gish W, Miller W, Myers EW, Lipman DJ. Basic local alignment search tool. J Mol Biol. 1990;215:403–10.
Kearse M, Moir R, Wilson A, Stones-Havas S, Cheung M, Sturrock S, et al. Geneious Basic: an integrated and extendable desktop software platform for the organization and analysis of sequence data. Bioinformatics. 2012;28:1647–9.
Laslett D, Canback B. ARWEN: a program to detect tRNA genes in metazoan mitochondrial nucleotide sequences. Bioinformatics. 2008;24:172–5.
Bernt M, Donath A, Juhling F, Externbrink F, Florentz C, Fritzsch G, et al. MITOS: Improved de novo metazoan mitochondrial genome annotation. Mol Phylogenet Evol. 2013;69:313–9.
Zhang D. MitoTool software: Secondary MitoTool software. Wuhan: Institute of Hydrobiology, Chinese Academy of Sciences; https://github.com/dongzhang0725/MitoTool. Accessed 24 May 2018
Hadley W. ggplot2: Elegant graphics for data analysis. J Stat Softw. 2010;35:65–88.
Benson G. Tandem repeats finder: a program to analyze DNA sequences. Nucleic Acids Res. 1999;27:573–80.
Zuker M. Mfold web server for nucleic acid folding and hybridization prediction. Nucleic Acids Res. 2003;31:3406–15.
Zhang Z, Li J, Zhao XQ, Wang J, Wong GKS, Yu J. KaKs_Calculator: calculating Ka and Ks through model selection and model averaging. Genomics Proteomics Bioinformatics. 2006;4:259–63.
Librado P, Rozas J. DnaSP v5: a software for comprehensive analysis of DNA polymorphism data. Bioinformatics. 2009;25:1451–2.
Zhang D. BioSuite software: Secondary BioSuite software. https://github.com/dongzhang0725/BioSuite. Accessed 24 May 2018.
Katoh K, Standley DM. MAFFT multiple sequence alignment software version 7: improvements in performance and usability. Mol Biol Evol. 2013;30:772–80.
Talavera G, Castresana J. Improvement of phylogenies after removing divergent and ambiguously aligned blocks from protein sequence alignments. Syst Biol. 2007;56:564–77.
Keane TM, Creevey CJ, Pentony MM, Naughton TJ, McLnerney JO. Assessment of methods for amino acid matrix selection and their use on empirical data shows that ad hoc assumptions for choice of matrix are not justified. BMC Evol Biol. 2006;6:29.
Silvestro D, Michalak I. raxmlGUI: a graphical front-end for RAxML. Org Divers Evol. 2011;12:335–7.
Ronquist F, Teslenko M, van der Mark P, Ayres DL, Darling A, Höhna S, et al. MrBayes 3.2: efficient Bayesian phylogenetic inference and model choice across a large model space. Syst Biol. 2012;61:539–42.
Letunic I, Bork P. Interactive tree of life (iTOL) v3: an online tool for the display and annotation of phylogenetic and other trees. Nucleic Acids Res. 2016;44:242–5.
Le TH, Blair D, McManus DP. Mitochondrial genomes of parasitic flatworms. Trends Parasitol. 2002;18:206–13.
Li WX, Zhang D, Boyce K, Xi BW, Zou H, Wu SG, et al. The complete mitochondrial DNA of three monozoic tapeworms in the Caryophyllidea: a mitogenomic perspective on the phylogeny of eucestodes. Parasit Vectors. 2017;10:314.
Guo A. Moniezia benedeni and Moniezia expansa are distinct cestode species based on complete mitochondrial genomes. Acta Trop. 2016;166:287–92.
Waeschenbach A, Webster BL, Littlewood DTJ. Adding resolution to ordinal level relationships of tapeworms (Platyhelminthes: Cestoda) with large fragments of mtDNA. Mol Phylogenet Evol. 2012;63:834–47.
Feng Y, Feng HL, Fang YH, Su YB. Characterization of the complete mitochondrial genome of Khawia sinensis belongs among platyhelminths, cestodes. Exp Parasitol. 2017;177:35–9.
Zhang X, Duan JY, Shi YL, Jiang P, Zeng J, Wang ZQ, et al. Comparative mitochondrial genomics among Spirometra (Cestoda: Diphyllobothriidae) and the molecular phylogeny of related tapeworms. Mol Phylogenet Evol. 2017;117:75–82.
Guo A. Characterization of the complete mitochondrial genome of the cloacal tapeworm Cloacotaenia megalops (Cestoda: Hymenolepididae). Parasit Vectors. 2016;9:490.
Liu GH, Li C, Li JY, Zhou DH, Xiong RC, Lin RQ, et al. Characterization of the complete mitochondrial genome sequence of Spirometra erinaceieuropaei (Cestoda: Diphyllobothriidae) from China. Int J Biol Sci. 2012;8:640–9.
Zhang X, Duan JY, Wang ZQ, Jiang P, Liu RD, Cui J. Using the small subunit of nuclear ribosomal DNA to reveal the phylogenetic position of the plerocercoid larvae of Spirometra tapeworms. Exp Parasitol. 2017;175:1–7.
Hebert PDN, Cywinska A, Ball SL, deWaard JR. Biological identifications through DNA barcodes. Proc R Soc Lond B. 2003;270:313–21.
Arizono N, Shedko M, Yamada M, Uchikawa R, Tegoshi T, Takeda K, et al. Mitochondrial DNA divergence in populations of the tapeworm Diphyllobothrium nihonkaiense and its phylogenetic relationship with Diphyllobothrium klebanovskii. Parasitol Int. 2009;58:22–8.
Badaraco JL, Ayala FJ, Bart JM, Gottstein B, Haag KL. Using mitochondrial and nuclear markers to evaluate the degree of genetic cohesion among Echinococcus populations. Exp Parasitol. 2008;119:453–9.
Bazsalovicsova E, Kralova-Hromadova I, Stefka J, Scholz T, Hanzelova V, Vavrova S, et al. Population study of Atractolytocestus huronensis (Cestoda: Caryophyllidea), an invasive parasite of common carp introduced to Europe: mitochondrial cox1 haplotypes and intragenomic ribosomal ITS2 variants. Parasitol Res. 2011;109:125–31.
Hu D, Song X, Wang N, Zhong X, Wang J, Liu T, et al. Molecular identification of Echinococcus granulosus on the Tibetan Plateau using mitochondrial DNA markers. Genet Mol Res. 2015;14:13915–23.
Huyse T, Buchmann K, Littlewood DT. The mitochondrial genome of Gyrodactylus derjavinoides (Platyhelminthes: Monogenea) - a mitogenomic approach for Gyrodactylus species and strain identification. Gene. 2008;417:27–34.
The authors would like to thank the editor and the two anonymous reviewers for the time they have invested into reviewing our manuscript.
This work was funded by the National Natural Science Foundation of China (31572658), the Major Scientific and Technological Innovation Project of Hubei Province (2015ABA045) and the Earmarked Fund for China Agriculture Research System (CARS-45-15).
Availability of data and materials
The datasets supporting the conclusions of this article are available in the GenBank international nucleotide sequence repository under accession numbers MF671696 and MF671697.
Tapeworms were collected from fish in accordance with the recommended guidelines for animal experimentation by the Chinese Association for Laboratory Animal Sciences.
The authors declare that they have no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Table S1. Primers used to amplify and sequence the mitochondrial genome of Digramma interrupta and Ligula intestinalis. (DOCX 15 kb)
Table S2. The list of cestode species and outgroups used for comparative mitogenomic and phylogenetic analyses, and accession number, A+T content and skewness of different elements of each mitogenome. (XLSX 19 kb)
Table S3. General statistics (length and codons) for mitochondrial protein-coding genes and rRNAs of 38 cestodes. Abbreviations of species name are the initials of genus and species name combined. (XLSX 21 kb)
Figure S1. Relative Synonymous Codon Usage (RSCU) of Digramma interrupta and Ligula intestinalis. Codon families are labelled on the x-axis. Values on the top of the bars denote amino acid usage. (PDF 37 kb)