Deep-sea chemosynthetic ecosystems are “oases” fueled by bacterial primary production (Van Dover et al., 2001; Dubilier et al., 2008), such as deep-sea hot vents, cold seeps, and organic falls like sunken woods and whale falls (Martin et al., 2008; He et al., 2023). Many species inhabiting these chemosynthesis-based ecosystems are specialists that cannot be found in other habitats (Mills and Harrison, 1998). Hydrothermal vent specialists are especially notable, as vents are often distributed like stepping stones along mid-ocean ridges species with dispersal most likely mediated by currents along the mid-ocean ridges (Tunnicliffe and Mary R. Fowler, 1996; Rogers et al., 2012), although other systems in settings like hot spots and back-arc basins are also present (Breusing et al., 2016; Xu et al., 2021; Brunner et al., 2022).

Neomphalida is an order of gastropod molluscs endemic to chemosynthetic ecosystems, forming the subclass Neomphaliones together with the cocculiniform gastropods in the order Cocculinida (Ponder et al., 2019; Uribe et al., 2022). Neomphalida includes a single superfamily Neomphaloidea with three families: Melanodrymiidae, Neomphalidae, and Peltospiridae. Among these three families, Melanodrymiidae has been found on sunken wood and hydrothermal vents, while living Neomphalidae and Peltospiridae are vent specialists (Mclean, 1989; Heß et al., 2008; Chen et al., 2015b; Chen et al., 2015c; Zhong et al., 2022). One genus, Retiskenea, is known from both recent and fossil hydrocarbon seeps as early as the Early Jurassic with unclear familial affinity, although it has been tentatively assigned to Neomphalidae (Warén and Bouchet, 2001; Kaim et al., 2014; Kaim, 2022). Neomphaloideans are morphologically and ecologically diverse, ranging from filter feeding to grazing to relying on endosymbionts (Sasaki et al., 2010). The latter is seen in the two coiled genera Chrysomallon and Gigantopelta (Chen et al., 2015a; Chen et al., 2015b; Sun et al., 2020), which have independently evolved (Chen et al., 2015c; Zhong et al., 2022) endosymbiosis in an enlarged and modified esophageal gland (Lan et al., 2021; Lan et al., 2022). To better understand the evolution of these traits, a robust phylogeny is required.

Reconstructing the internal phylogenetic relationships among different taxa within Neomphalida has proved challenging with different studies finding disparate topologies. For instance, multi-gene and mitogenome-based phylogenetic studies have placed the Scaly-foot Snail Chrysomallon squamiferum in Peltospiridae, which agrees with its morphology, but it has been recovered nested within Neomphalidae when only the mitochondrial COI gene was used (Chen et al., 2017). Monophyly of the three Neomphalida families was sometimes supported by molecular phylogeny (Heß et al., 2008) but not in others (Aktipis and Giribet, 2010; Chen et al., 2017; Zhong et al., 2022). These studies were performed either using a small number of genes or with a limited taxon coverage. For instance, in the only mitogenome-based molecular phylogeny study focusing on this group, Melanodrymiidae was unrepresented (Zhong et al., 2022). In this study, we expanded the mitogenome sampling to the whole order by adding a total of 11 mitogenomes from all three families. We are the first to include the five mitogenomes from Melanodrymiidae, with phylogenetic analyses performed using various methods, matrix, and models for robustness.

A total of 11 mitogenomes were assembled in Neomphalida in this study, and the mitogenomes of Melanodrymia laurelin, Melanodrymia galeronae, Melanodrymia brightae, Melanodrymia aurantiaca, Cyathermia naticoides, and Symmetromphalus regularis are circular. The mitogenomes ranged from 15,269 bp (Peltospira operculata) to 18,219 bp (Melanodrymia laurelin) in length ( Tables 1 , 2 ). All mitochondrial genomes contain exclusively 13 protein-coding genes, 22 tRNA genes, and two rRNA genes (12S and 16S), without any gene duplication or deletion.

