Skip to main content

Extensive gene rearrangements in the mitogenomes of congeneric annelid species and insights on the evolutionary history of the genus Ophryotrocha



Annelids are one the most speciose and ecologically diverse groups of metazoans. Although a significant effort has been recently invested in sequencing genomes of a wide array of metazoans, many orders and families within the phylum Annelida are still represented by a single specimen of a single species. The genus of interstitial annelids Ophryotrocha (Dorvilleidae, Errantia, Annelida) is among these neglected groups, despite its extensive use as model organism in numerous studies on the evolution of life history, physiological and ecological traits. To compensate for the paucity of genomic information in this genus, we here obtained novel complete mitochondrial genomes of six Ophryotrocha species using next generation sequencing. In addition, we investigated the evolution of the reproductive mode in the Ophryotrocha genus using a phylogeny based on two mitochondrial markers (COXI and 16S rDNA) and one nuclear fragment (Histone H3).


Surprisingly, gene order was not conserved among the six Ophryotrocha species investigated, and varied greatly as compared to those found in other annelid species within the class Errantia. The mitogenome phylogeny for the six Ophryotrocha species displayed a separation of gonochoric and hermaphroditic species. However, this separation was not observed in the phylogeny based on the COX1, 16S rDNA, and H3 genes. Parsimony and Bayesian ancestral trait reconstruction indicated that gonochorism was the most parsimonious ancestral reproductive mode in Ophryotrocha spp.


Our results highlight the remarkably high level of gene order variation among congeneric species, even in annelids. This encourages the need for additional mitogenome sequencing of annelid taxa in order to properly understand its mtDNA evolution, high biodiversity and phylogenetic relationships.


Mitochondrial DNA (hereafter mtDNA) has been invaluable in the field of molecular evolution and phylogenetics, and is still a widely used marker today [1, 2]. Maternal inheritance and near absence of recombination have popularised its use in many eukaryotes [3]. Most of metazoan mitochondrial genomes are circular molecules that typically include 13 protein coding genes (PCG), two ribosomal RNA genes, 22 transfer RNA genes and a control region ([4], but see [5, 6]). The mitochondrial gene content is almost invariant among species, but the gene order has been found to vary considerably across Metazoans (such as flatworms, molluscs and tunicates [7]), generating interest in using mitochondrial DNA (mtDNA) gene order for phylogenetic inference [2]. The advent of next generation sequencing has made it easier to obtain mitochondrial genomes even for classical non-model organisms. This enables the detection of gene rearrangements, as well as phylogenetic relationships among and within diverse phyla [8]. Gene order is known to vary extensively within the phyla Mollusca [9], Arthropoda [10] and Annelida [11, 12].

Although over two thousand metazoan genomes have been sequenced to date [13], entire phyla such as Annelida have been widely neglected, with many orders or families often represented by a single specimen of a single species [14]. With more than 17,000 species described, annelids are among the most speciose and ecologically important groups of metazoans [15,16,17]. The extraordinary morphological and ecological diversity of annelids is only comparable to that of crustaceans and molluscs [18]. Furthermore, their exceptional plasticity and adaptability have enabled them to colonise all domains, from marine and freshwater to terrestrial habitats, evolving a wide variety of life history strategies, reproductive modes and feeding habits [19]. Marine annelids have been extensively used as model organisms for the investigation of central questions in ecotoxicology, ecology, physiology, development and evolution, owing to the fact that they are relatively easy to maintain and rear under laboratory conditions (reviewed in [20,21,22]). For example, the bristle worm Capitella teleta Blake, Grassle & Eckelbarger, 2009 (former Capitella sp. I) and the Dumeril’s clam worm Platynereis dumerilii (Audouin & Milne Edwards, 1833) have been primary biological systems for developmental, evolutionary and neurobiological studies [23,24,25], and are both part of a long list of annelids used as indicators for biomonitoring and eco-toxicology tests (reviewed in [26, 27]).

For similar reasons, the interstitial worms in the genus Ophryotrocha Claparède & Mecznikow, 1869 have been widely used to investigate the eco-evolution of functional traits (e.g. [28, 29]) and reproductive strategies (e.g. [30,31,32]), and more recently, to investigate species’ ability to tolerate, respond and adapt to global changes [33,34,35,36,37,38]. Although this group of annelids has been the focus of several descriptive and experimental studies over several decades (Additional file 1), our understanding of their genomics is scarce when compared to other annelids and marine invertebrate genera, thus limiting its potential as model system [22]. As new areas in the sea are explored, new Ophryotrocha species are regularly discovered and described [39,40,41,42,43,44,45,46], and with them the need to clarify the phylogenetic relationships within this genus. For example, recent molecular evidence is suggesting the presence of complexes or lineages in species that were originally considered as independent taxonomic units [42, 45]. Moreover, based on mitochondrial fragments of 16S and cytochrome oxidase I (COXI) genes, and nuclear H3 genes, a clear separation between gonochoric and hermaphroditic species has been proposed by multiple studies, with the hermaphroditism considered as the plesiomophic state for this genus [40, 41]. Providing additional genomic information on Ophryotrocha would certainly increase its usefulness as an emerging interdisciplinary model.

In order to improve our knowledge on the evolutionary history of the genus Ophryotrocha, we: (i) characterised for the first time the complete mitochondrial genomes of six Ophryotrocha species (Ophryotrocha adherens Paavo, Bailey-Brock & Akesson, 2000; Ophryotrocha diadema Åkesson, 1976; Ophryotrocha japonica Paxton & Åkesson, 2010; Ophryotrocha labronica La Greca & Bacci, 1962; Ophryotrocha puerilis Claparède & Mecznikow, 1869 and Ophryotrocha robusta Paxton & Åkesson, 2010), (ii) compared the gene orders of the six species with those described for the main annelids’ taxonomic groups, (iii) used these novel mitogenomes to investigate the phylogenetic relationships among the six corresponding Ophryotrocha species and, finally, (iv) updating the Ophyrotrocha phylogeny to portray how reproductive mode may be linked to evolutionary history in this genus.



The number of pair-end reads obtained varied from 4,328,646 and 12,012,958 million (Table 1). We succeeded in circularising two mitochondrial genomes (O. diadema and O. puerilis). The other genomes were partial, but all 13 PCG were retrieved. The coverage for each mitochondrial genome or fragment varied between 65x (O. labronica) and 256x (O. adherens). The length of mtDNA genomes varied from 14,428 bp (O. robusta) to 15,941 bp (O. puerilis).

Table 1 Sequencing characteristics of the six Ophryotrocha species investigated

Genome organisation and features

Protein coding genes

We retrieved all the 13 PCG found in metazoans. All genes were coded on the plus strand, except for two ribosomal genes and some tRNAs of O. diadema that were found on the minus strand (Fig. 1, Additional files 2, 3, 4, 5, 6, and 7). Most of the genes started with the ATG codon and finished with the TAA stop codon (Table 2). Other alternative start codons were ATT, ATA, TTA and GTG. Ten codons were highly used in all species (Relative Synonymous Codon Usage (RSCU) > 1.25) and most of them were the type of NNU: UAU, GUU, ACU, AAU, UUU, AUA, CUC, AUU, GAU, GCU. Eleven codons were avoided in all species (RSCU < 0.75): UCG, UAC, GUC, ACG, CCG, UUA, CUA, AUC, GGC, UCC and GAC (Table 2). Valine (GUA) was mostly used in O. adherens, O. diadema and O. puerilis. Serine codons (AGU and AGG) were highly biased in O. japonica, O. labronica and O. robusta, while Threonine codon (ACA) was biased in O. labronica and O. puerilis. Regarding cysteine and glutamic acid, codon usage frequency (UGU and GAG) was higher in O. japonica, O. labronica and O. robusta. The AT content was high in the genome of all six Ophryotrocha species and ranged between 64.1% (O. japonica) and 71.5% (O. adherens). In particular, O. adherens, O. diadema and O. robusta harboured the highest AT content. Most of the PCG in each Ophryotrocha species harboured negative AT-skew (Table 3), suggesting a bias in T, with the exception of seven PCG in O. puerilis, which held a positive value of AT-skew. Positive values of AT-skew were observed in the ribosomal genes of five Ophryotrocha species, indicating a bias in A. In addition, negative AT-skew was observed for the ribosomal genes of O. diadema, indicating a bias in T. All PCG in O. adherens, O. diadema and O. puerilis had negative values of GC-skew suggesting a bias in C, whereas positive GC-skew was found for all PCG of O. japonica, O. labronica and O. robusta.

