The complete mitochondrial genomes of two freshwater snails provide new protein-coding gene rearrangement models and phylogenetic implications
Mu et al. Parasites & Vectors
The complete mitochondrial genomes of two freshwater snails provide new protein- coding gene rearrangement models and phylogenetic implications
Xidong Mu 0
Yexin Yang 0
Yi Liu 0
Du Luo 0
Meng Xu 0
Hui Wei 0
Dangen Gu 0
Hongmei Song 0
Yinchang Hu 0
0 Key Laboratory of Tropical&Subtropical Fishery Resource Application & Cultivation, Ministry of Agriculture, Pearl River Fisheries Research Institute, Chinese Academy of Fishery Sciences , Xingyu Road1, Guangzhou 510380 , China
Background: Mitochondrial (mt) genome sequences are widely used for species identification and to study the phylogenetic relationships among Gastropoda. However, to date, limited data are available as taxon sampling is narrow. In this study we sequenced the complete mt genomes of the freshwater gastropods Radix swinhoei (Lymnaeidae) and Planorbarius corneus (Planorbidae). Based on these sequences, we investigated the gene rearrangement in these two species and the relationships with respect to the ancestral gene order and assessed their phylogenetic relationships. Methods: The complete mt genomes of R. swinhoei and P. corneus were sequenced using Illumina-based paired-end sequencing and annotated by comparing the sequence information with that of related gastropod species. Putative models of mitochondrial gene rearrangements were predicted for both R. swinhoei and P. corneus, using Reishia clavigera mtDNA structure as the ancestral gene order. The phylogenetic relationships were inferred using thirteen protein sequences based on Maximum likelihood and Bayesian inference analyses. Results: The complete circular mt genome sequences of R. swinhoei and P. corneus were 14,241 bp and 13,687 bp in length, respectively. Comparison of the gene order demonstrated complex rearrangement events in Gastropoda, both for tRNA genes and protein-coding genes. The phylogenetic analyses showed that the family Lymnaeidae was more closely related to the family Planorbidae, consistent with previous classification. Nevertheless, due to the position recovered for R. swinhoei, the family Lymnaeidae was not monophyletic. Conclusion: This study provides the complete mt genomes of two freshwater snails, which will aid the development of useful molecular markers for epidemiological, ecological and phylogenetic studies. Additionally, the predicted models for mt gene rearrangement might provide novel insights into mt genome evolution in gastropods.
Pulmonate; Mitochondrial genome; Gene order; Phylogeny
The hyperdiverse pulmonate gastropods  contains the
medically important clade Hygrophila, which comprises
the freshwater families Acroloxidae, Chilinidae,
Planorbidae, Lymnaeidae and Physidae . Many of these
freshwater snails are intermediate hosts for flatworm
parasites and transmit infectious diseases of human and
veterinary importance such as fascioliasis, cercarial
dermatitis and schistosomiasis [3–5]. Accurate
identification of species and analysis of genetic variation within
populations is essential for studying molecular
epidemiology and controlling parasite infection. However,
previous studies suggest that pulmonate snails such as those
of the genera Radix and Planorbarius exhibit a great
diversity in shell morphology with extremely
homogeneous anatomical traits . Varying environmental
factors seem to affect the morphological features
resulting in variations and making it difficult to identify the
species accurately on the basis of external features.
© The Author(s). 2017 Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0
International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and
reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to
the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver
(http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
Additionally, the evolutionary relationships among
different molluscan classes and within some major clades
are still unclear, due to the limited taxon sampling [7, 8].
Owing to the unique features such as maternal
inheritance, lack of extensive recombination, a relatively high
evolutionary rate and abundantly available marker types
, mitochondrial (mt) genomes have been widely used
for species identification, population genetics and
evolutionary relationships studies in metazoans including
pulmonates [4, 5, 8]. In general, the metazoan mt genome is
a single circular DNA molecule of about 14–17 kb in size,
typically containing 37 genes : 13 protein-coding genes
(PCGs) [cytochrome c oxidase subunits 1–3 (cox1-cox3),
apocytochrome b (cytb), ATPsynthase subunits 6 and 8
(atp6 and atp8), NADHdehydrogenase subunits 1–6 and
4 L (nad1-6 and nad4L)], two ribosomal RNAs [small and
large subunit ribosomal RNA (rrnS and rrnL)] and 22
transfer RNA (tRNA). Typically, there are few non-coding
regions containing repeat elements or pseudogenes, which
are the usual source of size variations in metazoan
mtDNA. With recent methodological advances,
particularly the next generation sequencing technologies and the
associated cost reduction in DNA sequencing , a
growing number of complete mt genome sequences are
available for mollusks in the GenBank database.
Pulmonate snails like R. swinhoei and P. corneus are pathogen
carriers and research into their basic biology has medical
implications. Radix swinhoei serves as the major
intermediate host of pathogens such as Fasciola hepatica,
Trichobilharzia paoi, T. physellae, T. ocellata,
Echinostoma revolutum, E. hortense, Orientobilharzia
turkestanicum, Angiostrongylus cantonensis, Cercaria ohiensis,
Plagiorchis muris, Euparyphium ilocanum,
Echinoparyphium recurvatum, Diplostomum niedashui and D.
hupensis in China, Japan, Thailand, India and Vietnam [3,
11]. Planorbarius corneus is the dominant intermediate
host snail for the transmission of Prosthogonimus ovatus,
Apatemon gracilis, Hypoderaeum conoideum, Syngamus
trachea and Typhlocoelum sisowi worldwide [3, 12]. In this
study, we used the Illumina-based paired-end sequencing
 to report novel complete mt genomes of the
freshwater snails Radix swinhoei and Planorbarius corneus,
belonging to the families Lymnaeidae and Planorbidae,
respectively. Mt gene rearrangement in pulmonate
gastropods is of key interest to scientists from the perspective of
understanding evolution and genome diversification.
Further sequence analysis of the two snail species under
investigation revealed novel gene rearrangements involving both
protein-coding and tRNA genes. Together with other
published complete mt genomes of heterobranchs gastropods,
we reconstructed the phylogenetic relationship using the
amino acid sequences of the 13 protein-coding genes with
two different computational algorithms (maximum
likelihood and Bayesian inference analysis). These data would
provide valuable information not only for phylogenetic
studies but also for the development of useful genetic
markers for stock management and molecular
epidemiological studies of parasites.
Specimen collection and DNA extraction
One adult individual of each R. swinhoei and P. corneus
was collected from the Aquatic Invasive Risk Assessment
Center, Pearl River Fisheries Research Institute Chinese
Academy of Fishery Sciences (23°04′04.05″N, 113°13′
06.97″E) in Guangzhou, Guangdong Province, China. The
specimens were washed in physiological saline, identified
morphologically according to existing descriptions of
mollusc shape , fixed in 70% (v/v) ethanol and stored at
-20 °C. Total genomic DNA was isolated from each
species using approximately 30 mg of fresh foot tissue with
OMEGA EZNA Mollusc DNA kit following the
manufacturers’ instructions. Total DNA was eluted in sterile
deionized water and was stored at -20 °C.
Mitochondrial genomes sequencing, assembly and
Paired-end libraries (500 bp) using TruSeq DNA Sample
prep kit were prepared following the Illumina
instructions. The size-selected, adapter-modified DNA
fragments were PCR-amplified using PCR primers following
the protocol: polymerase activation (98 °C for 2 min)
followed by 10 cycles (denaturation at 98 °C for 30 s,
annealing at 65 °C for 30 s, and extension at 72 °C for
60 s) with a final, 4 min extension at 72 °C. DNA
libraries were purified by magnetic beads and quantified by
real time quantitative PCR (RT-PCR).
Sequencing using Hiseq 2500 plate resulted in 1.23 Gb
(R. swinhoei) and 1.44 Gb (P. corneus) high quality reads,
containing 18,322 reads and 4514 reads of mitochondrion,
respectively (Additional file 1: Table S1). Pair-End 100 bp
read length of Illumina reads were analyzed. Reads that
contained adapters were trimmed, and low quality reads
which have more than 3 “N” base were removed. The first
assembly using the chloroplast and mitochondrion
assemble (CMA) V1.1.1 software (Guangzhou SCGene Co., Ltd)
was based on overlap with the mt genomes of related
species and paired-end relationships. The assembled
complete mt genomes were tested for completeness and
preciseness through paired-end read mapping back to the
genome. For P. corneus, an uncertain region (nt 6275–
6731) was amplified using the primers (F: 5′-ATG TGG
GTT GTC AAT TAT CTG GT-3′; R: 5′-GCT ATA ACT
AAG CTA TTG GGC TC-3′). The PCR reactions were
prepared with a 40 μl total volume as follows: 20 μl 2×
Taq master Mix (100 μmol/l) (GC gene), 2.0 μl of each
primer (10 μmol/l), 15 μl ddH2O, and 1 μl DNA sample
(0.2 ~ 0.5 μmol/l). The following PCR cycle was used for
all fragment amplifications: an initial denaturation at
94 °C for 4 min; 35 cycles of 94 °C for 30 s
(denaturation), 55 °C for 30 s (annealing), and 72 °C for
2 min (extension); followed by a final extension at
72 °C for 10 min. PCR products were examined using
1% agarose gel electrophoresis to validate the
amplification efficiency and were sequenced using an ABI
377 (Applied Biosystems) automated DNA sequencer
using the same primers (one primer at a time) as that
for PCR. After de novo assembly and functional
annotation, 13 protein-coding genes, rDNA genes, and
tRNAs of mt genome were found, and compared with
the two known complete mt genome of species of the
family Lymnaeidae: Galba pervia (JN564796)  and
Radix balthica (HQ330989) . The related data
were deposited in the National Center for
Biotechnology Information (NCBI) Biosample databases with the
accession numbers: SRS781941 (R. swinhoei) and
SRS781939 (P. corneus). The two complete mt
genomes were deposited in the GenBank database:
KP279638 (R. swinhoei) and KP279639 (P. corneus).
The complete mt genome sequence of R. swinhoei and
P. corneus was aligned with other pulmonate complete
mt sequences obtained from GenBank (Additional file 1:
Table S2) by Clustal ×1.83 , using default
parameters, and following the guidelines of the SeaView
software . Codon usage and nucleotide composition
statistics were computed using Molecular Evolutionary
Genetics Analysis (MEGA) 6 . The protein-coding
regions were identified by Basic Local Alignment Search
Tool (BLAST; blastn, tblastx) from National Center for
Biotechnology Information (NCBI) database and by
using the MITOS WebServer BETA
(http://bloodymary.bioinf.uni-leipzig.de/mitos/index.py) . The transfer
RNA genes were annotated using tRNA scan-SE v.1.21
(http://lowelab.ucsc.edu/tRNA scan -SE) with Search
Mode = “Eufind tRNA- Cove”, Genetic Code =
“Invertebrate Mito”, and Cove score cut-off = 0.1, and the
software ARWEN (http://126.96.36.199/ARWEN/) . The
map of the species was visualized using the Genome Vx
online tool (http://wolfe.ucd.ie/GenomeVx/) . Repeat
sequences were found using Spectral Repeat Finder v1.1
. Strand asymmetry was calculated using the
formulas: AT skew = (A-T)/(A + T) and GC skew = (G-C)/(G +
C) . Codon usage and building block distributions
were determined gene-wise for all protein-coding genes,
and merged using MEGA6.06 and statistical package R.
Statistical analyses of distribution and codon usage
heatmaps were generated using package R as well . The
stem-loop secondary structures of the non-coding
regions were predicted using the default parameters
under RNA folding option in the Mfold Server (http://
www.bioinfo.rpi.edu/applications/mfold/) . To
conduct pair-wise comparison of the mt gene order of R.
swinhoei and P. corneus with that of Reishia clavigera
(name currently accepted for Thais clavigera ) as the
standard gene pattern of molluscan mt genomes , we
used CREx the program . CREx is an efficient
software suite which could analyze complex genome
rearrangements scenarios in the gene order of a pair of taxa
and determine the most parsimonious steps required for
the rearrangement. In terms of rearrangement
mechanism, the software can handle transpositions, reverse
transpositions, reversals, and tandem duplication -
random loss (TDRL) events among others. The analysis was
performed by applying the common interval parameters
for distance measurement and by using only
proteincoding and ribosomal RNA genes. The more variable
tRNAs were excluded from the analysis. Linear mt
genome comparison of R. swinhoei, P. corneus and related
species was performed using EasyFig2.1 (BLASTn, default
setting) . The graphical map was visualized with the
CGView Comparison Tool (CCT) .
Phylogenetic analyses were performed using Bayesian
inference (BI) and maximum likelihood (ML) methods.
For the best-fit models of evolution for the amino acid
sequence datasets (13 protein-coding genes) was selected
ProtTest 2.4 . BI was performed on combined
database using MrBayes version 3.1.2 for 10,000,000
generations with a random starting tree. Four independent
Markov chains were simultaneously run for ten million
generations with a heating scheme (temp = 0.2) .
Trees were sampled every 100 generations (sample-freq
= 100) with the first 25% of the generations were
discarded as ‘burn-in’ and the remaining generations were
used to compute the consensus tree. Stationarity was
considered to be reached when the average standard
deviation of split frequencies was below 0.01. ML analyses
were conducted using PhyML 3.0 with 1000
bootstrapping based on the MtArt + I + G model . The
phylogenetic trees were drawn using the Evolview (http://
Results and discussion
Structural features of the mitocondrial genome
The complete mt genomes of R. swinhoei and P. corneus
are 14,271 bp and 13,687 bp in length, respectively. The
mt genome length of the two species of snail are
comparable to that of other sequences of pulmonates
(Additional file 1: Table S2). Mt genome of both R.
swinhoei and P. corneus is a circular double-stranded DNA
molecule, containing a total of 37 genes typically found
in metazoans. These 37 genes belong to the following
categories: 13 PCGs (cox1-3, nad1-6, nad4L, atp6, atp 8
Fig. 1 Gene maps and organization of the complete mitochondrial genomes of Radix swinhoei (a) and Planorbarius corneus (c). The outer and the
inner circles represent the positive (H-strand) and the negative (L-strand) strand, respectively. The tRNA genes are named using single-letter amino acid
abbreviations as shown in the Tables (b, d). aIndicates that the gene is encoded by H or L strand; b Intergenic nucleotides, indicates the number of
nucleotides separating a gene from the one upstream of it; negative numbers indicate an overlap between the adjacent genes; cIncomplete termination
codon, which might probably be extended by post-transcriptional adenylation; dRepresent the presence of repeat masker
and cytb), 1 rrnS, 1 rrnL and 22 tRNAs (Fig. 1). A high
variation in nucleotide composition of pulmonate mt
genomes has been reported [2, 4–6, 9, 14]. The variation
of overall A + T content ranges from 54.76% (Ovatella
vulcani) to 77.0% (Succinea putris), with an average value of
65.5% (Additional file 1: Table S2). The A + T content of
R. swinhoei is 69.45% and of P. corneus is 72.66%,
corresponding well with that of related species. The high A + T
content is also reflected in the individual PCGs, with the
values especially higher for nad6 gene (77.1%) for R.
swinhoei, and cox2 gene (77.4%) for P. corneus (Adittional file
1: Table S3; Additional file 1: Table S4).
There are small variations in the AT- and GC-skews in
pulmonate mt genomes. Due to the strand asymmetry
(strand compositional bias) , the AT skews of the
whole mt genome ranges from -0.210 (Siphonaria gigas)
to -0.073 (Albinaria caerulea), while the GC skew values
are between 0.047 (P. corneus) and 0.215 (S. gigas)
(Additional file 1: Table S1). As with other pulmonate species,
a similar AT and GC skews were detected in the mt
genomes of both R. swinhoei and P. corneus (Additional
file 1: Figure S1). Interestingly, the mt genome AT and
GC skew values are similar between the two snails here
studied. However, individual PCGs showed different and
variable AT and GC skews in the R. swinhoei and P.
corneus mt genomes (Additional file 1: Figure S2). AT skews
were negative for most of the protein-coding genes
except for rrnS (0.1) and rrnL (0.001) in R. swinhoei, and
for nad2 (0.07) in P. corneus. On the other hand,
GCskews were positive generally, with negative values for
atp6 gene (-0.1) in R. swinhoei, and for four genes, atp6
(-0.02), cox3 (-0.21), rrnS (-0.15); rrnL (-0.01) in P.
corneus (Additional file 1: Table S3, Table S4, Figure S3).
The nucleotide composition bias and skew may be
caused by the selection-mutation-drift equilibrium of
molecular evolution .
To better visualize the nucleotide identity in
pulmonate mt genomes, we generated the graphical identity
map using the CGView comparison tool (CCT). Taking
gene order into account, an easy to track comparable
graphic gene identity map was generated (Fig. 2). High
conservation in the cox genes was observed with cox1
showing the highest similarity among several pulmonate
species. On the other hand, nad genes were most
variable with nad4L showing the maximal variation.
Pairwise comparisons of the concatenated amino acid
sequences revealed the overall mt genome similarity of
40.3–91.8% among pulmonate snails. This variation can
be attributed to the rapid rate of mtDNA evolution.
Protein-coding genes (PCGs) and codon usage patterns
The full set of 13 PCGs, typically found in pulmonate
species, were identified in the mt genome of R. swinhoei
and P. corneus. Inferred initiation and termination
codons from each protein-coding gene are shown in Fig. 1.
Of the 13 PCGs, ten genes were found to initiate with
the ATN codon, whereas three started with TTG in both
R. swinhoei and P. corneus. These data are in accordance
with previous findings from different gastropod clades
[2, 4, 9, 14]. Most of PCGs were inferred to use TAA/
TAG as stop codons except for T (P. corneus: cox1 and
cox2), which frequently occurred in protein-coding genes
of most gastropod mt genomes [5, 6, 8, 13]. The
incomplete stop codon was thought to be complemented via
post-transcription alpolyadenylation .
The various codon families and the Relative
Synonymous Codon Usage (RSCU) for PCGs in R.
swinhoei and P. corneus and their related species are
summarized in Fig. 3. The total number of codons
for all protein-coding genes in the mt genome of R.
swinhoei and P. corneus were found to be 3443 and
3595, respectively (Additional file 1: Figure S4, Table
S5, and Table S6). These numbers are distinctly small
in comparison to that for G. pervia (3655) . A bias
towards T-rich codon was observed in the
proteincoding genes, which may be attributed to the high
percentage of Thymine in the mt genome of R.
swinhoei and P. corneus. TTA coding for Leucine
(Leu) was the most frequently used codon (203 times)
in R. swinhoei (Additional file 1: Table S5) while TTT
coding phenylalanine (336 times) was maximally
represented in P. corneus (Additional file 1: Table S6).
The codon families also exhibited a difference in their
usage pattern among the two species, with Leu being
the most frequent amino acid in R. swinhoei (16.94%)
and Ser in P. corneus (11.10%) (Additional file 1:
Figure S4). We also found a similarity in codon usage of
R. swinhoei and P. corneus to that of the suborder
Hygrophila. Although no major difference in codon
usage was observed, the codons varied between
different species and different genes. For example, a minor
variation in the frequency of the codon TTA was
observed among the investigated species.
Ribosomal and transfer RNA genes, and non-coding
Similar to most of the other pulmonate mt genomes [4, 5],
the location of rrnL is between tRNA-Val (V) and
tRNALeu (L1), while that of rrnS is between tRNA-Glu (E) and
tRNA-Met (M) in both R. swinhoei and P. corneus mt
genomes (Fig. 1). The length of rrnL and rrnS in the R.
swinhoei mt genome is 1370 bp and 827 bp and in P. corneus
mt genome is 884 bp and 647 bp, respectively. The A + T
contents of the rrnL and rrnS of both R. swinhoei (72.8 and
70.4%) and P. corneus (70.8 and 73.3%) were lower
compared to that of G. pervia (rrnL: 74.93%; rrnS: 72.09%) .
Additionally, sequence alignment of R. swinhoei and P.
corneus demonstrated sequence similarities for rrnL (67.6%)
and rrnS (64.7%).
Both R. swinhoei and P. corneus contained 22 tRNA
genes, ranging in size from 49 bp for both tRNA-Cys
and tRNA-Glu (in P. corneus) to 75 bp for tRNA-Val (in
R. swinhoei) with variations mainly arising from
differences in stem and loop sizes of dihydrouridine (DHU)
and TΨC. Most of the tRNA genes were predicted to
have the typical cloverleaf secondary structure, except
for tRNA-Gly (G), tRNA-Ser (S1) (AGN) and tRNA-Ser
(S2) (UCN) in R. swinhoei (Fig. 4a), and tRNA-Cys (C),
Fig. 2 Graphical map showing nucleotide identity among the mitochondrial genes from Radix swinhoei, Planorbarius corneus and other pulmonate
species. Gene specific identity was obtained by BLAST searches. The map is visualized by using the CGView comparison tool (CCT), which arranges
BLAST result in an order where the sequence that is most similar to the reference (in this case R. swinhoei), is placed closer to the outer edge of
tRNA-Ser (S1) (AGN) and tRNA-Ser (S2) (UCN) in P.
corneus (Fig. 4b). The tRNA-Gly (G) of R. swinhoei as well as
tRNA-Cys (C) of P. corneus harbor a simple loop in the
TΨC stem, while tRNA-Ser (S1) (AGN) and tRNA-Ser
(S2) (UCN) of both R. swinhoei and P. corneus harbor a
simple loop in the dihydrouridine (DHU) arm.
Furthermore, tRNA rearrangements were predicted to occur in
the R. swinhoei and P. corneus mt genomes. Such similar
tRNA rearrangements have been reported in multiple
divergent taxa, like G. pervia  and R. balthica .
As in most pulmonate snail species, both mt genomes
contained a number of unassigned nucleotides, with the
number ranging from 220 in P. corneus (1.6% of the
genome) to 294 in R. swinhoei (2.1% of the genome). There
are more than 30 non-coding regions throughout R.
swinhoei (49, 45 and 36 bp in length) and P. corneus
(112, 71 and 42 bp in length). The longest non-coding
region (49 bp) in R. swinhoei, located between cox3 and
tRNA-Ile gene, has a high A + T content (89.8%) and
two stem-loop secondary structures (Fig. 4c), whereas
the longest non-coding region in P. corneus (112 bp) lies
between tRNA-Val and rrnL gene with a high A + T
content (89.3%) and three stem-loop secondary structures
(Fig. 4d). Although the functions of most of these
nonFig. 3 Codon usage distributions of all merged protein-coding genes among various pulmonate species. Color key: Color value around 1 (green region)
implies that a codon was observed about 1 time in 60 codons indicating a uniform distribution. Higher (red) and lower (blue) values indicate deviations
of uniform distribution (as factors). This map also helps us visualize the implicit variability in amino acid distribution. a Hierarchical clustering (average
linkage method) of codon patterns (y-axis) approximately groups together the examined families. Closely related species (i.e., from the same tribe)
generally have highly similar codon usage. Clustering of codon frequencies (x-axis) facilitates the visual identification of deviations. b The evaluation of
relative synonymous codon usage
coding regions remain unclear, the longest regions from
both the species most likely are “putative control
regions” owing to their sequence length and the presence
of characteristic stem loop structure.
Comparison of mitochondrial gene order
The gene order of the two mt genomes under
investigation were compared to each other and to that of R.
clavigera as a representative of the ancestral mollusc gene
Fig. 4 (See legend on next page.)
(See figure on previous page.)
Fig. 4 Analysis of possible secondary structure of mitochondrial tRNA genes and non-coding regions from Radix swinhoei and Planorbarius corneus. a,
b Putative secondary structures of three representative tRNA genes identified in the mitochondrial genomes of R. swinhoei (a) and P. corneus (b); c, d
Stem-loop secondary structures of two non-coding regions in the mt genomes of R. swinhoei (49 bp and 45 bp) (c) and P. corneus (112 bp and 71 bp)
(d). Bars indicate Watson-Crick base pairings; dots indicate canonical base pairing between G and U nucleotides in RNA
order within gastropods . The gene order of the mt
genomes of R. swinhoei and P. corneus is significantly
different, with positional mismatches for both
proteincoding and tRNA genes (Additional file 1: Figure S5).
These findings are consistent with studies demonstrating
high diversity in gene arrangements in pulmonate species,
such as G. pervia  and R. balthica . Some of the
studies provide novel insights into the mechanism of
mtDNA rearrangement in certain gastropod species
[14, 37, 38]. Our results further support the view that
tRNAs are involved in more frequent rearrangements
than protein-coding genes and ribosomal DNA in
metazoan mt genomes . By comparing closely
related species with different gene orders, the mt gene
rearrangements could be usually explained by four
possible models (the recombination model, tandem
duplication and random loss (TDRL) model, tandem
duplication and non-random loss (TDNL) model and
tRNA miss-priming model) based on three
mechanisms of gene arrangement (shuffling, translocation,
and inversion) [39, 40]. However, for all
rearrangement scenarios, tRNAs were not compared due to
their higher variability in location.
We furthermore used the software CREx to reflect the
genomic rearrangement history of R. swinhoei and P.
corneus. The results effectively presumed that the gene
rearrangement of R. swinhoei was postulated as follows.
The first step was three times of continuous reversal: a
reversal of 14 PCG genes except for nad5, a reversal of
“nad4L-nad4” and a reversal of
“cox2-cox1-nad1-nad2rrnL” and “Cytb” (see three TDRL models in Fig. 5).
Then, we performed the same analysis to presume the
gene rearrangement form R. clavigera to P. corneus. The
result showed that the first step was reverse
transposition of cox2 gene. In the second of step, there were at
least four putative reversals, including two reversals of
14 PCG genes except for nad5, a reversal of 13 PCG
genes except for nad5 and nad6, a reversal of
“nad3cox3- cox2-atp6- atp8-cox1-nad2- nad1-rrnL-rrnS”. The
third step included two TDRL. At last, there was a
transposition of rrnL gene.
The phylogenetic relationships of pulmonate species
based on concatenated amino acid sequence datasets
using BI and ML analyses were reconstructed. The two
50% consensus trees had a similar topology with
wellsupported branches for major clades (Fig. 6). In the tree,
all Panpulmonate species were clustered with high
statistical support. Among the families represented by more
than one species, the Helicidae, Planorbidae,
Siphonariidae and Onchidiidae were recovered as monophyletic,
while the Ellobiidae and Lymnaeidae were
nonmonophyletic due the recovered position of Myosotella
myosotis and R. swinhoei, respectively. Some authors also
recovered Ellobiidae as paraphylelic using the complete
mt genomes [41–43] or partial genes , mainly due to
the position of Pedipes pedipes or M.myosotis.
Nevertheless, other phylogenies, such as those of Dayrat et al. 
and Romero et al.  using nuclear and mt genes
rendered Ellobiidae as monophyletic. Meanwhile, the
taxonomic position of R. swinhoei should be revised. As
pointed out by Lawton et al.  molecular
identification was the only reliable method to identify Radix
species and other Lymnaeidae since shell and other
anatomical features are morphologically plastic and most
species share morphological characters as a result of
convergent adaptations to shared limnic environments.
The results also revealed that the families Lymnaeidae
and Planorbidae are closely related with high statistical
support, and the data obtained basically agreed with
those of previous phylogenetic analyses based on
complete mt genomes [4, 41–43, 47]. The taxa
Sacoglossa (families Volvatellidae and Placobranchidae) and
Siphonarioidea (family Siphonariidae) were recovered as
sister clades, indicating closely relationships (already
noted by Grande et al. ). Likewise, the Trimusculidae
and Ellobiidae also showed sister-group relationships,
previously pointed out by White et al. .
Additional complete mt genomes are needed for
pulmonate snails (especially from missing lineages) in order
to resolve the phylogenetic framework of this diverse
group of gastropods to further understanding its evolution
and obtain new information about the detailed processes
and mechanisms of mt genome rearrangements.
In this study, the characterization of the complete mt
genomes of R. swinhoei and P. corneus, both encoding
all the thirty-seven genes typical for pulmonates,
revealed considerable interspecific differences in length
and sequence composition. These genes were arranged
in the same order as that of the proposed ancestral
gastropod. However, the most remarkable feature was
that, unlike in other pulmonate snails, a novel gene
Fig. 5 Putative gene rearrangement events in gastropod mitochondrial genomes. Reishia clavigera was used as the standard ancestral gene pattern.
Rearranged genes are indicated by different colors. Only the rearrangements for protein-coding and rRNA genes are taken into consideration. Four types
of putative rearrangement events are possible in this context: reversal, reverse transposition, transposition and tandem duplication and random loss
(TDRL). In the step-by-step scheme, the intermediate statuses are used to show different types of gene rearrangement events and the rearrangement
process. a Putative gene rearrangement events from R. clavigera to Radix swinhoei mt genome. b Putative gene rearrangement events from R. clavigera
to Planorbarius corneus mt genome
Fig. 6 Phylogenomic analysis of concatenated proteins in mitochondrial genomes. The concatenated amino acid sequences of 13 protein-coding
genes were analyzed utilizing Bayesian analysis (BI) and maximum likelihood (ML), using Micromelo undata and Nautilus macromphalus as the outgroup
arrangement was observed. This study also provides an
idea about novel mt genetic markers for species
identification and population genetics of freshwater pulmonates
and also has implications for the diagnosis, prevention
and control of Fasciola spp. infection in hosts.
Additional file 1: Table S1. Summary of Radix swinhoei and
Planorbarius corneus using Illumina sequencing. Table S2. General
characteristics of the mitochondrial genomes of various members of
pulmonate gastropods. Table S3. Nucleotide composition and AT- and
GC-skews of the mitochondrial protein-coding and ribosomal RNA genes
in the complete Radix swinhoei mt genome. Table S4. Nucleotide
composition and AT- and GC-skews of the mitochondrial protein-coding
and ribosomal RNA genes in the complete Planorbarius corneus mt
genome. Table S5. Codon usage of Radix swinhoei mitochondrial
protein-coding genes. Table S6. Codon usage of Planorbarius corneus
mitochondrial protein-coding genes. Figure S1. Comparison of AT and
GC skews among the 30 pulmonate species in Table S2. Circles and
triangles separately represent AT and GC skews of the complete
mitochondrial genomes. Figure S2. Gene specific strand composition of
mitochondrial genome in Radix swinhoei and Planorbarius corneus. Figure S3.
Graphical summary of nucleotide composition across complete
mitochondrial genomes. Figure S4. Percentage of synonymous codon usage
for each amino acid in Radix swinhoei and Planorbarius corneus
mitochondrial protein-coding genes. Figure S5. Linear comparison of the
gene distribution pattern of mitochondrial genomes between Radix
swinhoei and Planorbarius corneus. (DOC 1430 kb)
ATP6 and ATP8: ATPase subunits 6 and 8; cox1-cox3: cytochrome c oxidase
subunits1-3; Cytb: cytochrome b; ND1-ND6 and ND4L: NADH dehydrogenase
subunits 1-6 and 4 L; NGS: Next-generation sequencing; PCG: Protein-coding
gene; rRNA: ribosomal RNA; rrnL: Large rRNA subunit (gene); rrnS: Small rRNA
subunit (gene); trn: Transfer RNA
Availability of data and materials
The datasets generated during the current study were submitted to the
National Center for Biotechnology Information (NCBI) Biosample databases
with the accession numbers: SRS781941 (R. swinhoei) and SRS781939 (P.
corneus). The two complete mt genomes were deposited in the GenBank
database: KP279638 (R. swinhoei) and KP279639 (P. corneus).
XDM designed the study, analyzed the data, and drafted the whole
manuscript. YXY, YL, DL, MX analyzed the bioinformatic data and participated
in the manuscript revision. HW, GEG and HMS collected samples, assisted
with data analysis. XM and YH co-designed the experiments and obtained
the funds. All authors read and approved the final manuscript.
This study did not involve the use of endangered or protected species, and
was carried out in strict accordance with the recommendations in the Guide
for the Care and Use of Laboratory Animals of Pearl River Fisheries Research
Institute, Chinese Academy of Fishery Sciences. All experiments were
conducted maintaining current China laws. The protocol was approved by
the Committee on the Ethics of Animal Experiments of Pearl River Fisheries
Research Institute, Chinese Academy of Fishery Sciences.
1. Jörger KM , Stöger I , Kano Y , Fukuda H , Knebelsberger T , Schrödl M. On the origin of Acochlidia and other enigmatic euthyneuran gastropods, with implications for the systematics of Heterobranchia . BMC Evol Biol . 2010 ; 10 : 323 .
2. Dayrat B , Conrad M , Balayan S , White TR , Albrecht C , Golding R. Phylogenetic relationships and evolution of pulmonate gastropods (Mollusca): new insights from increased taxon sampling . Mol Phylogenet Evol . 2011 ; 59 ( 2 ): 425 - 37 .
3. Liu YY , Zhang WZ , Wang YX , Wang EY . Chinese economics animal: freshwater snail . Sci Press (China). 1979 ; 1 : 54 - 6 .
4. Liu GH , Wang SY , Huang WY , Zhao GH , Wei SJ , Song HQ , et al. The complete mitochondrial genome of Galba pervia (Gastropoda: Mollusca), an intermediate host snail of Fasciola spp . PLoS One . 2012 ; 7 ( 7 ): e42172 .
5. Nolan JR , Bergthorsson U , Adema CM . Physella acuta: atypical mitochondrial gene order among panpulmonates (Gastropoda) . J Moll Stud . 2014 ; 80 ( 4 ): 388 - 99 .
6. Pfenninger M , Cordellier M , Streit B. Comparing the efficacy of morphologic and DNA-based taxonomy in the freshwater gastropod genus Radix (Basommatophora, Pulmonata) . BMC Evol Biol . 2006 ; 6 ( 1 ): 100 .
7. Zou S , Li Q , Kong L. Additional gene data and increased sampling give new insights into the phylogenetic relationships of Neogastropoda, within the caenogastropod phylogenetic framework . Mol Phylogenet Evol . 2011 ; 61 ( 2 ): 425 - 35 .
8. Kocot KM , Halanych KM , Krug PJ . Phylogenomics supports Panpulmonata: opisthobranch paraphyly and key evolutionary steps in a major radiation of gastropod molluscs . Mol Phylogenet Evol . 2013 ; 69 ( 3 ): 764 - 71 .
9. Wolstenholme DR . Animal mitochondrial DNA: structure and evolution . Int Rev Cytol . 1992 ; 141 : 173 - 216 .
10. Rasmussen DA , Noor MA. What can you do with 0.1× genome coverage? A case study based on a genome survey of the scuttle fly Megaselia scalaris (Phoridae) . BMC Genomics . 2009 ; 10 : 382 .
11. Guo YH , Wang CM , Luo J , He HX . Intermediate host of main parasites: mollusks distributed in Beijing region . Chin J Vector Bio Control . 2009 ; 20 ( 5 ): 449 - 54 .
12. Chen YX , Zhang W , Tian M , Chen YJ , Wang WL , Zhang NG . Investigation of snails transmitting parasites diseases in Yunnan Province . J Pathogen Biol . 2009 ; 4 ( 3 ): 211 - 4 .
13. Williams ST , Foster PG , Littlewood DT . The complete mitochondrial genome of a turbinid vetigastropod from MiSeq Illumina sequencing of genomic DNA and steps towards a resolved gastropod phylogeny . Gene . 2014 ; 533 ( 1 ): 38 - 47 .
14. Feldmeyer B , Hoffmeier K , Pfenninger M. The complete mitochondrial genome of Radix balthica (Pulmonata, Basommatophora), obtained by low coverage shot gun next generation sequencing . Mol Phylogenet Evol . 2010 ; 57 ( 3 ): 1329 - 33 .
15. Thompson JD , Gibson TJ , Plewniak F , Jeanmougin F , Higgins DG . The CLUSTAL_ X windows interface: flexible strategies for multiple sequence alignment aided by quality analysis tools . Nucleic Acids Res . 1997 ; 25 ( 24 ): 4876 - 82 .
16. Gouy M , Guindon S , Gascuel O. SeaView version 4: a multiplatform graphical user interface for sequence alignment and phylogenetic tree building . Mol Biol Evol . 2010 ; 27 ( 2 ): 221 - 4 .
17. Tamura K , Stecher G , Peterson D , Filipski A , Kumar S. MEGA6: Molecular evolutionary genetics analysis version 6.0. Mol Biol Evol . 2013 ; 30 ( 12 ): 2725 - 9 .
18. Bernt M , Donath A , Jühling F , Externbrink F , Florentz C , Fritzsch G , et al. MITOS: improved de novo metazoan mitochondrial genome annotation . Mol Phylogenet Evol . 2013 ; 69 ( 2 ): 313 - 9 .
19. Laslett D , Canback B. ARWEN: a program to detect tRNA genes in metazoan mitochondrial nucleotide sequences . Bioinformatics . 2008 ; 24 ( 2 ): 172 - 5 .
20. Conant GC , Wolfe KH . Genome Vx: simple web-based creation of editable circular chromosome maps . Bioinformatics . 2008 ; 24 ( 6 ): 861 - 2 .
21. Sharma D , Issac B , Raghava GP , Ramaswamy R. Spectral repeat finder (SRF): identification of repetitive sequences using Fourier transformation . Bioinformatics . 2004 ; 20 ( 9 ): 1405 - 12 .
22. Perna NT , Kocher TD . Patterns of nucleotide composition at fourfold degenerate sites of animal mitochondrial genomes . J Mol Evol . 1995 ; 41 ( 3 ): 353 - 8 .
23. Charif D , Lobry JR . SeqinR 1.0-2: a contributed package to the R project for statistical computing devoted to biological sequences retrieval and analysis . In: Bastolla U, Porto M , Roman HE , Vendruscolo M, editors. Structural approaches to sequence evolution: molecules , networks, populations. New York : Springer Verlag ; 2007 . p. 207 - 32 .
24. Zuker M. Mfold web server for nucleic acid folding and hybridization prediction . Nucleic Acids Res . 2003 ; 31 ( 13 ): 3406 - 15 .
25. Houart R : Reishia clavigera (Küster, 1860 ). In: MolluscaBase ( 2015 ). WoRMS database : http://www.marinespecies. org/aphia.php?p=taxdetails&id=397003.
26. Ki JS , Lee YM , Jung SO , Horiguchi T , Cho HS , Lee JS . Mitochondrial genome of Thais clavigera (Mollusca: Gastropoda): affirmation of the conserved, ancestral gene pattern within the mollusks . Mol Phylogenet Evol . 2010 ; 54 ( 3 ): 1016 - 20 .
27. Bernt M , Merkle D , Ramsch K , Fritzsch G , Perseke M , Bernhard D , et al. CREx: inferring genomic rearrangements based on common intervals . Bioinformatics . 2007 ; 23 ( 21 ): 2957 - 8 .
28. Sullivan MJ , Petty NK , Beatson SA . Easyfig: a genome comparison visualizer . Bioinformatics . 2011 ; 27 ( 7 ): 1009 - 10 .
29. Grant JR , Arantes AS , Stothard P. Comparing thousands of circular genomes using the CGView Comparison Tool . BMC Genomics . 2012 ; 13 : 202 .
30. Abascal F , Zardoya R , Posada D. ProtTest: selection of best-fit models of protein evolution . Bioinformatics . 2005 ; 21 ( 9 ): 2104 - 5 .
31. Ronquist F , Huelsenbeck JP . MrBayes 3: Bayesian phylogenetic inference under mixed models . Bioinformatics . 2003 ; 19 ( 12 ): 1572 - 4 .
32. Guindon S , Gascuel O. A simple, fast, and accurate algorithm to estimate large phylogenies by maximum likelihood . Syst Biol . 2003 ; 52 ( 5 ): 696 - 704 .
33. Zhang H , Gao S , Lercher MJ , Hu S , Chen WH. EvolView, an online tool for visualizing, annotating and managing phylogenetic trees . Nucleic Acids Res . 2012 ; 40 (Web Server issue): W569 - 72 .
34. McLean MJ , Wolfe KH , Devine KM . Base composition skews, replication orientation, and gene orientation in 12 prokaryote genomes . J Mol Evol . 1998 ; 47 ( 6 ): 691 - 6 .
35. San Mauro D , Gower DJ , Zardoya R , Wilkinson M. A hotspot of gene order rearrangement by tandem duplication and random loss in the vertebrate mitochondrial genome . Mol Biol Evol . 2006 ; 23 ( 1 ): 227 - 34 .
36. Anderson S , Bankier AT , Barrell BG , de Bruijn MHL , Coulson AR , Drouin J , et al. Sequence and organization of the human mitochondrial genome . Nature . 1981 ; 290 ( 5806 ): 457 - 65 .
37. Grande C , Templado J , Zardoya R. Evolution of gastropod mitochondrial genome arrangements . BMC Evol Biol . 2008 ; 8 ( 1 ): 61 .
38. Yamazaki N , Ueshima R , Terrett JA , Yokobori S , Kaifu M , Segawa R , et al. Evolution of pulmonate gastropod mitochondrial genomes: comparisons of gene organizations of Euhadra, Cepaea and Albinaria and implications of unusual tRNA secondary structures . Genetics . 1997 ; 145 ( 3 ): 749 - 58 .
39. Macey JR , Larson A , Ananjeva NB , Fang Z , Papenfuss TJ . Two novel gene orders and the role of light-strand replication in rearrangement of the vertebrate mitochondrial genome . Mol Biol Evol . 1997 ; 14 ( 1 ): 91 - 104 .
40. Lunt DH , Hyman BC . Animal mitochondrial DNA recombination . Nature . 1997 ; 387 ( 6630 ): 247 .
41. Gaitan-Espitia JD , Nespolo RF , Opazo JC . The complete mitochondrial genome of the land snail Cornu aspersum (Helicidae: Mollusca): intra-specific divergence of protein-coding genes and phylogenetic considerations within Euthyneura . PLoS One . 2013 ; 8 ( 6 ): e67299 .
42. Price MR , Forsman ZH , Knapp I , Hadfield MG , Toonen RJ . The complete mitochondrial genome of Achatinella mustelina (Gastropoda: Pulmonata: Stylommatophora) . Mitochondrial DNA Part B: Resour . 2016 ; 1 ( 1 ): 175 - 7 .
43. Romero PE , Weigand AM , Pfenninger M. Positive selection on panpulmonate mitogenomes provide new clues on adaptations to terrestrial life . BMC Evol Biol . 2016 ; 16 : 164 .
44. Klussmann-Kolb A , Dinapoli A , Kuhn K , Streit B , Albrecht C. From sea to land and beyond - new insights into the evolution of euthyneuran Gastropoda (Mollusca) . BMC Evol Biol . 2008 ; 8 : 57 .
45. Romero PE , Pfenninger M , Kano Y , Klussmann-Kolb A. Molecular phylogeny of the Ellobiidae (Gastropoda: Panpulmonata) supports independent terrestrial invasions . Mol Phylogenet Evol . 2016 ; 97 : 43 - 54 .
46. Lawton SP , Lim RM , Dukes JP , Kett SM , Cook RT , Walker AJ , et al. Unravelling the riddle of Radix: DNA barcoding for species identification of freshwater snail intermediate hosts of zoonotic digeneans and estimating their interpopulation evolutionary relationships . Infect Genet Evol . 2015 ; 35 : 63 - 74 .
47. White TR , Conrad MM , Tseng R , Balayan S , Golding R , Martins AM , et al. Ten new complete mitochondrial genomes of pulmonates (Mollusca: Gastropoda) and their impact on phylogenetic relationships . BMC Evol Biol . 2011 ; 11 : 295 .