By applying the concatenated and coalescent methods of phylogenetic tree construction, we obtained a total of 13 phylogenetic trees, comprising two different tree topologies with conflicting phylogenetic signals. In topology 1 ( Figure 1A ), Melanodrymiidae is sister to (Neomphalidae + Peltospiridae), supported by the following methods: 13PCGs with MFP, C60, RAxML, PB, ASTRAL − MP, and PCGs with BP. In topology 2, however, Peltospiridae was recovered sister to (Neomphalidae + Melanodrymiidae), supported by the following methods: 13PCGs with MFP, RAxML, ASTRAL – MP and PCGs + rRNA with MFP, RAxML, BP, ASTRAL − MP ( Figure 1B ). The three Neomphalida families were monophyletic with maximal support across all methods, and the whole order Neomphalida.

To further test the phylogenetic signals of these two topologies on all the partitioned genes, 15 nucleotide sequences of the 13 PCGs + 2 rRNAs and 13 PCGs were separately tested on the two trees. For the 13 PCGs and two rRNAs, the gene tree from 10 of them supported topology 2 and the remaining five supported topology 1. Meanwhile, for the 13 PCGs, seven supported topology 1 and six supported topology 2 ( Figure 2 ).

For intra-family level relationships, the phylogenetic positions were relatively stable in both topologies. The lineage with contentious position was the clade of (Dracogyra + Gigantopelta), which was found to be either sister to ((Nodopelta + Peltospira) + Symmetriapelta) in topology 1 or to (Nodopelta + Peltospira) in topology 2. In addition, the phylogenetic position of Melanodrymia telperion was variable, being either sister to (Melanodrymia aurantiaca +Melanodrymia brightae) or to all other Melanodrymia species included ( Figure 1 ).

The gene order of each mitogenome was analyzed. Overall, the gene order within each family was consistent, and they shared the identical gene order of atp6–atp8–trnD–cox2–nad5–trnH–nad4–nad4L–trnT–trnS2 ( Figure 3 ). Comparing the gene order of each family with the hypothetical ancestral gene order of the whole Gastropoda, the gene order of both Neomphalidae + Peltospiridae was more similar with less translocation and transversion, while the gene order of Melanodrymiidae was more deviated with the translocation of cob–trnE–trnF–trnR–trnC–trnG–cox3–trnK–trnI–nad3–trnS1–nad2. Both Neomphalida and Cocculinida have a shared transversion of atp6–atp8–trnD–cox2 when compared with the hypothetical ancestral gene order.

The mitogenome GC content, mitogenome size, and substitution rate of the three families were also calculated ( Figure 4 ). Melanodrymiidae had the lowest GC content, which significantly differed from Peltospiridae and Neomphalidae (p < 0.01) ( Figure 4A ). There is a significant difference in amino acid substitution rates among the three families, with Peltospiridae exhibiting the highest value compared to Neomphalidae (p < 0.01) ( Figure 4B ). For the nucleotide substitution rate, Neomphalidae has the lowest value compared with other two families and differs significantly from Melanodrymiidae (p < 0.05) ( Figure 4C ). Melanodrymiidae was found to have the largest mitogenome size among the three families (p < 0.05) ( Figure 4D ). Among the three families, Peltospiridae showed a much wider range of GC content and substitution rates of nucleotide and amino acids. For the GC content, C. squamiferum had the highest (33.84%), while S. wareni had the lowest (27.81%). For the substitution rates of nucleotide and amino acids, G. aegis showed the lowest substitution rates (1.03 and 0.59, respectively), while P. operculata had the highest rates (1.13 and 0.74, respectively).

As one of the six subclasses in Gastropoda, Neomphaliones is relatively understudied compared to the other five, largely due to sample availability since species in this clade are mainly from deep-sea chemosynthetic ecosystems. By including 11 new mitogenomes and a total of 21 mitogenomes covering all three known families in Neomophlida, we tested their phylogenetic positions with a combination of different methods and matrixes, gene order analysis, GC content, and evolutionary rate analysis. As of now, 26 genera are included in Neomphalida, and our analyses cover 13 genera plus an undescribed genus (Neomphalidae gen. et sp. Hatoma sensu Zhong et al., 2022).