Fig. 1
figure 1

Mitochondrial genome of the six Ophryotrocha species. Colors represent gene complexes. Lightblue = ATP, Pink = cytochrome oxidase, lightpink = Cytochrome b, Purple = NADH, orange = ribosomal subunits, Blue = tRNA, Grey = Non coding region (NCR). The orientation of the arrows represents the orientation of the gene

Table 2 Codon usage in the 13 PCG for each of the six Ophryotrocha species investigated
Table 3 AT content and AT/GC skew of the mitochondrial genomes of the six Ophryotrocha species investigated

The highest variability in gene length was observed for NAD4, whose size ranged from 1239 bp (O. labronica) to 1341 bp (O. diadema and O. adherens). Cytochrome b in O. adherens was the longest compared to the other Ophryotrocha species. Several overlaps and small non-coding regions between tRNA and genes, or between two adjacent genes, were observed in all species and showed sizes that ranged from a few to hundreds of bases (Additional files 2, 3, 4, 5, 6, and 7).

rRNA and tRNA

The length of the small ribosomal unit varied between 736 bp (O. japonica) and 794 bp (O. adherens), while that of the large ribosomal subunit ranged between 1064 bp (O. labronica) and 1140 bp (O. puerilis). The lowest GC content was observed in O adherens for both small and large ribosomal units (24 and 20%, respectively). Ophryotrocha japonica harboured the highest GC content for both ribosomal units: 32% for 12S and 30% for 16S.

In O. labronica and O. diadema, we retrieved 22 tRNA genes based on their secondary leaf structure and their anticodon (see Additional files 8, 9, 10, 11, 12 and 13). In the remaining Ophryotrocha species, Valine and Glutamate were not found, while Aspartate was missing only in O. japonica. In O. labronica, tRNA-Pro was not completely recovered.

Gene rearrangement

The gene order varied widely among Ophryotrocha species (Fig. 2). The higher level of similarity between gene order was found between O. robusta and O puerilis, while the highest dissimilarity was observed between O. japonica and the two species O. adherens and O. diadema (Additional file 14). Ophryotrocha japonica and O. labronica, showed the same gene order, and differed from O. robusta by a single transposition of several genes from COX2 to NAD5. Ophryotrocha robusta and O. puerilis differed by a transposition of NAD2. Ophryotrocha adherens and O. robusta differed by a transposition of ATP8. Ophryotrocha robusta and O. diadema differed by two transpositions of NAD5 and the block of NAD1, NAD2, NAD3 and a reversion of the two ribosomal genes. Ophryotrocha puerilis and O. adherens differed by two transpositions of ATP8 and NAD2. Ophryotrocha adherens and O. diadema differed by three transpositions of NAD5, ATP8 and the block of NAD1, NAD2 and NAD3 genes, and a reversion of the ribosomal genes. Ophryotrocha adherens and O. japonica differed by a tandem duplication random loss of ATP8 and a transposition of the 16S, NAD1, NAD3. The gene order of O. puerilis was the most similar to that found in Eunicidae and Pleistoannelid (number of common interval = 91, Additional file 15), differing only by the position of NAD2 (Fig. 2). Ophryotrocha diadema and the Ampharetidae gene order were similar and differed by the position of NAD1, NAD2, NAD3 and the reversal of the two ribosomal genes. The gene order pattern of Pleistoannelida was similar to the one found in O. adherens and differed only by the transposition of ATP8. Ophryotrocha robusta had a gene order very similar to the Pleistoannelida differing only by a single transposition between NAD4 and NAD4L. The Errantia pattern with O. japonica and O. labronica gene order differed by the transposition of the COX2, ATP8, COX3, NAD6, CYTB, ATP6 and NAD5 regions. Ophryotrocha robusta gene order did not differ from that of Errantia. Magelonidae was closest to the putative bilateral pattern differing by several reversals of NAD3, NAD4L, NAD4, NAD5, NAD6, CYTB, ribosomal genes and NAD1. It also differed from the Lophotrochozoa pattern by several transpositions of COX3, NAD6 and CYTB.

Fig. 2
figure 2

PCG gene order and ribosomal RNAs of the Ophryotrocha species studied and other annelids. Since O. labronica has the same PCG order as O. japonica, only the gene order of the latter was displayed. Each gene is represented by a specific color. The up and down blocks represent the position on the plus or the minus strand of the gene


Mitochondrial phylogeny of the six studied Ophryotrocha species

The concatenation of the amino-acid sequences resulted in a fragment of 4098 residues. After removing the poorly aligned positions, 2242 residues were kept for phylogenetic analyses. The concatenation of the nucleotide sequences resulted in a fragment of 12,297 bp. After removing the poorly aligned positions, 6726 bp were kept for phylogenetic analyses. Bayesian and maximum likelihood amino-acid and nucleotide phylogenies were mostly congruent (Fig. 3). The six Ophryotrocha species clustered into two groups: one including the gonochoric species O. labronica, O. japonica and O. robusta and a second group including the hermaphroditic species O. adherens, O. puerilis and O. diadema. Only the position of O. puerilis within the hermaphrodite clade changes between ML and BI phylogenies.

Fig. 3
figure 3

Phylogenetic trees based on concatenated mitochondrial PCG nucleotides under the K3P model with invariant gamma sites. a: ML tree. The number represents the aLTR obtained by maximum likelihood. Only LTR values greater than 70 are shown. b: BI tree. The number represents the posterior probability (pp). Only posterior probabilities greater than 0.7 are shown

Ophryotrocha genus phylogeny and ancestral state reconstruction of the reproductive mode

Both Bayesian inference (BI) and Maximum likelihood (ML) trees were mostly congruent and two major clades were found, i.e. clade 1 and 2, which displayed high support (aLTR = 100, pp. = 0.99) (Fig. 4 and Fig. 5). In both phylogenies, clade 1 was composed mostly of gonochoric species, except for the hermaphroditic species O. diadema #3 (EF464534). Clade 2 contained species with all three types of reproductive mode: such as O. adherens and species of the O. puerilis group, which are respectively simultaneous and protandrous hermaphroditic species, and some gonochoric species such as Iphitime cuenoti, Iphitime paguri, Ophryotrocha geryonicola and Ophryotrocha globopalpata. A third minor clade, clade 3, containing O. nauarchus, O. globopalpata, O. flabella, and O. longidentata, was identified in both phylogenies. Within clade 3, only the reproductive mode of Ophryotrocha globopalpata is known to be gonochorism. In the BI phylogeny, the species O. diadema EF464534, O. permanni EF464535 and Ophryotrocha sp. Benidorm were grouped together, but not included in a specific clade as observed in the ML phylogeny (aLTR = 83).

Fig. 4
figure 4

Maximum likelihood phylogeny of the genus Ophryotrocha based on the fragment of COXI, 16S and H3 under the GTR model with invariant gamma sites and ancestral state reconstruction. Values represent the aLTR support and the posterior probabilities, respectively. Only values greater than 70 are shown. * represents the species we used in this study. Branch color represents the reproductive mode for all the species and their common ancestors: light blue (protandrous hermaphrodite), dark blue (simultaneous hermaphrodite), pink (gonochoric) and grey for unknown mode

Fig. 5
figure 5

Bayesian phylogeny of the genus Ophryotrocha based on the fragments of COXI, 16S and H3 under the GTR model with invariant gamma sites and ancestral state reconstruction. Values represent the aLTR support and the posterior probabilities, respectively. Only posterior probabilities greater than 0.7 are shown. * represents the species we used in this study. Colors of the branch represent the reproductive mode for all the species and their common ancestors: light blue (protandrous hermaphrodite), dark blue (simultaneous hermaphrodite), pink (gonochoric) and grey for unknown mode

