Introduction
Human echinostomiasis is a global zoonotic foodborne trematodiasis caused by flukes within the Echinostoma revolutum group, and despite its worldwide distribution, it is a particular public health problem in South East Asia (Chai, Reference Chai, Fried and Toledo2009; Toledo and Esteban, Reference Toledo and Esteban2016). Echinostoma revolutum (Fröhlich, 1802) Rudolphi, 1809, is a member of the family Echinostomatidae (Platyhelminthes: Echinostomata), and the ‘E. revolutum’ group is characterized by the ‘37-collar-spines’ found on the cercariae (Kostadinova, Reference Kostadinova, Gibson, Jones and Bray2005; Georgieva et al., Reference Georgieva, Faltýnková, Brown, Blasco-Costa, Soldánová, Sitko, Scholz and Kostadinova2014). There are nine Echinostoma species within the E. revolutum group including Echinostoma caproni, Echinostoma echinatum, Echinostoma friedi, Echinostoma jurini, Echinostoma miyagawai, Echinostoma paraensei, Echinostoma parvocirrus, E. revolutum and Echinostoma trivolvis; while in other Echinostomatidae species the number of collar spines may vary, such as 25–29 on Echinostoma hortense, 43 on Echinostoma malayanum, 41–45 on Hypoderaeum conoideum and 43–50 on Echinoparyphium recurvatum (Chai, Reference Chai, Fried and Toledo2009; Saijuntha et al., Reference Saijuntha, Sithithaworn, Duenngai, Kiatsopit, Andrews and Petney2011a). The similarity of these species within the E. revolutum complex usually required additional identification approaches for their discrimination, mostly enzymatic and molecular techniques (Saijuntha et al., Reference Saijuntha, Sithithaworn, Duenngai, Kiatsopit, Andrews and Petney2011a, Reference Saijuntha, Tantrawatpan, Sithithaworn, Andrews and Petney2011c; Georgieva et al., Reference Georgieva, Faltýnková, Brown, Blasco-Costa, Soldánová, Sitko, Scholz and Kostadinova2014; Tkach et al., Reference Tkach, Kudlai and Kostadinova2016).
The taxonomic status of E. revolutum is still controversial although recently a number of molecular studies have identified the parasite to be a highly cosmopolitan species comprising of several distinct geographical lineages corresponding to parasite populations with European, American, and Southeast Asian origins (Saijuntha et al., Reference Saijuntha, Sithithaworn, Duenngai, Kiatsopit, Andrews and Petney2011a; Georgieva et al., Reference Georgieva, Faltýnková, Brown, Blasco-Costa, Soldánová, Sitko, Scholz and Kostadinova2014; Faltýnková et al., Reference Faltýnková, Georgieva, Soldánová and Kostadinova2015; Nagataki et al., Reference Nagataki, Tantrawatpan, Agatsuma, Sugiura, Duenngai, Sithithaworn, Andrews, Petney and Saijuntha2015). The taxonomic identification and the phylogenetic assessment of each species within the ‘E. revolutum’ group and as well between member taxa in the family Echinostomatidae require accurate genomic data. Many attempts of interspecific clarification for the echinostomatids, particularly for those within the ‘37-collar-spined’ taxa have relied predominantly on tenuous morphological features (Georgieva et al., Reference Georgieva, Faltýnková, Brown, Blasco-Costa, Soldánová, Sitko, Scholz and Kostadinova2014; Faltýnková et al., Reference Faltýnková, Georgieva, Soldánová and Kostadinova2015; Nagataki et al., Reference Nagataki, Tantrawatpan, Agatsuma, Sugiura, Duenngai, Sithithaworn, Andrews, Petney and Saijuntha2015; Tkach et al., Reference Tkach, Kudlai and Kostadinova2016). However, by using single 28S ribosomal DNA, limited short mitochondrial DNA sequences (mtDNA) or a combination of both, new cryptic echinostome species and the systematic relationships within and between members within the Echinostomatidae have been revealed as well as their association with the other families in the superfamily Echinostomatoidea (Platyhelminthes: Echinostomata) (Olson et al., Reference Olson, Cribb, Tkach, Bray and Littlewood2003; Georgieva et al., Reference Georgieva, Selbach, Faltýnková, Soldánová, Sures, Skírnisson and Kostadinova2013, Reference Georgieva, Faltýnková, Brown, Blasco-Costa, Soldánová, Sitko, Scholz and Kostadinova2014; Nagataki et al., Reference Nagataki, Tantrawatpan, Agatsuma, Sugiura, Duenngai, Sithithaworn, Andrews, Petney and Saijuntha2015; Tkach et al., Reference Tkach, Kudlai and Kostadinova2016). However, in order to provide a detailed account of current species and to taxonomically validate echinostomes more effectively, it has been argued that genomic analyses could provide insights into the fine scale inter-relationships between echinostome species (Detwiler et al., Reference Detwiler, Bos and Minchella2010; Faltýnková et al., Reference Faltýnková, Georgieva, Soldánová and Kostadinova2015; Gordy and Hanington, Reference Gordy and Hanington2019). In fact, the analyses of complete mitochondrial genomes to perform taxonomic and phylogenetic analyses of other members of the Echinostomata, as well as other trematode species, have been widely used and have provided not only a deeper understanding of the evolutionary relationships within and between trematode families but have also provided essential molecular markers for population genetics and diagnostics, crucial for modern epidemiological studies (Wey-Fabrizius et al., Reference Wey-Fabrizius, Podsiadlowski, Herlyn and Hankeln2013; Georgieva et al., Reference Georgieva, Faltýnková, Brown, Blasco-Costa, Soldánová, Sitko, Scholz and Kostadinova2014; Faltýnková et al., Reference Faltýnková, Georgieva, Soldánová and Kostadinova2015).
However, many morphologically similar species, and particularly, for those of the ‘collar-spined’ Echinostoma spp. in the Echinostomatoidea lack complete mitochondrial genomic data. Currently, only four of the nine species of the ‘E. revolutum’ group, including E. caproni, E. paraensei, E. miyagawai, E. hortense (Saijuntha et al., Reference Saijuntha, Tantrawatpan, Sithithaworn, Andrews and Petney2011c), and a few species within the Echinostomata suborder have complete mitochondrial genomes available (Yang et al., Reference Yang, Gasser, Koehler, Wang, Zhu, Chen, Feng, Hu and Fang2015; Fu et al., Reference Fu, Jin, Li and Liu2019; Suleman et al., Reference Suleman, Heneberg, Zhou, Muhammad, Zhu and Ma2019; Li et al., Reference Li, Ma, Lv, Hu, Qiu, Chang and Wang2019a, Reference Li, Qiu, Zeng, Diao, Chang, Gao, Zhang and Wang2019b).
This current study determined the complete mitochondrial genome sequence of E. revolutum and correlatively characterized its genomic features and compared them with those previously reported in the superfamily Echinostomatoidea. A phylogeny for members of families in the suborders Echinostomata, Opisthorchiata, Troglotremata, Pronocephalata, and Xiphidiata is provided.
Materials and methods
Samples, DNA extraction and species identification
Adult E. revolutum flukes were obtained from the intestines of the naturally infected domestic ducks from abattoirs in Khon Kaen province, Thailand. The flukes were thoroughly washed in physiological saline and morphologically identified based on the size of the body and circumoral disc, the appearance of testes and the presence of ‘37-collar spines’ around the head (Miliotis and Bier, Reference Miliotis and Bier2003; Georgieva et al., Reference Georgieva, Faltýnková, Brown, Blasco-Costa, Soldánová, Sitko, Scholz and Kostadinova2014). The worms were individually fixed in 70% (v/v) ethanol and stored at −20°C until use. Subsequently, species were confirmed by molecular phylogenetic analyses using nuclear ITS-1, mitochondrial cox1 and nad1 markers (Saijuntha et al., Reference Saijuntha, Sithithaworn, Duenngai, Kiatsopit, Andrews and Petney2011a, Reference Saijuntha, Tantrawatpan, Sithithaworn, Andrews and Petney2011b; Nagataki et al., Reference Nagataki, Tantrawatpan, Agatsuma, Sugiura, Duenngai, Sithithaworn, Andrews, Petney and Saijuntha2015).
Total genomic DNA was extracted from individual worms using the DNA extraction kit (QIAGEN, Hilden, Germany) following the manufacturer's protocol. The E. revolutum-species used for mitochondrial sequencing in this study belonged to the nad1-based E. revolutum-Eurasian lineage (Nagataki et al., Reference Nagataki, Tantrawatpan, Agatsuma, Sugiura, Duenngai, Sithithaworn, Andrews, Petney and Saijuntha2015).
PCR strategies for obtaining the complete mitochondrial genome
The first, initial specific primer pairs (ERE1F/ERE2R; ERE3F/ERE4R; ERE5F/ERE6R) designed based on the conserved nucleotide sequences aligned by those E. revolutum-mt sequences, cox1, nad1, rrnS (12S), respectively, available in GenBank and others, namely platyhelminth-universal primers (TRECOBF; TRECOBR; GLYF; GLYR) previously described in Le et al. (Reference Le, Nguyen, Nguyen, Doan, Agatsuma and Blair2019) were used. They were paired to bind on the target regions for amplification of long PCR of 4.0–7.5 kb or short of <4.0 kb overlapping fragments. The sequence data obtained were used to design further E. revolutum-specific primers (Table 1).
Primers for amplification and sequencing fragments of the mitochondrial genome of Echinostoma revolutum
NCR, non-coding region.
All reagents and kits used in this study were from Thermo Fisher Scientific Inc. (Waltham, MA, USA), including Phusion for long, and Dream Taq PCR Master Kits for short PCRs. PCRs were prepared in 50 μL volume with the addition of DMSO to 1.5%, and performed in an MJ PTC-100 Thermal Cycler. Long PCRs were conducted with initial denaturation at 98°C for 30 s, followed by 35 cycles, each consisting of denaturation step for 30 s at 98°C, annealing/extension step at 72°C for 6–8 min and final extension at 72°C for 10 min (in some cases, at 68°C). Short PCRs were started at 95°C for 5 min, followed by 35 cycles consisting of denaturation for 1 min at 94°C, annealing at 52°C for 1 min, extension at 72°C for 2–5 min and a final extension at 72°C for 7 or 10 min. A negative (no-DNA) control was included in some cases. The PCR products (5 − 10 μL of each) were examined on a 1% agarose gel, stained with ethidium bromide and visualized under UV light (Wealtec, Sparks, NV, USA). The primer-walking sequencing was applied until the complete sequence for the whole fragment, and the overlapping assembly was used to complete the mitochondrial genome.
Characterization of mitogenomic features
Protein-encoding genes (PCGs) were identified by alignment with the available mt genomes of other Echinostoma trematode species and ATG/GTG as start and TAA/TAG as stop codons were used to define gene boundaries. PCGs were translated using the echinoderm/flatworm mitochondrial genetic code: translation Table 9 in GenBank. Nucleotide and codon composition were analysed with MEGA 7.0 (Kumar et al., Reference Kumar, Stecher and Tamura2016) and codon usage for all PCGs was determined with the program GENE INFINITY (Codon Usage: http://www.geneinfinity.org/sms/sms_codonusage.html). Nucleotide percentage (%) for comparison of individual/concatenated PCGs and mitochondrial ribosomal genes (MRGs) between E. revolutum and 14 representative members of the superfamily Echinostomatoidea (Table 2) was determined by using GENEDOC 2.7 for alignment, Gblocks 0.91b (Castresana, Reference Castresana2000) (online accession at http://molevol.cmima.csic.es/castresana/Gblocks_server.html) for picking the best quality block (10,112 bp) and MEGA 7.0 for percentage estimation.
Locations of genes and other features in the complete mitochondrial genome of Echinostoma revolutum (17,030 bp) (GenBank: MN496162)
bp, base pair; aa, amino acid; start, start codon; stop, stop codon; Int. seq., intergenic sequence (+, number of nucleotides before the start of the following gene; −, number of nucleotides overlapping with the following gene); LRU, long repeat unit; SRU, short repeat unit; IntS, internal spacer sequence between LRU7 and SRU4; unique seq, nucleotide sequence between SRU1 and cox3.
a tRNAs lacking DHU-arm.
The transfer RNA genes (tRNA or trn) were identified using tRNAscan-SE 1.21 program (www.genetics.wustl.edu/eddy/tRNAscan-SE/) (Lowe and Eddy, Reference Lowe and Eddy1997); ARWEN at http://mbio-serv2.mbioekol.lu.se/ARWEN/ (Laslett and Canback, Reference Laslett and Canback2008) for finding the best final tRNA sequences and secondary structures. Any tRNAs not detected by these programs were found by inspection of the sequences, based on the alignment with sequences of other trematode and by their potential formation of tRNA configuration. The ribosomal 16S (rrnL) and 12S (rrnS) RNA genes were recognized as described in Le et al. (Reference Le, Nguyen, Nguyen, Doan, Agatsuma and Blair2019) in the region located between tRNAThr (trnT) and cox2 separated by tRNACys (trnC), respectively.
The nucleotide composition, AT and GC content for concatenated 12 PCGs (not excluding the overlapped sequences between nad4L and nad4), two MRGs and complete mt genome for 15 members of the Echinostomatoidea were determined by MEGA7.0, and the AT and GC skewness values (from −1 to +1) calculated according to the formula by Perna and Kocher (Reference Perna and Kocher1995) [AT skew = (A + T)/(A–T) and GC skew = (G + C)/(G–C)].
The non-coding region (NCR) was determined by the recognition of boundaries between tRNAGlu (trnE) and cox3. Tandem Repeat Finder v3.01 (Benson, Reference Benson1999) was used to detect repeat units (RUs) in the NCR of mitogenome of E. revolutum in this study and other Echinostoma spp. and digeneans which were not available in GenBank or not previously reported.
Phylogenetic analyses
Concatenated amino acid sequences of the 12 PCGs of E. revolutum and 44 species from 13 families [i.e., Echinostomatidae, Fasciolidae, Himasthlidae, Echinochasmidae, Cyclocoelidae, Paramphistomidae, Gastrothylacidae, Notocotylidae, Troglotrematidae/(Paragonimidae), Heterophyidae, Opisthorchiidae, Diclocoeliidae, and Schistosomatidae] in the superfamilies of Echinostomatoidea, Paramphistomoidea, Pronocephaloidea, Troglotrematoidea, Opisthorchioidea, and Gorgoderoidea were aligned for phylogenetic analysis. The sequence of Schistosoma haematobium (Digenea: Schistosomatidae) was chosen as an outgroup (Littlewood et al., Reference Littlewood, Lockyer, Webster, Johnston and Le2006). The alignment was constructed by GENEDOC2.7, confirmed by MAFFT 7.122 (Katoh and Standley, Reference Katoh and Standley2013) and finalized by Gblocks 0.91b. The final alignment block of 2,993–3,025 amino acids without poorly aligned regions was picked out for phylogenetic analysis. Tree was constructed using maximum likelihood by MEGA 7.0 with a bootstrap of 1000 replications. The substitution model with the best score according to the Bayesian information criterion was the Jones, Taylor and Thornton +F+G+I model, with residue frequencies estimated from the data (+F), rate variation along the length of the alignment (+G) and allowing for a proportion of invariant sites (+I).
Results
Gene organization and genomic features
The complete mitochondrial genome of E. revolutum was shown to be 17,030 bp in size (GenBank accession no. MN496162) (Fig. 1). As common in other trematodes, the E. revolutum mitogenome has one-direction transcription, similar gene organization and content with the exception of African Schistosoma spp. It comprises of 12 protein coding genes (atp6, cox1-3, cytb, nad1-6, nad4L), two ribosomal RNA (rrnL and rrnS) and 22 transfer RNA genes (tRNA or trn) similar to those of common digeneans (Table 2).
A schematic drawing of a circular map of the mitochondrial genome of Echinostoma revolutum (GenBank: MN496162). Protein-coding and ribosomal large and small subunit genes are abbreviated according to our previous publications (Le et al., Reference Le, Nguyen, Nguyen, Doan, Dung and Blair2016, Reference Le, Nguyen, Nguyen, Doan, Agatsuma and Blair2019). The transfer RNA genes (tRNA) are marked with three letter-amino acid abbreviations (see: Table 2). The non-coding region (NCR) located between tRNAGlu and cox3 consists of seven long (LRU1–7), four short repeat units (SRU1–4), and internal spacer sequence (IntS) between LRU7 and SRU4.
Echinostoma revolutum has typical mt structural features of the platyhelminths and does not contain atp8 and has the overlapped region between nad4L and nad4 genes by 40 bp (Table 2). Five protein-coding genes used GTG (nad4L, nad2, nad1, cox1, nad5) and other seven used ATG as start codons; and seven genes used TAG and five used TAA for termination. Boundaries between cytb and nad4L, between tRNAAsp and nad1, from tRNAThr to rrnS (12S), covering rrnL (16S), tRNACys genes, and between repeats in the NCR are continuous whilst there are large intergenic spacers of 33 or 30 bp between other genes (cox1 and tRNATh; and tRNAVal and tRNAAla), respectively.
The mt genome of E. revolutum encodes 22 transfer RNAs, ranging from 60 (tRNAS1(AGN)) to 71 nucleotides (tRNAHis). Twenty have common ‘cloverleaf’ folding into secondary structures with the complete four arms but two for Serine, tRNAS1(AGN) and tRNAS2(UCN), possess special forms missing DHU-arms (Table 2; SFig. 1). Two ribosomal RNA genes, rrnL (977 bp) and rrnS (756 bp long), are located between the tRNAThr and cox2, separated by tRNACys. The order of the mitochondrial DNA block of [cox1-tRNAThr-rrnL-tRNACys-rrnS-cox2-nad6] is highly conserved in all the trematodes, including E. miyagawai, Ech. japonicus, Fas. magna, F. hepatica, F. gigantica, and Asian Schistosoma spp. (Le et al., Reference Le, Blair and McManus2001, Reference Le, Blair and McManus2002, Reference Le, Nguyen, Nguyen, Doan, Dung and Blair2016; Liu et al., Reference Liu, Gasser, Young, Song, Ai and Zhu2014; Ma et al., Reference Ma, He, Liu, Leontovyč, Kašný and Zhu2016; Fu et al., Reference Fu, Jin, Li and Liu2019; Li et al., Reference Li, Qiu, Zeng, Diao, Chang, Gao, Zhang and Wang2019b).
Base composition and comparative analyses
The base composition was A (18.81%), T (47.40%), G (23.50) and C (10.29% in the mt genome of E. revolutum and the A + T content was 62.21% for PCGs and their skewness values were −0.46 for A + T and 0.391 for G + C, respectively). MRGs showed a similar percentage of overall A + T (62.73%) and G + C (37.27%) but their skewness values were considerably different (−0.179/A + T; and 0.275/G + C) due to the bias use of A over T in PCGs than in MRGs (Table 3).
Base composition and skewness value for the mitochondrial protein-coding (PCGs) and mitoribosomal genes (MRGs) of 15 members of the superfamily Echinostomatoidea
The divergence rate (%) inferred from the nucleotide pairwise comparison of 12 individual mitochondrial protein-coding and two ribosomal genes between E. revolutum and 14 members of Echinostomatoidea indicated that the rate was the lowest level of divergence between E. revolutum and E. miyagawai (8.99%/nad4L–18.4%/nad4; 6.63%/rrnS–8.93%/rrnL), and in average, 14.89%/PCGs for protein-coding genes and 8.29%/MRGs for ribosomal genes, respectively.
The highest nucleotide sequence divergence between E. revolutum and Echinostomatoidea trematodes was 39.5% in comparison with Tracheophilus cymbius (Cyclocoelidae) and 38.16% for Ech. japonicus (Echinochasmidae) for PCGs (Table 4). Overall, the nucleotide sequence of E. revolutum in each gene differed from 6.63%/rrnS/(E. miyagawai) to 59.89%/nad5/(T. cymbius). Within the Echinostomatidae, the interspecific variation does not exceed 37%, as seen between atp6 genes of E. revolutum and H. conoideum.
Nucleotide comparison for divergence rate (%) of individual and concatenated protein-coding (PCGs) and mitoribosomal genes (MRGs) between Echinostoma revolutum and 14 representative members of the superfamily Echinostomatoidea (Platyhelminthes: Echinostomata)
Ecap: Echinostoma caproni; Emiy: E. miyagawai; Epar: E. paraensei; EJM: Echinostoma sp. JM-2019; Asuf: Artyfechinostomum sufrartyfex; Hcon: Hypoderaeum conoideum; Ejap: Echinochasmus japonicus; Fhep: Fasciola hepatica; Fgig: F. gigantica; Fjac: Fascioloides jacksoni; Fmag: Fas. magna; Fbus: Fasciolopsis buski; AWAK: Acanthoparyphium sp. WAK-2018; Tcym: Tracheophilus cymbius. Rows of data subjected to discussion in the text were background colour shaded, including the lowest divergence between E. revolutum and E. miyagawa (family: Echinostomatidae) and the highest rate between E. revolutum and Ech. japonicus (Echinochasmidae), and between E. revolutum and T. cymbius (Cyclocoelidae).
The codon usage in mtDNAs of all the Echinostomatidae trematodes (E. revolutum; E. caproni; E. miyagawai; E. paraensei; Echinostoma sp. JM-2019; A. sufrartyfex; H. conoideum) is biased to the use of TTT (for Phenylanine), TTG (for Leucine) and GTT (for Valine). Multiple Thymine (T) in use in these codons facilitates the mostly used frequency (from 5.96% GTT/Val in H. conoideum to 10.65% TTT/Phe in E. caproni). The least frequently used codons, comprising mostly G and C, are CGC (for Arginine), AAC (for Asparagine) and GAC (for Aspartic acid), ranging from one to two (0.03–0.06%) to six to seven (0.18–0.21%) were noted. Clear bias was seen to the use of TAG (7–12 codons) for termination of 12 PCGs rather than TAA (0–5) in mt PCG genes of all the eight echinostomids (STable 1).
Polymorphism featured by repeat units in non-coding regions of Echinostoma spp.
The NCR of E. revolutum was identified by recognition of boundary of tRNAGlu (trnE) and cox3 gene, which is of 3,549 bp in length, perhaps the longest in the mt genomes of the echinostomid flatworm sever fully sequenced to date (Tables 2 and 5). The NCR of this species possesses seven long, identical RUs (LRU1 to LRU7, 317 bp each) and four short, identical RUs (SRU1 to SRU4, 207 bp each) tandemly arranged after each other (Tables 2 and 5; Fig. 1; GenBank: MN496162). Between LRU7 and SRU4, there is a linking region of an internal spacer sequence of 377 bp which contained 188 bp, partial of LRU (designated as IntS-half1) and 189 bp, partial of SRU (IntS-half2). A unique sequence region of 130 nucleotides continuously occurs between SRU1 and cox3 (Table 2; Fig. 1).
Number and type of the repetitive sequences in non-coding regions (NCR) of 15 representative members of the superfamily Echinostomatoidea indicating high polymorphism and interspecific/intergeneric variation
a Non-coding region in E. paraensei not fully sequenced.
Tandem RUs were also found in E. miyagawai (two RUs, 319 bp each), in E. paraensei (at least three RUs, 206 bp each in the partially sequenced NCR), in Echinostoma sp. JM-2019 (five LRUs, 245 bp each and two SRUs, 166 bp each), in A. sufrartyfex (two RUs, 144 bp each) which is variable in numbers and length (Table 5). The size of the mt genome differed among echinostomes and digeneans; this is due to the variable length of their NCRs rich in multiple RUs (Table 3).
Phylogenetic analysis
The topology of the phylogenetic tree of taxonomic relationship indicated clear positions of five suborders, including Echinostomata, Pronocephalata, Troglotremata, Opisthorchiata, and Xiphidiata where E. revolutum, grouped in a monophyletic subclade as a sister taxa to E. miyagawai and paraphyletic to the other echinostomatids in the Echinostomatidae (Fig. 2). Monophyly of Echinostomatidae and Fasciolidae clearly resolved with respect to Echinochasmidae, Himasthlidae, and Cyclocoelidae; these were rendered paraphyletic in the suborder Echinostomata (Fig. 2). The high nodal bootstrap values well supported clear taxonomic relationships of the ‘E. revolutum’ group in the Echinostomatoidea and this seemed to be in the paraphyletic position with all the other superfamilies, Pronocephalata, Troglotremata, Opisthorchiata, and Xiphidiata in the digenean order Plagiorchiida.
A maximum likelihood phylogenetic tree showing the position of Echinostoma revolutum (diamond symbol) based on the analysis of concatenated amino acid sequence data for the 12 mitochondrial proteins of 45 digenean species/strains. Thirteen families [Echinostomatidae, Fasciolidae, Himasthlidae, Echinochasmidae, Cyclocoelidae, Paramphistomidae, Gastrothylacidae, Notocotylidae, Troglotrematidae/(Paragonimidae), Heterophyidae, Opisthorchiidae, Diclocoeliidae, and Schistosomatidae] belonging to six superfamilies indicated by arrows, Echinostomatoidea (ECH), Paramphistomoidea (PAR), Pronocephaloidea (PRO), Troglotrematoidea (TRO), Opisthorchioidea (OPI), and Gorgoderoidea (GOR), are represented. Schistosoma haematobium (Platyhelminthes: Schistosomatidae) is included as an outgroup. Nodal support values evaluated using 1000 bootstrap resamplings are shown on each branch. The scale bar represents the number of substitutions per site. Accession numbers are given for each species/strains and country name (in bracket) of their origin (where available) at the end of each sequence.
Discussion
The complete mitochondrial genome of E. revolutum (Fröhlich, 1802) Rudolphi, 1809, was 17,030 bp in size; the longest of all the Echinostomatoidea to date sequenced, although the mitogenome of E. paraensei (KT008005) was claimed longer, 20,298 bp, but some of 5,600 nucleotides were only of estimation (Liu et al., Reference Liu, Gasser, Young, Song, Ai and Zhu2014; Yang et al., Reference Yang, Gasser, Koehler, Wang, Zhu, Chen, Feng, Hu and Fang2015; Le et al., Reference Le, Nguyen, Nguyen, Doan, Dung and Blair2016; Ma et al., Reference Ma, He, Liu, Leontovyč, Kašný and Zhu2016, Reference Ma, Sun, He, Liu, Ai, Chen and Zhu2017; Fu et al., Reference Fu, Jin, Li and Liu2019; Suleman et al., Reference Suleman, Heneberg, Zhou, Muhammad, Zhu and Ma2019; Li et al., Reference Li, Ma, Lv, Hu, Qiu, Chang and Wang2019a, Reference Li, Qiu, Zeng, Diao, Chang, Gao, Zhang and Wang2019b).
The length of the mt genome of E. revolutum seemed to be one of the longest among trematodes fully obtained to date, shorter than the estimated, partially sequenced congener E. paraensei, but was slightly longer than other echinostomids, including H. conoideum (Yang et al., Reference Yang, Gasser, Koehler, Wang, Zhu, Chen, Feng, Hu and Fang2015), E. miyagawai (Fu et al., Reference Fu, Jin, Li and Liu2019; Li et al., Reference Li, Qiu, Zeng, Diao, Chang, Gao, Zhang and Wang2019b), Ech. japonicus (Le et al., Reference Le, Nguyen, Nguyen, Doan, Dung and Blair2016) and two cyclocoelids, Uvitellina sp. and T. cymbius (Suleman et al., Reference Suleman, Heneberg, Zhou, Muhammad, Zhu and Ma2019; Li et al., Reference Li, Ma, Lv, Hu, Qiu, Chang and Wang2019a). It was considerably longer than many fasciolids, such as Fa. buski (GenBank: KX169163) (Ma et al., Reference Ma, Sun, He, Liu, Ai, Chen and Zhu2017), F. gigantica (KF543342), F. hepatica (AF216697), Fasciola/Fascioloides jacksoni (KX787886) and Fas. magna (KU060148) (Liu et al., Reference Liu, Gasser, Young, Song, Ai and Zhu2014; Ma et al., Reference Ma, He, Liu, Leontovyč, Kašný and Zhu2016).
The tRNAs which were lacking DHU-arm for Serine in E. revolutum are usually found in many digenean mitogenomes, i.e. Echinococcus granulosus, F. hepatica (Le et al., Reference Le, Blair and McManus2001, Reference Le, Blair and McManus2002), Fas. magna (Ma et al., Reference Ma, He, Liu, Leontovyč, Kašný and Zhu2016), E. miyagawai (Li et al., Reference Li, Qiu, Zeng, Diao, Chang, Gao, Zhang and Wang2019b) and Fas. jacksoni (KX787886).
The gene organization, comparative description of genomic features with other members of Echinostomatidae, particularly, with E. miyagawai isolates (from Hunan and Helongjang of China) and those of the digenean Echinostomata were presented. In mtDNA sequence of E. revolutum, the nucleotide usage clearly biased to AT, and thus, constituting their negative skewness. Skewness values for A + T are consistent with those of E. miyagawai, E. paraensei (Echinostomatidae) and Acanthoparyphium sp. WAK-2018 (Himasthlidae), considerably higher than all of the members of Fasciolidae, slightly higher than other echinostomatids (E. caproni, Echinostoma sp. JM-2019, A. sufrartyfex and H. conoideum), echinochasmid (Ech. japonicus) but lower than the cyclocoelid T. cymbius. The G + C content and skewness of E. revolutum seemed to be of the lowest (GC skew = 0.391) among all species studied here (Table 3).
Echinostoma revolutum and E. miyagawai shared more common genomic features than others in the genus Echinostoma and family Echinostomatidae. The pattern of the usage of ATG/GTG start and TAG/TAA stop codons, the AT composition bias, the negative AT-skewness, and the most for Phe/Leu/Val and the least for Arg/Asn/Asp codons in E. revolutum were usual and similar to the members of Echinostoma and digenean trematodes.
The presence of 11 tandem repeats in NCR (GenBank: MN496162) made the NCR of E. revolutum longer and more complex relative to other echinostomatids. The repetitive sequence richness in NCR was a typical genomic feature commonly seen in a number of species, specifying high-level polymorphism in Echinostomata and other digeneans (Table 5) (Le et al., Reference Le, Blair and McManus2001, Reference Le, Nguyen, Nguyen, Doan, Dung and Blair2016, Reference Le, Nguyen, Nguyen, Doan, Agatsuma and Blair2019; Liu et al., Reference Liu, Gasser, Young, Song, Ai and Zhu2014; Ma et al., Reference Ma, He, Liu, Leontovyč, Kašný and Zhu2016, Reference Ma, Sun, He, Liu, Ai, Chen and Zhu2017; Fu et al., Reference Fu, Jin, Li and Liu2019). For some of Echinostoma spp. which had their complete mitogenomes fully sequenced to date, the number of RUs was fewer or absent, and the length of the NCR was less than those of E. revolutum.
The actual size of the mitogenomes of other echinostomids may have been an underestimation in some of the original individuals sampled as several of the repeat elements may not have been considered or incorporated in the initial analyses (Yang et al., Reference Yang, Gasser, Koehler, Wang, Zhu, Chen, Feng, Hu and Fang2015; Fu et al., Reference Fu, Jin, Li and Liu2019; and GenBank: MH212284; KY548763; AP017706) (Table 5), as a result of missing a part of the region containing more RUs. The missing part of the NCR may be the result of an inaccurate PCR experiment that was carried out without verification (Kinkar et al., Reference Kinkar, Korhonen, Cai, Gauci, Lightowlers, Saarma, Jenkins, Li, Li, Young and Gasser2019; Oey et al., Reference Oey, Zakrzewski, Narain, Devi, Agatsuma, Nawaratna, Gobert, Jones, Ragan, McManus and Krause2019). In E. revolutum, the NCR was successfully amplified and accurately sequenced from a number of the verified PCR products and the RUs were confirmed to occur in the expanded NCR giving its complete mtDNA sequence as the second longest among members of the Echinostomatidae. Such repetitive regions have also occurred in the mtDNA of, for example, E. granulosus G1 with the addition of a 4.4 kb tandem repeat region consisting ten RUs (Kinkar et al., Reference Kinkar, Korhonen, Cai, Gauci, Lightowlers, Saarma, Jenkins, Li, Li, Young and Gasser2019), or Paragonimus westermani from the Arunachal Pradesh State (India), with the full mtDNA of 20.3 kb comprising of a long repetitive region in the isolate of the East Siang district (Oey et al., Reference Oey, Zakrzewski, Narain, Devi, Agatsuma, Nawaratna, Gobert, Jones, Ragan, McManus and Krause2019) instead of 14,965 bp in the isolate of the Changlang District (Biswal et al., Reference Biswal, Chatterjee, Bhattacharya and Tandon2014). However, it should be noted that the length and number of repeats are genetically variable between geographical isolates of a trematode species, as seen in P. ohirai and P. westermani (Le et al., Reference Le, Nguyen, Nguyen, Doan, Agatsuma and Blair2019; Oey et al., Reference Oey, Zakrzewski, Narain, Devi, Agatsuma, Nawaratna, Gobert, Jones, Ragan, McManus and Krause2019) and there are no quantity of repeats in individuals to be considered fixed. In many other taxa of trematodes reported to date, for example, Ech. japonicus, F. hepatica, Fa. buski, Fas. magna, P. ohirai, repetitive units either of long or short sequences and even various quantity within a species, frequently occurred and commonly found (Le et al., Reference Le, Blair and McManus2001, Reference Le, Nguyen, Nguyen, Doan, Dung and Blair2016, Reference Le, Nguyen, Nguyen, Doan, Agatsuma and Blair2019; Liu et al., Reference Liu, Gasser, Young, Song, Ai and Zhu2014; Ma et al., Reference Ma, He, Liu, Leontovyč, Kašný and Zhu2016, Reference Ma, Sun, He, Liu, Ai, Chen and Zhu2017; Fu et al., Reference Fu, Jin, Li and Liu2019; Li et al., Reference Li, Qiu, Zeng, Diao, Chang, Gao, Zhang and Wang2019b). Interestingly, none of the repetitive units was found in E. caproni, H. conoideum and T. cymbius (Yang et al., Reference Yang, Gasser, Koehler, Wang, Zhu, Chen, Feng, Hu and Fang2015; Li et al., Reference Li, Ma, Lv, Hu, Qiu, Chang and Wang2019a). The occurrence of repetitive sequences in tandem order in many species certainly is one of the most interesting genomic features, specifying the high-level polymorphism in the NCRs of digenean trematodes. Also, although tandem repeats are common in eukaryotic mitochondrial genomes, their functional role is still not completely understood. However, they do appear to have an accelerated rate of evolution and some involvements in the regulation of the mtDNA coding region (Lunt et al., Reference Lunt, Whipple and Hyman1998; Gemayel et al., Reference Gemayel, Vinces, Legendre and Verstrepen2010).
The phylogenetic tree presented in this study indicated the precise placement of E. revolutum in the Echinostomatidae matched closely the relationships described in previous studies using nuclear ribosomal sequences (Olson et al., Reference Olson, Cribb, Tkach, Bray and Littlewood2003; Tkach et al., Reference Tkach, Kudlai and Kostadinova2016). The Fasciolidae and Echinostomatidae are always sister groups within the Echinostomata, as are the Heterophyidae and Opisthorchiidae (within the Opisthorchiata). The echinostomatid species in the tree were also clustered well in the phylogenetic studies by Li et al. (Reference Li, Qiu, Zeng, Diao, Chang, Gao, Zhang and Wang2019b) and Fu et al. (Reference Fu, Jin, Li and Liu2019) using the complete mitochondrial genome sequences with one exception. The exception was, in their studies, the closeness of E. myiagawai and E. paraensei (sister groups); however, in this current study, E. myiagawai is closely associated with E. revolutum. This discrepancy of echinostomatid relationships was explained by the lack of mtDNA of E. revolutum for comparative analysis at that time. In our present study, the echinostomid relationship was also resolved, that the ‘37-collar-spine E. revolutum’ group members, E. revolutum, E. myiagawai, E. caproni and E. paraensei, were clustered together indicating their genetically close relationships, rather than other Echinostoma species, A. sufrartyfex, Echinostoma sp. JM-2019 and H. conoideum (Fig. 2). This relatedness of the ‘revolutum’ group is reflected by the very low divergence rate (%) of individual and concatenated PCGs and MRGs between E. revolutum and E. myiagawai/E. caproni, which varied within the least, 6.63% and the highest, 20.03%, compared to the rate of more than 20% in all cases of other echinostomid species (Table 4).
Conclusion
The fully annotated mitogenome of E. revolutum and comparative description of mitogenomic features of echinostomids in the present study provide well-supported resolution of relationships of the ‘revolutum’ group and the Echinostomata in the relation of other suborders in Plagiorchiida (Platyhelminthes: Digenea). The characterization revealed the taxonomic and phylogenetic relationships of E. revolutum to the echinostomatid species and other members in Echinostomatoidea. Molecular analyses of recently available mitogenomic sequences from Echinostomatidae, Himasthlidae and Cyclocoeliidae and comparisons of genetic features have emphasized the ‘revolutum’ group to be complex, but phylogenetic analysis has confirmed monophyly Echinostomatidae and Fasciolidae. Data from this species and additional Echinostoma spp. will be useful for clarification and reappraisal of the complex echinostome group and for the use in the field of molecular taxonomic, diagnostic, systematic, epidemiological, phylogenetic and population studies of trematodes.
Supplementary material
The supplementary material for this article can be found at https://doi.org/10.1017/S0031182020000128
Acknowledgements
We express our thanks to colleagues and technicians for contribution to our laboratory work. We would also like to express our thanks to the external reviewers for their extremely informative and constructive comments.
Author contributions
Thanh Hoa Le conceived the study, analysed the final data, prepared figures and tables and wrote the manuscript; Linh Thi Khanh Pham, Huong Thi Thanh Doan and Xuyen Thi Kim Le conducted laboratory works and sequence analyses. Weerachai Saijuntha collected, molecularly identified and provided specimens. R.P.V. Jayanthe Rajapakse and Weerachai Saijuntha reviewed the drafts. Scott P. Lawton and Thanh Hoa Le completed and approved the manuscript. All authors read and approved the final manuscript.
Financial support
This research is funded by the Vietnam National Foundation for Science and Technology Development (NAFOSTED) under grant number 108.02-2017.09 (PI: ThanhHoa Le).
Conflict of interest
None.
Ethical standards
None.