Our phylogenetic analysis has recovered the monophyly of the three Neomphalida families across all methods, indicating a good consistency between molecular phylogeny and morphology at the family level. This result is in line with a former phylogenetic tree analysis (Heß et al., 2008) but with better support on the nodes and improved taxon coverage. Nevertheless, the relationships among these three families remain unsettled even with complete mitogenome sequencing. By checking the phylogenetic signals on each partition of both nucleotide and amino acid sequences of our two conflicting topologies recovered, we show that 10 nucleotide partitions support topology 2 versus five supporting topology 1, while seven amino acid partitions supported topology 1 versus six supporting topology 2. There shall be very little effect of the third codon saturations on the molecular phylogeny, since the third codons was replaced by degenerated nucleotides. It is possible that the phylogenetic signal in the nucleotide data is lost in the amino acid dataset after translation for Neomphalida, and overall, it seems that topology 2 is better supported by our results. To further resolve the phylogenetic positions of these three families, a more comprehensive study using transcriptome or genome sequencing is needed. However, these methods require high-quality samples, which are usually lacking for many small and rare species in this order. In addition, more taxon sampling is also required, especially sunken wood melanodrymiid species like Leptogyra and Leptogyropsis with some plesiomorphic characters (Heß et al., 2008).

Our topology 2 is different from a former study (Heß et al., 2008), which included only 11 species with two partial genes (COI and histone H3). Based on their phylogenetic result and that some melanodrymiids show plesiomorphic characters inhabit the likely ancestral sunken wood habitat, the authors deduced a scenario that sunken wood could have served as an ecological stepping stone to hot vent. From our phylogenetic results, it seems equally likely that the subclass Neomphaliones diverged from non-chemosynthetic deep seafloor according to habitats, first splitting into the sunken-wood lineage (i.e., Cocculinida) and the hot-vent lineage (i.e., Neomphalida). Then, there could have been a secondary habitat shift from hydrothermal vents to sunken wood in Melanodrymiidae. This secondary habitat shift could be a plausible reason of the elevated nucleotide and amino acid substitution rate and the deviated mitogenome gene order in Melanodrymiidae when compared to Neomphalidae. This habitat transition scenario is complex and also differs from other deep-sea lineages, such as bathymodioline mussels, which took “wooden steps” to the vents and seeps (Lorion et al., 2013)—warranting future phylogenetic analyses using genomic-level data to test.

Species in Peltospiridae have a higher amino acid substitution rate (p < 0.05) and nucleotide substitution rate (though not very significant) comparing with the other two families, which may be related to its wider species distributions. By mapping the species distribution in combination with the phylogenetic analysis, it is clear that the Neomphalida species included are mostly from the eastern Pacific, especially Melanodrymiidae and Neomphalidae ( Figure 5 ) (Mclean, 1981). This is consistent with the distribution of known species of these two families, with Neomphalidae limited to the Pacific and Melanodrymiidae being largely Pacific except four species of Leptogyra from the Atlantic (Bush, 1897; Warén and Bouchet, 2001). For Peltospiridae, species in this family are also mostly found in the Pacific, but with species also distributing in the Indian Ocean, Atlantic Ocean, and Southern Ocean (East Scotia Ridge). Our tree suggests that peltospirids have colonized the Indian Ocean hot vents at least three times in historical dispersal. The earliest event was the clade leading up to Chrysomallon (Zhong et al., 2022), followed by the (Gigantopelta aegis + Dracogyra subfusca) clade and then Lirapex politus, in agreement with previous COI analyses (Chen et al., 2017). These three Indian Ocean colonization events may be linked to major tectonic events in the Indian Ocean (McElhinny, 1970), but the limited fossil records in the whole order means that it is difficult to carry out accurate molecular clock analysis with sufficient calibration nodes. Our results, combined with overall distribution pattern of Neomphalida, point to an eastern Pacific origin of this group.

The data presented in the study are deposited in the NCBI repository (https://www.ncbi.nlm.nih.gov/), under BioProject number PRJNA1048888.

The manuscript presents research on animals that do not require ethical approval for their study. No experiment on live animals were carried out in this study.

LZ: Data curation, Formal analysis, Methodology, Software, Validation, Visualization, Writing – original draft, Writing – review & editing. XG: Methodology, Software, Visualization, Writing – review & editing. CC: Data curation, Methodology, Visualization, Writing – review & editing. XH: Data curation, Methodology, Writing – review & editing. YQ: Methodology, Software, Writing – review & editing. JS: Conceptualization, Funding acquisition, Investigation, Project administration, Supervision, Writing – original draft, Writing – review & editing.