In the phylogeny of the Ophryotrocha genus, we identified 25 gonochoric species, seven protandrous hermaphrodite species and 11 simultaneous hermaphrodite species. All ancestral reproductive mode analyses suggest gonochorism as the reproductive mode for the last common ancestor of Ophryotrocha spp. (Figs. 4, 5, additional file 17). However, for the protandrous group containing all the Ophryotrocha puerilis species and the gonochoric species Ophryotrocha eutrophila, the support of the ancestral state was low (pp < 0.7), suggesting that the reproductive mode for the ancestor is not fully supported.


Our work is the first to report gene rearrangements involving protein coding genes (PCG) among annelid species belonging to the same genus. Also, in contrast to other known annelid mitochondrial genomes, we show that Ophryotrocha diadema possesses the ribosomal genes encoded on the minus strand. Finally, the ancestral state reconstruction of the reproductive mode shows that gonochorism is the plesiomorphic condition in the Ophryotrocha genus, not hermaphroditism as previously hypothesized. Below we discuss in detail these major findings, and their implications for understanding the molecular evolution of the mitochondrial genome within the genus Ophryotrocha.

Genome organisation and features of the six Ophryotrocha species

We provide here the complete mitochondrial genome sequences for six Ophryotrocha species. The features of these mitogenomes did not differ from those reported from other annelid species. Their AT-content, AT and GC-skew, the length of the ribosomal genes and the initiative codons of the PCG are congruent with those reported for other annelids [18, 47,48,49,50,51,52,53,54] and metazoans [55]. GC content and GC-skew have an influence on the usage of codon and proteins [50, 56, 57]. According to their base compositions, codons ending with C or G are avoided in all Ophryotrocha species, which is a trend observed for other annelid species [52].

Surprisingly, the gene order of the mitochondrial genome was not conserved as the positions of PCG and tRNA differed among the six Ophryotrocha species. Mitochondrial gene order tends to be highly conserved among metazoans [58] except in some groups such as molluscs [59, 60] and c [61]. Most frequent changes in gene order involve tRNA, while rearrangements between rRNAs or PCG are rarer [9]. Interestingly, we found rearrangements involving both tRNAs, PCG and ribosomal genes within congeneric species. Gene rearrangements have already been observed in other annelid families. Specimens from the Diurodrilidae [62], Chaetopteridae [12] and Syllidae [11] have shown different gene orders. However, previous findings suggest that these rearrangements are confined to basal lineages and should not be found in more recent lineages such as the Errantia and Sedentaria [12]. Our results do not provide support for this hypothesis, and are instead in line with findings by Ocerguera-Figueroa et al. [63], who reported gene rearrangement in more recent lineages. Several hypotheses have been proposed to explain these rearrangements: the role of tRNAs as mobile elements within the mitochondrial genome [55], intra-mitochondrial recombination [64], or the influence of the oxidative stress [65]. The latter hypothesis would explain why certain clades with species experiencing greater levels of oxidative stress are more prone to rearrangements compared to clades with species experiencing lower levels of oxidative stress [66]. As a consequence, it remains to be further tested whether annelid species with greater degrees of gene rearrangement are also those species that exhibit higher oxidative stress levels.

Ribosomal genes on the minus strand

In addition, we report the presence of ribosomal genes on the minus strand. Usually, all annelid mitochondrial genes are encoded on the same strand except in the tubeworm Owenia fusiformis, the magelonid Magelona mirabilis (Johnston, 1865), and the ragworm Laeonereis culveri (Webster, 1879), which have one or two tRNA located on the minus strand [12, 67]. To explain this, it was hypothesised that, in the last common ancestor of annelids, all the genes were encoded on the same strand by chance wherein a ratchet effect took place, eliminating all the transcriptional elements and preventing the translocation of genes to the other strand [1, 68]. Weigert et al. [12] proposed two hypotheses to explain the translocation of genes between strands in annelids. The first hypothesis implies that in the last common ancestor of annelids, a single strand first encoded all genes and then underwent an inversion-transposition of the tRNAs and their transcriptional elements. The second hypothesis implies that the last common ancestor of annelids still had transcription signals on both strands and that these signals were kept in the basal lineages and lost in the more recent ones. However, O. diadema and Laeonereis culveri [67], two species found in most recent lineages, possess genes encoded on both strands, suggesting either (1) that a more recent ancestor still possessed transcript signals on both strands or (2) that tRNAs were transposed with transcript elements.

Indeed, sequencing strategies and/or annotation methodologies can influence the results in terms of length of the mitogenome sequence recovered or in the gene length [69]. Some studies have documented the presence of dissimilarities in the assembly and annotation of mitochondrial genomes that are not always associated with low coverage [69, 70]. In addition, not all the species are necessarily sensitive to this problem, which complicates the identification of a uniform approach to obtaining reliable mitochondrial genomes [69]. Since gene order is not affected by the methodology used, our results can be reasonably considered not to be an artefact.

Phylogenetic relationships within the genus Ophryotrocha using mitogenomes and traditional marker

Based on our first mitogenome phylogeny of the genus Ophryotrocha, the two main clusters differed in their sexual strategies: one group including only gonochoric species (O. japonica, O. labronica and O. robusta), and the other including only hermaphroditic species (O. adherens, O. diadema and O. puerilis). A segregation based on the reproductive mode has been reported in previous phylogenetic studies [40, 41]. Genomic features differed between these two groups. The hermaphroditic species harboured negative AT-skew and GC-skew values in most of the PCG and ribosomal genes, as reported for other mitogenomes of annelids [50]. Gonochoric species had positive GC-Skew and negative AT-skew for each PCG. Codon usage was also different between the two-sexual modes: to code glutamic acid, the gonochoric species use the GAG codon, whereas the hermaphrodites use the GAA. Similar differences in codon usage remain to be confirmed in additional hermaphroditic and gonochoric Ophryotrocha species.

Based on the phylogeny of the Ophryotrocha genus we obtained, two main clades were identified: clade 1 that contains most of gonochoric species and clade 2 mostly composed of hermaphrodite species, but also including species with gonochoric strategies. Lending a closer look within the clade, reveal that a separation between reproductive modes was not as clear as it was reported in previous phylogenies including only a few species [40]. In addition, several groups were identified as O. labronica, this suggesting the existence of issues related to taxonomic identification. In contrast to Dahlgren et al. [40] and Heggoy et al. [41], we observed another clade within the hermaphrodite group. The separation of this clade was not linked to the reproductive mode of the species, although most of them are lacking information on their reproductive mode. Wilkund et al. [42] suggested a separation according to their habitat. As we do not possess this information for all species, it was not possible to comprehensively test for this hypothesis.

In addition, it is worth noticing the presence of different lineages with the same name in our phylogenies. This observation confirm the presence of cryptic species within groups that were originally described as independent species, as hypothesized for O. labronica lineages and for the O. puerilis complex, this confirming what reported in recent studies (e.g. [42, 43, 45]). In some cases, some gene sequences with different names corresponded to the same species, as for O. obscura and O. sanya (previously O. vellae by Paxton et al. [71]). In certain cases, the wrong taxonomic identification of the species may have occurred, prompting the incorrect association to a given reproductive mode, such as for sequences labeled as O. diadema #3 and O. permanni #2. Finally, we confirmed that Iphitime species do cluster in the Ophryotrocha genus [39, 41, 43, 45, 72, 73]. However, in our phylogeny, the Iphitime species are closely related to O. adherens, O. puerilis and O. socialis, whereas in the phylogeny from Heggoy et al. [40] the Iphitime specimen was closer to O. gracilis and O. hartmanni. Contrary to Taboada et al. [45], we found O. clava closely related to O. jiaolongi and the Iphitime species. This result was also observed by Zhang et al. [39].

As reported in previous studies, three different reproductive modes are present in the Ophryotrocha genus, i.e. gonochorism, simultaneous hermaphroditism and proterandrous hermaphroditism: the latter found in only one species, O. puerilis [31]. Contrary to other phyla, such as Mollusca and Arthrpoda, the Annelida phylum is known for its high diversity in reproductive modes, even within families [74]. The relative simplicity of their reproductive system and the relaxed morpho-physiological constraints to the evolution of alternative reproductive strategies seem to have favoured this remarkable variation of reproductive modes even among congeneric species [75]. Our results suggest that transitions between reproductive modes, i.e. gonochorism to hermaphroditism and vice versa, seem to have occurred multiple times along the evolutionary history of this genus. Indeed, both gonochoric and hermaphroditic species of the Ophryotrocha genus show some degree of sexual lability in the population that can potentially favour the expression of alternative reproductive strategies [76]. In particular, the presence of sexual phenotypes in the gonochoric species (i.e. pure male, male with a few oocytes, pure female, and female with a few sperm) is considered a vestigial trait of an ancestral hermaphroditic state [77]. Interestingly, our study shows that the ancestral reproductive strategy of this genus is most likely gonochorism, not hermaphroditism as previously reported [40, 41, 78]. Moreover, while clade 1 is composed of only gonochoristic species, clade 2 contained all modes of reproduction. These findings are in line with a recent study suggesting that transitions from gonochorism to hermaphroditism are more common than the reverse in many animal taxa. Factors promoting hermaphroditism in gonochoristic animals, such as low densities and inbreeding depression, are in fact more widespread or create stronger selection pressures than the conditions promoting gonochorism in hermaphroditic animals, such as high density and reproductive assurance [79].

The evolution of hermaphroditic strategies in this genus has been generally explained as an adaptation to conditions of low density, stabilized by poor mate search efficiency and high costs of searching [32]. In the Mediterranean Sea, the gonochoristic species O. labronica and O. japonica are indeed present in greater densities compared to the other hermaphrodite species [80]. In addition, Prevedelli et al. [31] demonstrated that gonochoristic and hermaphroditic species differed in a number of life-history traits. These differences confer to the former higher fitness and demographic advantages over the latter in conditions of high mate-search efficiency. However, the lack of information on environmental population densities and on life history for the majority of the Ophryotrocha species comprising our phylogeny prevents this hypothesis from being formally tested.

Finally, we cannot completely discard the idea that the ancestral reproductive mode could have been somewhat an “intermediate” one between gonochorism and hermaphroditism, given the documented wide variety of sexual phenotypes found in some species of these genus expressing these two extant forms of reproduction [76, 77]. Further research is therefore needed to better understand the phylogenetic relationship among Ophryotrocha species as emerging models for the investigation of evolutionary global change biology [37, 38].


The descriptions of unique gene rearrangements within the Ophryotrocha genus are remarkable as they suggest that mitochondrial genomes in this taxonomic group are highly dynamic, signalling that gene rearrangements can occur more rapidly than previously thought. The within-genus PCG rearrangement refutes the idea that the gene order is conserved among the Errantia, although further studies are required to determine the mechanisms involved. The use of next generation sequencing techniques on Ophryotrocha has revealed the significant potential of these species as model organisms for studying evolutionary history within this genus. Moreover, this study displays the remarkable level of genetic diversification in annelids found even among closely related species. It highlights the need to increase the taxonomic representation in future phylogenetic studies for a more accurate understanding of this phylum’s diversity. Finally, developing an in-depth genomic understanding of the Ophryotrocha genus will help further the investigation of both evolution of life-history traits and the emerging field of evolutionary global change biology [37, 38, 81, 82].


Specimen collection, genomic extractions and sequencing

Specimens of the six Ophryotrocha species originated from individuals collected in 2008 (O. robusta, O. diadema, O. adherens) in the harbour of Porto Empedocle (Italy; 37°17′4″N, 13°31′3″E) and in 2015 (O. labronica, O. japonica, O. puerilis) in the harbour of La Spezia (Italy, 44° 6′ 24″ N, 9° 49′ 45″ E). Species identification was performed using morphological and reproductive traits as described by Simonini et al. [80]. For each species, three samples (each from a single breeding pair) of 30 individuals were used for mtDNA extraction. Genomic extractions were performed using QIAGEN’s DNeasy Blood and Tissue Kit with RNAse according to the manufacturer’s protocol. Quantification of DNA was done with Quant-iT™ PicoGreen® dsDNA Assay Kit (Invitrogen™, Canada following the manufacter’s protocol. Genomic DNA (500 ng) was mechanically fragmented for 40 s. using a covaris M220 (Covaris, Woburn MA, USA) with default settings. Fragmented DNA was transferred to PCR tubes and library synthesis was performed with the NEB Next Ultra II (New England Biolabs) according to manufacturer’s instructions. TruSeq HT adapters (Illumina, SanDiego, CA, USA) were used to barcode the samples. The libraries were quantified and pooled using an equimolar ratio and sequenced on an Illumina MiSeq 300 base-pair (bp) paired-end run (600 cycle, v3 kit) at the Plateforme d’Analyses Génomiques of the Institut de Biologie Intégrative et des Systèmes (Laval University, Quebec, Canada).

Assembly and annotations

The quality of the sequencing was assessed with Fastqc [83], and adapters were removed with Trimmomatic [84] available in [85]. Default parameters were used to retrieve mitochondrial genomes using the perl script Novoplasty2.7.2 [86]. Briefly, based on a seed-and-extend algorithm, mitochondrial genome is retrieved from the whole genome sequencing data, using a related or distant single seed sequence [87]. For each species, we used a fragment of COX1 as seed: JQ310756.1 for O. adherens, JQ310758.1 for O. diadema, EF46454.1 for O. japonica, KF305814.1 for O. labronica, EF 464544.1 for O. puerilis, EF464547.1 for O. robusta. As no mitochondrial genome from close relatives of each species was available, no references were used for the assembly. Annotation was performed with MITOS2 Web Server [87], verified with ORF finder [88] and ARWEN [89], and visualized using GeSeq [90]. Determination of the A + T content of protein-coding genes, tRNA genes, rRNA genes and the RSCU was performed with DAMBE 6 [91]. All the mitochondrion genomes have been deposited in Genbank under the following accession number: MT737360 (O. robusta), MT737361 (O. labronica), MT737362 (O. japonica), MT737363 (O. adherence), MT737364 (O. diadema), MT737365 (O. puerilis).

Gene rearrangement

CREx [92] was used to examine the possible scenarios of rearrangement between pairs of complete genomes. Briefly, this method compares two genomes and determines the most parsimonious scenario that has led to the observed rearrangement accounting for duplication reversals, transpositions or other events. As the tRNA order is not conserved among Annelida species [12], we kept only the PCG order to infer the possible scenarios of rearrangement. We examined all the rearrangement scenarios among Ophryotrocha species and other gene orders known for annelids obtained from Lavrov and Lang [93], Mwinyi et al. [54] and Weigert et al. [12]. Results were visualized with the R-package genoPlotR [94].

Phylogeny reconstruction

Phylogenies were assessed using maximum likelihood and Bayesian inference. First, PCG from whole mitochondrial sequences were used to reconstruct phylogenies for the six Ophryotrocha species considered. Each PCG was separately aligned and then concatenated using seaview4 [95]. Poorly aligned regions were removed using Gblocks [96] with default parameters. This step was performed on both nucleotide and amino-acid sequences of PCG. Two species were used as outgroup in the phylogenies: Marphysa sanguinea (NC_023124.1) and O. fusiformis (NC_028712.1). Secondly, in order to investigate the relationship among Ophryotrocha species, we built a phylogeny based on two mitochondrial fragments, the COI gene and the 16S, and a nuclear fragment from Histone 3 (H3). First, we retrieved the Histone 3 fragment for each Ophryotrocha species used in this study using a BLASTn search available in Galaxy [85] with a H3 fragment from a close relative available in Genbank. Sequences of H3 have been deposited in Genbank (MT733538-MT733543). In order to determine if the three genes could be concatenated for phylogenetic analyses, we tested the congruence among the three distance matrices genes using the congruence among distance matrices approach (CADM) developed by [97, 98] available in ape R-package [99]. Briefly, the CADM tests for the presence of incongruency among all the distance matrices. The significance of the test was performed with a 1000 permutations. All three distance matrices were obtained in MEGA5 [100]. No incongruence was detected among the matrices (W = 0.65, p = 0.001) and all the sequences were concatenated. A total of 71 sequences representing 46 species from Genbank (listed in Additional file 16) and the six-species used in this study were aligned for each gene with MAFFT version 7 [101] available at . All the COXI aligned sequences were trimmed to 400 bp length and all the H3 sequences were trimmed to 298 bp length. All three genes were concatenated. We used a combination of close and distant species as outgroups: Eunice pennata, Eunice norvegica, Protodorvillea gracilis, Dorvillea erucaeformis. For each phylogeny, the most adapted evolutionary model was determined based on Bayesian information criterion as implemented in W-IQ TREE [102]. Maximum likelihood (ML) trees were subsequently generated using IQ-Tree [103] and branch support estimated using the Shimodaira-Hasegawa approximate likelihood ratio test (aLTR), as described in Anisimova and Gascuel [104]. For each dataset, Bayesian inferences (BI) were performed on two runs until convergence was reached under the appropriate evolutionary model. Tree sampling was done every 1000 generations with a burn-in of 25%. Bayesian reconstructions were performed using MrBayes 3.2.6 [105] available on the CIPRES gateway [106]. The posterior probabilities (pp) were obtained for the 50% majority-rule consensus tree. Strong support of branches was considered when pp. ≥ 0.95, whereas intermediate support was defined with pp. values between 0.85 and 0.94.

Ancestral state reconstruction of reproductive trait in Ophryotrocha

Information about the reproductive mode of each species (gonochoric, simultaneous hermaphrodite and protandrous hermaphrodite) was retrieved from the literature (Additional file 16) and all species coded accordingly. When no information on the reproductive mode was found, species were marked as unknown. Parsimony ancestral trait reconstruction was performed on both ML and BI phylogeny in order to find the ancestral states that minimize the number of changes according to the phylogeny. Analysis was conducted in Mesquite v3.4 [107].

In addition to the maximum parsimony reconstruction, we also used a method based on Bayesian MCMC sampling methods to reconstruct the ancestral reproductive mode of the genus Ophryotrocha as implemented in BEAST2 [108]. In particular, we used the Bayesian phylogeny previously obtained and estimate the posterior probability of the state for each ancestor for each node of the tree. We conducted 10,000,000 iterations and sampled parameters every 10,000 generations. The posterior distribution was first verified with Tracer1.7 [109]. We discarded the first 25% samples of states obtained as burn-in as implemented in TreeAnnotaor from the BEAST2 package.

Availability of data and materials

All mitochondrial sequences and the Histone 3 fragment were deposited under the accession number MT737360-MT737365 and MT733538-MT733543 in the Genbank database. Alignments and script were available as Supplementary file. The accession number of all sequences of COI, 16S and H3 used in this study are listed in the additional file 16 of the Supplementary Material and are available in Genbank database. Mitogenome of Marphysa sanguinea (NC_023124.1) and O. fusiformis (NC_028712.1), used for phylogenies reconstruction, are available in Genbank. The following accession numbers: JQ310756.1, JQ310758.1, KF305814.1, EF 464544.1 and EF464547.1 are available in Genbank.





Protein Coding Gene


Cytochrome oxidase I


Relative Synonymous Codon Usage


Bayesian inference


Maximum likelihood


Approximate likelihood ratio test


Posterior probabilities


  1. Boore JL. Animal mitochondrial genomes. Nucleic Acids Res. 1999;27(8):1767–80.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  2. Snel B, Huynen MA, Dutilh BE. Genome trees and the nature of genome evolution. Annu Rev Microbiol. 2005;59:191–209.

    Article  CAS  PubMed  Google Scholar 

  3. Avise JC. Phylogeography: the history and formation of species. Cambridge: Harvard University press; 2000.

  4. Wolstenholme DR. Animal mitochondrial DNA: structure and evolution. Int Rev Cytol. 1992;141:173–216.

    Article  CAS  PubMed  Google Scholar 

  5. Stampar SN, Broe MB, Macrander J, Reitzel AM, Brugler MR, Daly M. Linear mitochondrial genome in Anthozoa (Cnidaria): a case study in Ceriantharia. Sci Rep. 2019;9(1):1–12.

    Article  CAS  Google Scholar 

  6. Lavrov DV, Adamski M, Chevaldonne P, Adamska M. Extensive mitochondrial mRNA editing and unusual mitochondrial genome organization in calcaronean sponges. Curr Biol. 2016;26(1):86–92.

    Article  CAS  PubMed  Google Scholar 

  7. Podsiadlowski L, Mwinyi A, Lesný P, Bartolomaeus T. Mitochondrial gene order in Metazoa–theme and variations. In: Wägele JW, Bartolomaeus T, editors. Deep Metazoan Phylogeny: the backbone of the tree of life. Berlin: Walter De Gruyter GmbH; 2014. p. 459–72.

    Google Scholar 

  8. Richter S, Schwarz F, Hering L, Böggemann M, Bleidorn C. The utility of genome skimming for phylogenomic analyses as demonstrated for glycerid relationships (Annelida, Glyceridae). Genome Biol Evol. 2015;7(12):3443–62.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  9. Rawlings TA, Collins TM, Bieler R. A major mitochondrial gene rearrangement among closely related species. Mol Biol Evol. 2001;18(8):1604–9.

    Article  CAS  PubMed  Google Scholar 

  10. Hickerson MJ, Cunningham CW. Dramatic mitochondrial gene rearrangements in the hermit crab Pagurus longicarpus (Crustacea, Anomura). Mol Biol Evol. 2000;17(4):639–44.

    Article  CAS  PubMed  Google Scholar 

  11. Aguado MT, Richter S, Sontowski R, Golombek A, Struck TH, Bleidorn C. Syllidae mitochondrial gene order is unusually variable for Annelida. Gene. 2016;594(1):89–96.

    Article  CAS  PubMed  Google Scholar 

  12. Weigert A, Golombek A, Gerth M, Schwarz F, Struck TH, Bleidorn C. Evolution of mitochondrial gene order in Annelida. Mol Phylogenet Evol. 2016;94:196–206.

    Article  CAS  PubMed  Google Scholar 

  13. Organelle Genome Resource. Accessed 20 May 2019.

  14. Zhang Y, Sun J, Rouse GW, Wiklund H, Pleijel F, Watanabe HK, Chen C, Qian PY, Qiu JW. Phylogeny, evolution and mitochondrial gene order rearrangement in scale worms (Aphroditiformia, Annelida). Mol Phylogenet Evol. 2018;125:220–31.

    Article  CAS  PubMed  Google Scholar 

  15. McHugh D. Molecular phylogeny of the Annelida. Can J Zool. 2000;78(11):1873–84.

    Article  CAS  Google Scholar 

  16. Brusca RC, Brusca GJ. Invertebrates. Basingstoke; 2003.

    Google Scholar 

  17. Zhang ZQ. Animal biodiversity: an update of classification and diversity in 2013. Zootaxa. 2013;3703(1):5–11.

    Article  Google Scholar 

  18. Schmidt-Rhaesa A, Harzsch S, Purschke G. Structure and evolution of invertebrate nervous systems. Oxford and New York: Oxford University Press; 2015.

  19. Giangrande A. Polychaete reproductive patterns, life cycles and life histories: an overview. Oceanogr Mar Biol Ann Rev. 1997;35:323–86.

    Google Scholar 

  20. Weigert A, Bleidorn C. Current status of annelid phylogeny. Org Divers Evol. 2016;16(2):345–62.

    Article  Google Scholar 

  21. Shain DH. Annelids in modern biology. Hoboken: Wiley-Blackwell; 2009.

  22. Halanych KM, Borda E. Developing models for Lophotrochozoan and Annelid biology. In: Shain DH, editor. Annelids in Modern Biology. Hoboken: Wiley-Blackwell; 2009. p. 1-12.

  23. Fischer AH, Henrich T, Arendt D. The normal development of Platynereis dumerilii (Nereididae, Annelida). Front Zool. 2010;7(1):31.

    Article  PubMed  PubMed Central  Google Scholar 

  24. Arendt D, Technau U, Wittbrodt J. Evolution of the bilaterian larval foregut. Nature. 2001;4409:81–5.

    Article  CAS  Google Scholar 

  25. Simakov O, Marletaz F, Cho SJ, Edsinger-Gonzales E, Havlak P, Hellsten U, Kuo DH, Larsson T, Lv J, Arendt D, et al. Insights into bilaterian evolution from three spiralian genomes. Nature. 2013;493:526–31.

    Article  CAS  PubMed  Google Scholar 

  26. Libralato G, Minetto D, Lofrano G, Guida M, Carotenuto M, Aliberti F, Barbara C, Notarnicola M. Toxicity assessment within the application of in situ contaminated sediment remediation technologies: a review. Sci Total Environ. 2018;621:85–94.

    Article  CAS  PubMed  Google Scholar 

  27. Diaz V, Reish D. Polychaetes in Environmental Studies. In: Shain DH, editor. Annelids in Modern Biology. Hoboken: Wiley-Blackwell; 2009. Chapter 11..

  28. Thornhill DJ, Dahlgren TG, Halanych K. The evolution and ecology of Ophryotrocha (Dorvilleidae, Eunicida). In: Shain DH, editor. Annelids as Model Systems in the Biological Sciences. Hoboken: Wiley-Blackwell; 2009. p. 242-256.

  29. Massamba-N’Siala G, Calosi P, Bilton DT, Prevedelli D, Simonini R. Life-history and thermal tolerance traits display different thermal plasticities and relationships with temperature in the marine polychaete Ophryotrocha labronica La Greca and Bacci (Dorvilleidae). J Exp Mar Biol Ecol. 2012;438:109–17.

    Article  Google Scholar 

  30. Massamba-N’Siala G, Prevedelli D, Simonini R. Trans-generational plasticity in physiological thermal tolerance is modulated by maternal pre-reproductive environment in the polychaete Ophryotrocha labronica. J Exp Biol. 2014;217:2004–12.

    Article  PubMed  Google Scholar 

  31. Prevedelli D, Massamba N’Siala G, Simonini R. Gonochorism vs hermaphroditism: relationship between life history and fitness in three species of Ophryotrocha (Polychaeta: Dorvilleidae) with different form of sexuality. J An Ecol. 2006;75:203–12.

    Article  Google Scholar 

  32. Schleicherova D, Lorenzi MC, Sella G, Michiels NK. Gender expression and group size: a test in a hermaphroditic and a gonochoric congeneric species of Ophryotrocha (Polychaeta). J Exp Biol. 2010;213(9):1586–90.

    Article  CAS  PubMed  Google Scholar 

  33. Lorenzi MC, Meconcelli S, Sella G. Social recognition in Annelids and the evolution of social recognition and cognitive abilities by sexual selection. In: Social Recognition in Invertebrates. Cham: Springer; 2015. p. 1–15.

    Google Scholar 

  34. Calosi P, De Wit P, Thor P, Dupont S. Will life find a way? Evolution of marine species under global change. Evol Appl. 2016;9:1035–42.

    Article  PubMed  PubMed Central  Google Scholar 

  35. Chakravarti LJ, Jarrold MD, Gibbin EM, Christen F, Massamba-N’Siala G, Blier PU, Calosi P. Can trans-generational experiments be used to enhance species resilience to ocean warming and acidification? Evol Appl. 2016;9(9):1133–46.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  36. Gibbin EM, Massamba-N’Siala G, Chakravarti LJ, Jarrold MD, Calosi P. The evolution of phenotypic plasticity under global change. Sci Rep. 2017;7(1):17253.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  37. Jarrold MD, Chakravarti LJ, Gibbin EM, Christen F, Massamba-N'Siala G, Blier PU, Calosi P. Life-history trade-offs and limitations associated with phenotypic adaptation under future ocean warming and elevated salinity. Philos Trans Royal Soc B. 2019;374(1768):20180428.

    Article  Google Scholar 

  38. Thibault C, Massamba-N’Siala G, Noisette F, Vermandele F, Babin M, Calosi P. Within- and trans-generational responses to ocean warming and acidification are highly divergent in two congeneric species of marine annelids. Mar Biol. 2020;167(4):1-17.

  39. Zhang DS, Zhou YD, Wang CS, Rouse GW. A new species of Ophryotrocha (Annelida, Eunicida, Dorvilleidae) from hydrothermal vents on the southwest Indian ridge. ZooKeys. 2017;687:1.

    Article  Google Scholar 

  40. Dahlgren TG, Åkesson B, Schander C, Halanych KM, Sundberg P. Molecular phylogeny of the model annelid Ophryotrocha. Biol Bull. 2001;201(2):193–203.

    Article  CAS  PubMed  Google Scholar 

  41. Heggøy KK, Schander C, Åkesson B. The phylogeny of the annelid genus Ophryotrocha (Dorvilleidae). Mar Biol Res. 2007;3(6):412–20.

    Article  Google Scholar 

  42. Wiklund H, Glover AG, Dahlgren TG. Three new species of Ophryotrocha (Annelida: Dorvilleidae) from a whale-fall in the north-East Atlantic. Zootaxa. 2009;2228:43–56.

    Article  Google Scholar 

  43. Wiklund H, Altamira IV, Glover AG, Smith CR, Baco AR, Dahlgren TG. Systematics and biodiversity of Ophryotrocha (Annelida, Dorvilleidae) with descriptions of six new species from deep-sea whale-fall and wood-fall habitats in the north-East Pacific. Syst Biodivers. 2012;10(2):243–59.

    Article  Google Scholar 

  44. Paxton H, Wiklund H, Alexander F, Taboada S. Is the Antarctic Ophryotrocha orensanzi (Annelida: Dorvilleidae) a circumpolar non-specialized opportunist? Syst Biodivers. 2017;15(2):105–14.

    Article  Google Scholar 

  45. Taboada S, Wiklund H, Glover AG, Dahlgren TG, Cristobo J, Avila C. Two new Antarctic Ophryotrocha (Annelida: Dorvilleidae) described from shallow-water whale bones. Polar Biol. 2013;36(7):1031–45.

    Article  Google Scholar 

  46. Salvo F, Dufour SC, Hamoutene D, Parrish CC. Lipid classes and fatty acids in Ophryotrocha cyclops, a Dorvilleid from Newfoundland aquaculture sites. PLoS One. 2015;10(8):e0136772.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  47. Okimoto R, Macfarlane JL, Wolstenholme DR. Evidence for the frequent use of TTG as the translation initiation codon of mitochondrial protein genes in the nematodes, Ascaris suum and Caenorhabditis elegans. Nucleic Acids Res. 1990;18(20):6113–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  48. Boore JL, Brown WM. Mitochondrial genomes of Galathealinum, Helobdella, and Platynereis: sequence and gene arrangement comparisons indicate that Pogonophora is not a phylum and Annelida and Arthropoda are not sister taxa. Mol Biol Evol. 2000;17(1):87–106.

    Article  CAS  PubMed  Google Scholar 

  49. Bleidorn C, Vogt L, Bartolomaeus T. New insights into polychaete phylogeny (Annelida) inferred from 18S rDNA sequences. Mol Phylogen Evol. 2003;29(2):279–88.

    Article  CAS  Google Scholar 

  50. Zhong M, Struck TH, Halanych KM. Phylogenetic information from three mitochondrial genomes of Terebelliformia (Annelida) worms and duplication of the methionine tRNA. Gene. 2008;416:11–21.

    Article  CAS  PubMed  Google Scholar 

  51. Shen X, Ma X, Ren J, Zhao F. A close phylogenetic relationship between Sipuncula and Annelida evidenced from the complete mitochondrial genome sequence of Phascolosoma esculenta. BMC Genomics. 2009;10(1):136.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  52. Kim H, Kim HJ, Lee YH. The complete mitochondrial genome of the marine polychaete: Hediste diadroma (Phyllodocida, Nereididae). Mitochondrial DNA Part B. 2016;1(1):822–3.

    Article  PubMed  PubMed Central  Google Scholar 

  53. Xu Y, Nie J, Hou J, Xiao L, Lv P. Complete mitochondrial genome of Hirudo nipponia (Annelida, Hirudinea). Mitochondrial DNA Part A. 2016;27(1):257–8.

    Article  CAS  Google Scholar 

  54. Mwinyi A, Meyer A, Bleidorn C, Lieb B, Bartolomaeus T, Podsiadlowski L. Mitochondrial genome sequence and gene order of Sipunculus nudus give additional support for an inclusion of Sipuncula into Annelida. BMC Genomics. 2009;10(1):27.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  55. Saccone C, De Giorgi C, Gissi C, Pesole G, Reyes A. Evolutionary genomics in Metazoa: the mitochondrial DNA as a model system. Gene. 1999;238(1):195–209.

    Article  CAS  PubMed  Google Scholar 

  56. Knight RD, Freeland SJ, Landweber LF. A simple model based on mutation and selection explains trends in codon and amino-acid usage and GC composition within and across genomes. Genome Biol. 2001;2(4):research0010–1.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  57. Min XJ, Hickey DA. DNA asymmetric strand bias affects the amino acid composition of mitochondrial proteins. DNA Res. 2007;14(5):201–6.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  58. Boore JL, Brown WM. Big trees from little genomes: mitochondrial gene order as a phylogenetic tool. Curr Opin Genet Dev. 1998;8(6):668–74.

    Article  CAS  PubMed  Google Scholar 

  59. Bernt M, Braband A, Middendorf M, Misof B, Rota-Stabelli O, Stadler PF. Bioinformatics methods for the comparative analysis of metazoan mitochondrial genome sequences. Mol Phylogenet Evol. 2013;69(2):320–7.

    Article  CAS  PubMed  Google Scholar 

  60. Osca D, Irisarri I, Todt C, Grande C, Zardoya R. The complete mitochondrial genome of Scutopus ventrolineatus (Mollusca: Chaetodermomorpha) supports the Aculifera hypothesis. BMC Evolut Biol. 2014;14(1):197.

    Google Scholar 

  61. Black WC 4th, Roehrdanz RL. Mitochondrial gene order is not conserved in arthropods: Prostriate and metastriate tick mitochondrial genomes. Mol Biol Evol. 1998;15(12):1772–85.

    Article  CAS  PubMed  Google Scholar 

  62. Golombek A, Tobergte S, Nesnidal MP, Purschke G, Struck TH. Mitochondrial genomes to the rescue – Diurodrilidae in the myzostomid trap. Mol Phylogenet Evol. 2013;68:312–26.

    Article  CAS  PubMed  Google Scholar 

  63. Oceguera-Figueroa A, Manzano-Marin A, Kvist S, Moya A, Siddall ME, Latorre A. Comparative mitogenomics of leeches (Annelida: Clitellata): Genome conservation and Placobdella-specific trnD gene duplication. PLoS One. 2016;11(5):e0155441.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  64. Luo YJ, Satoh N, Endo K. Mitochondrial gene order variation in the brachiopod Lingula anatina and its implications for mitochondrial evolution in lophotrochozoans. Mar Genomics. 2015;24:31–40.

    Article  CAS  PubMed  Google Scholar 

  65. Kajander OA, Rovio AT, Majamaa K, Poulton J, Spelbrink JN, Holt IJ, Karhunen PJ, Jacobs HT. Human mtDNA sublimons resemble rearranged mitochondrial genomes found in pathological states. Hum Mol Gen. 2000;9(19):2821–35.

    Article  CAS  PubMed  Google Scholar 

  66. Dowton M, Campbell NJ. Intramitochondrial recombination–is it why some mitochondrial genes sleep around? Trends Ecol Evol. 2001;16(6):269–71.

    Article  CAS  PubMed  Google Scholar 

  67. Seixas VC, Paiva PC, de Moraes Russo CA. Complete mitochondrial genomes are not necessarily more informative than individual mitochondrial genes to recover a well-established annelid phylogeny. Gene ReP. 2016;5:10–7.

    Article  Google Scholar 

  68. Vallès Y, Boore JL. Lophotrochozoan mitochondrial genomes. Integr Comp Biol. 2006;46(4):544–57.

    Article  PubMed  Google Scholar 

  69. Timbó RV, Togawa RC, Costa MM, Andow DA, Paula DP. Mitogenome sequence accuracy using different elucidation methods. PLoS One. 2017;12(6):e0179971.

    Article  CAS  Google Scholar 

  70. Sun X, Wang Y, Chen P, Wang H, Lu L, Ye Z, Wu Y, Bu W, Xie Q. Biased heteroplasmy within the mitogenomic sequences of Gigantometra gigas revealed by sanger and high-throughput methods. Syst Zool. 2018;43(4):356–86.

    Google Scholar 

  71. Paxton H, Åkesson B. The Ophryotrocha labronica group (Annelida: Dorvilleidae) - with the description of seven new species. Zootaxa. 2010;2713(1):1–24.

    Article  Google Scholar 

  72. Orensanz JM. The eunicemorph polychaete annelids from Antarctic and Subantarctic seas, biology of the Antarctic seas XXI. Antarctic Res Series. 1990;52:1–183.

    Article  Google Scholar 

  73. Høisæter T, Samuelsen TJ. Taxonomic and biological notes on a species of Iphitime (Polychaeta, Eunicida) associated with Pagurus prideauxfrom Western Norway. Mar Biol Res. 2006;2:333–54.

    Article  Google Scholar 

  74. Wilson WH. Sexual reproductive modes in polychaetes: classification and diversity. Bull Mar Sci. 1991;48(2):500–16.

    Google Scholar 

  75. Da] es RP. Annelids. Oxford: Hutchinson Library; 1967. p. 200.

    Google Scholar 

  76. Meconcelli S, Lorenzi MC, Sella G. Labile sex expression and the evolution of dioecy in Ophryotrocha Polychaete worms. EvolBiol. 2015;42(1):42–53.

    Google Scholar 

  77. Lorenzi MC, Sella G. In between breeding systems: neither dioecy nor 512 androdioecy explains sexual polymorphism in functionally dioecious worms. 513. Integr Comp Biol. 2013;53(4):689–700.

    Article  PubMed  Google Scholar 

  78. Pleijel F, Eide R. The phylogeny of Ophryotrocha (Dorvilleidae: Eunicida: Polychaeta). J Nat Hist. 1996;30:647–59.

    Article  Google Scholar 

  79. Sasson DA, Ryan JF. A reconstruction of sexual modes throughout animal evolution. BMC Evol Biol. 2017;17(1):242.

    Article  PubMed  PubMed Central  Google Scholar 

  80. Simonini R, Massamba-N’Siala G, Grandi V, Prevedelli D. Distribution of the genus Ophryotrocha (Polychaeta) in Italy: new reports and comments on the biogeography of Mediterranean species. Vie Milieu. 2009;59(1):79–88.

    Google Scholar 

  81. Munday PL, Warner RR, Monro K, Pandolfi JM, Marshall DJ. Predicting evolutionary responses to climate change in the sea. Ecol Lett. 2013;16(12):1488–500.

    Article  PubMed  Google Scholar 

  82. Reusch TBH. Climate change in the oceans: evolutionary versus phenotypically plastic responses of marine animals and plants. Evol Appl. 2014;7:104–22.

    Article  PubMed  Google Scholar 

  83. Andrews S. FastQCA, quality control tool for high throughput sequence data. 2010. Accessed May 2019.

    Google Scholar 

  84. Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics. 2014;30(15):2114–20.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  85. Afgan E, Baker D, Batut B, Van Den Beek M, Bouvier D, Čech M, et al. The galaxy platform for accessible, reproducible and collaborative biomedical analyses: 2018 update. Nucleic Acids Res. 2018;46(W1):W537–44.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  86. Dierckxsens N, Mardulyn P, Smits G. NOVOPlasty: de novo assembly of organelle genomes from whole genome data. Nucleic Acids Res. 2016;45(4):e18.

    PubMed Central  Google Scholar 

  87. Bernt M, Donath A, Jühling F, Externbrink F, Florentz C, Fritzsch G, Pütze J, Middendorf M, Stadler PF. MITOS: improved de novo metazoan mitochondrial genome annotation. Mol Phylogenet Evol. 2013;69(2):313–9

    Article  PubMed  Google Scholar 

  88. Rombel IT, Sykes KF, Rayner S, Johnston SA. ORF-FINDER: a vector for high-throughput gene identification. Gene. 2002;282(1–2):33–41.

    Article  CAS  PubMed  Google Scholar 

  89. Laslett D, Canbäck B. ARWEN: a program to detect tRNA genes in metazoan mitochondrial nucleotide sequences. Bioinformatics. 2007;24(2):172–5.

    Article  PubMed  CAS  Google Scholar 

  90. Tillich M, Lehwark P, Pellizzer T, Ulbricht-Jones ES, Fischer A, Bock R, Greiner S. GeSeq–versatile and accurate annotation of organelle genomes. Nucleic Acids Res. 2017;45(W1):W6–11.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  91. Xia X. DAMBE6: new tools for microbial genomics, phylogenetics, and molecular evolution. J Hered. 2017;108(4):431–7.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  92. Bernt M, Merkle D, Ramsch K, Fritzsch G, Perseke M, Bernhard D, Schlegel M, Stadler PF, Middendorf M. CREx: inferring genomic rearrangements based on common intervals. Bioinformatics. 2007;23(21):2957–8

    Article  CAS  PubMed  Google Scholar 

  93. Lavrov DV, Lang BF. Poriferan mtDNA and animal phylogeny based on mitochondrial gene arrangements. Syst Biol. 2005;54(4):651–9.

    Article  PubMed  Google Scholar 

  94. Guy L, Roat Kultima J, Andersson SGE. genoPlotR: comparative gene and genome visualization in R. Bioinformatics. 2010;26(18):2334–5.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  95. Gouy M, Guindon S, Gascuel O. SeaView version 4: a multiplatform graphical user interface for sequence alignment and phylogenetic tree building. Mol Biol Evol. 2009;27(2):221–4.

    Article  PubMed  CAS  Google Scholar 

  96. Castresana J. GBLOCKS: selection of conserved blocks from multiple alignments for their use in phylogenetic analysis. Version 0.91 b. Copyrighted by Castresana J. EMBL; 2002.

    Google Scholar 

  97. Legendre P, Lapointe F-J. Assessing congruence among distance matrices: single-malt scotch whiskies revisited. Australian and New Zealand J Statistics. 2004;46:615–29.

    Article  Google Scholar 

  98. Campbell V, Legendre P, Lapointe F-J. Assessing congruence among ultrametric distance matrices. J Classif. 2009;26:103–17.

    Article  Google Scholar 

  99. Paradis E, Claude J, Strimmer K. APE: analyses of phylogenetics and evolution in R language. Bioinformatics. 2004;20:289sific.

    Article  Google Scholar 

  100. Tamura K, Peterson D, Peterson N, Stecher G, Nei M, Kumar S. MEGA5: molecular evolutionary genetics analysis using maximum likelihood, evolutionary distance, and maximum parsimony methods. Mol Biol Evol. 2011;28:2731–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  101. Katoh K, Rozewicki J, Yamada KD. MAFFT online service: multiple sequence alignment, interactive sequence choice and visualization. Brief Bioinform. 2019;20(4):1160–6.

    Article  CAS  PubMed  Google Scholar 

  102. Trifinopoulos J, Nguyen LT, von Haeseler A, Minh BQ. W-IQ-TREE: a fast online phylogenetic tool for maximum likelihood analysis. Nucleic Acids Res. 2016;44(W1):W232–5

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  103. Nguyen LT, Schmidt HA, von Haeseler A, Minh BQ. IQ-TREE: a fast and effective stochastic algorithm for estimating maximum-likelihood phylogenies. Mol bioland evol. 2014;32(1):268–74.

    Article  CAS  Google Scholar 

  104. Anisimova M, Gascuel O. Approximate likelihood-ratio test for branches: a fast, accurate, and powerful alternative. Syst Biol. 2006;55(4):539–52.

    Article  PubMed  Google Scholar 

  105. Ronquist F, Teslenko M, van der Mark P, Ayres DL, Darling A, Höhna S, Larget B, Liu L, Suchard MA, Huelsenbeck JP. MrBayes 3.2: efficient Bayesian phylogenetic inference and model choice across a large model space. Syst Biol. 2012;61(3):539–42.

    Article  PubMed  PubMed Central  Google Scholar 

  106. Miller MA, Pfeiffer W, Schwartz T. The CIPRES science gateway: a community resource for phylogenetic analyses. In: Proceedings of the 2011 TeraGrid conference: extreme digital discovery. ACM; 2011. p. 41.

    Google Scholar 

  107. Maddison WP, Maddison DR. Mesquite: a modular system for evolutionary analysis; 2018.

    Google Scholar 

  108. Bouckaert R, Heled J, Kühnert D, Vaughan T, Wu C-H, Xie D, et al. BEAST 2: a software platform for Bayesian evolutionary analysis. PLoS Comput Biol. 2014;10(4):e1003537 pmid:24722319.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  109. Rambaut A, Drummond AJ, Xie D, Baele G and Suchard MA. Posterior summarisation in Bayesian phylogenetics using Tracer 1.7. Syst Biol. 2018;67(5):901-4.

Download references


We would like to thank Marie-Hélène Carignan for her help in the laboratory and Dr. Brian Boyle from IBIS at Laval University (Quebec City, QC, Canada) for his help in designing the work and in carrying out the analyses of assembly.


This work was realized with the financial support of (i) a UQAR BORÉAS collaborative research program and a QCBS Seed Grant (#56 entitled « Mitogenomics of polychaete worms of the Ophryotrocha genus: an emerging model for multigenerational experiments in marine invertebrates ») granted to AT, GMN, FV, FD and PC by the University of Quebec in Rimouski; (ii) a Natural Sciences and Engineering Research Council of Canada Discovery Program grant (RGPIN- 06349) awarded to FD; (iii) a Natural Sciences and Engineering Research Council of Canada Discovery Program grant (RGPIN- 2020-05627), the Programme Établissement de Nouveaux Chercheurs Universitaires of the Fonds de Recherche du Québec—Nature et Technologies (No.199173), and the Fonds de recherche de l’Université du Québec à Rimouski, all awarded to PC. GMN was supported by a Marie Skłodowska-Curie Action (Individual Fellowship, Proposal Number: 659359). FV and PC are members of Quebec Ocean FRQNT-funded research excellence networks. AT, FV, FD and PC are members of QCBS FRQNT-funded research excellence networks.

Author information

Authors and Affiliations



AT, GMN, FV, FD and PC conceptualized and designed the work and analyses to be carried out. NB carried out the laboratory analyses. AT carried out the bioinformatic analyses with support from MM. AT wrote the first draft of this manuscript, which was further developed with the support of GMN, FD and PC. All authors have contributed and gave final approval for publication to the last version of this MS.

Corresponding author

Correspondence to Piero Calosi.

Ethics declarations

Ethics approval and consent to participate

All applicable international, national and/or institutional guidelines for the care and use of animals were followed.

Consent for publication

Not applicable.

Competing interests

The authors declare no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1.

Main literature on the genus Ophryotrocha.

Additional file 2.

Genome annotation for Ophryotrocha adherens.

Additional file 3.

Genome annotation for Ophryotrocha diadema.

Additional file 4.

Genome annotation of Ophryotrocha japonica.

Additional file 5.

Genome annotation of Ophryotrocha labronica.

Additional file 6.

Genome annotation of Ophryotrocha puerilis.

Additional file 7.

Genome annotation of Ophryotrocha robusta.

Additional file 8.

tRNA structures in O. adherens.

Additional file 9.

tRNA structures in O. diadema.

Additional file 10.

tRNA structures in O. japonica.

Additional file 11.

tRNA structures in O. labronica.

Additional file 12.

tRNA structures in O. puerilis.

Additional file 13.

tRNA structures in O. robusta.

Additional file 14.

Matrix of the comparison between gene order of the six Ophryotrocha species investigated.

Additional file 15.

Matrix of the comparison between gene order of Ophryotrocha and other annelids.

Additional file 16.

List of Ophryotrocha species used for the phylogeny.

Additional file 17.

Ancestral reproductive mode reconstruction based on Bayesian posterior probabilities.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Tempestini, A., Massamba-N’Siala, G., Vermandele, F. et al. Extensive gene rearrangements in the mitogenomes of congeneric annelid species and insights on the evolutionary history of the genus Ophryotrocha. BMC Genomics 21, 815 (2020).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI:


  • Molecular phylogeny
  • Dorvilleidae
  • Mitogenome
  • Next generation sequencing
  • Model species
  • Reproductive mode