Three Cases Revealing Remarkable Genetic Similarity Between Vent-Endemic Rimicaris Shrimps Across Distant Geographic Regions: Toward a New Conservation Perspective

1. Introduction

Deep-sea hydrothermal vents, which are found along mid-ocean ridges and back-arc basins, are characterized by darkness, temperatures over 300°C, and chemically enriched fluids (Ramirez-Llodra et al., 2007; Vrijenhoek, 2010; Sogin et al., 2021). These environments host unique ecosystems sustained by chemosynthetic microbiomes as primary producers, supporting high biodiversity and considerable biomass (McArthur and Tunnicliffe, 1998; Van Dover et al., 2002; Fisher et al., 2007; Moalic et al., 2012; Tunnicliffe et al., 2024). Notably, these ecosystems are shaped by evolutionary forces that promote symbiotic adaptations between chemosynthetic bacteria and vent invertebrates, contributing to the high levels of endemism observed among vent-inhabiting invertebrates (Won et al., 2003; Cavanaugh et al., 2006; Dubilier et al., 2008; Mullineaux et al., 2018; Breusing et al., 2022, 2023). While often regarded as remote and isolated, these ecosystems are attracting growing commercial interest due to the potential of deep-sea mineral resources.

Hydrothermally active vents are both biologically rich and geochemically significant, and are also sites of polymetallic sulfide deposits containing valuable metals such as copper, zinc, gold, silver, and rare earth elements (German et al., 2016; Van Dover et al., 2018; Fuchs et al., 2019). As interest in deep-sea mineral extraction grows, these habitats face increasing threats from human activities, particularly commercial deep-sea mining. Despite their ecological and economic importance, our understanding of how such disturbances will affect hydrothermal vent biodiversity remains limited, and uncertainties persist regarding species connectivity, gene flow barriers, dispersal capacity, and the potential for population recovery following habitat degradation (Van Dover, 2014; Mullineaux et al., 2020; van der Most et al., 2023). Expanding our knowledge in these areas is essential for developing effective conservation strategies prior to the commencement of large-scale exploitation. Elucidating the genetic structures within and between dominant vent species is a critical step toward informed management and effective conservation of hydrothermal vent ecosystems.

Among the dominant members of vent communities, the caridean family Alvinocarididae is a key indicator of genetic connectivity and ecosystem resilience (Teixeira et al., 2013; Sun et al., 2018; Dai et al., 2025). This family is one of the most abundant crustacean taxa in hydrothermal vent ecosystems worldwide, comprising 35 species in five genera: Alvinocaris, Rimicaris, Mirocaris, Nautilocaris, and Keldyshicaris (Komai and Segonzac, 2005; Martin and Haney, 2005; WoRMS Editorial Board, 2025). Rimicaris currently includes 15 valid species, all of which are geographically restricted to specific vent regions (Nye et al., 2012; Vereshchaka et al., 2015; Methou et al., 2024b). To understand the apparent endemism of Rimicaris species, previous studies have examined their trophic strategies and symbiotic associations with chemoautotrophic microbial communities (Assié, 2016; Apremont et al., 2018; Lee et al., 2021; Methou et al., 2023, 2024b). This genus also exhibits substantial morphological diversity and has recently undergone major taxonomic revisions, with six formerly distinct genera now synonymized under Rimicaris (Vereshchaka et al., 2015; Methou et al., 2024a). Despite the morphological and ecological diversity within Rimicaris, our pilot study revealed unexpectedly high genetic similarity between some species from geographically distant ocean basins, indicating that their gene flow and dispersal patterns may be more complex than previously understood. However, the mechanisms driving these patterns, particularly those related to genetic divergence, larval dispersal, and migration, remain poorly understood.

The COI barcode is an approximately 700-bp region of the mitochondrial cytochrome c oxidase I (COI) gene that has been widely used for species identification and population genetic studies in metazoans due to its high substitution rate and the availability of universal primers (Hebert et al., 2003a). However, in cases where closely related congeners have recently diverged from a common ancestor, COI barcoding may fail to distinguish species due to insufficient sequence divergence (Shearer and Coffroth, 2008; Matzen da Silva et al., 2011; Ranasinghe et al., 2022). To address these limitations, researchers have increasingly adopted multi-locus approaches such as incorporating additional mitochondrial and nuclear gene regions or analyzing whole mitogenomes or nuclear genomes (Hinsinger et al., 2015; Chagas et al., 2020; Nymoen et al., 2024).

In this study, we first confirmed three cases of paired Rimicaris species exhibiting high genetic affinity based on COI barcode sequences. Then, we assessed the sequence similarity and genetic connectivity within each species pair. To validate these affinity data at broader genetic levels, we conducted comparative analyses using the mitochondrial 16S rRNA (16S) gene, nuclear histone 3 (H3) gene, and all 13 mitochondrial protein-coding genes (PCGs). Our findings are anticipated to have implications relating to genetic connectivity, potential migration pathways, and speciation processes.

2. Materials and Methods

2.1. Ethics Approval

The Korea Institute of Ocean Science and Technology obtained permission to collect vent fauna, including shrimps, from hydrothermal vent regions in the Southwestern Pacific Ocean located within the Exclusive Economic Zones of Fiji and Tonga. Approval was granted by the Ministry of Land and Natural Resources of the Republic of Fiji and the Ministry of Lands, Survey and Natural Resources of the Kingdom of Tonga. Shrimps from the Manus Basin, which were loaned to Duke University, were collected with permission from the Government of Papua New Guinea.

2.2. Vent Shrimp Sampling and Identification

Alvinocaridid shrimp specimens were collected using suction samplers mounted on remotely operated vehicles from nine vent sites in the southwestern Pacific (SWP) and two vent sites in the northern Central Indian Ridge (nCIR) (Figure 1). On board the research vessel, all specimens were immediately preserved in 95% ethanol or stored at –80°C for genetic analysis. The specimens were identified based on morphological characteristics and COI barcodes (Hebert et al., 2003a; Komai and Tsuchida, 2015). Detailed information on the specimens is provided in Table S1.

2.3. DNA Extraction, Partial Gene and Mitogenome Sequencing, and Sequence Data Preprocessing

A small amount of muscle tissue was dissected from a pereopod of each specimen for DNA extraction. Total genomic DNA was extracted using the QIAamp Fast DNA Tissue kit (Qiagen, Hilden, Germany).

Partial sequences of the COI, 16S, and H3 genes were amplified using published universal primers (Table S2). Polymerase chain reaction (PCR) amplification was performed in a total volume of 50 μL containing 1 μL of genomic DNA, 4 μL of dNTP mixture (2.5 mM each), 1 μL of each primer (10 pmol), 5 μL of 10× Ex Taq Buffer (Mg²⁺ plus), and 1.25 U of Takara Ex Taq DNA Polymerase (Takara Bio, Kusatsu, Japan), with an initial denaturation at 94°C for 2 min, followed by 35 cycles of denaturation at 95°C for 10 s, primer annealing at 48°C for COI, 46°C for 16S rRNA, and 50°C for H3 for 30 s, and extension at 72°C for 1 min, and a final 5-min extension at 72°C. The PCR products were sequenced by Macrogen (Seoul, Korea) using an ABI 3730xl Analyzer (Applied Biosystems, Waltham, MA, USA) with BigDye Terminator v3.1 Cycle Sequencing Kits (Applied Biosystems). Newly obtained sequences were trimmed, annotated, and aligned with Geneious Prime v2023.0.1 (Biomatters, Auckland, New Zealand) and adjusted manually by visual inspection.

For mitogenome sequencing, mitochondrial DNA was amplified using the REPLI-g Mitochondrial DNA Kit (Qiagen). Libraries were prepared with a TruSeq Nano DNA Kit (Illumina, San Diego, CA, USA) and short-read sequencing was performed using the Illumina HiSeq 4000 platform at Macrogen. The mitogenome was assembled using NOVOPlasty v4.3.1 (Dierckxsens et al., 2017), annotated using MITOS2 (Bernt et al., 2013), and curated manually in Geneious Prime v2023.0.1 (Biomatters).

The newly generated sequences were registered in GenBank (Table S3).

2.4. Tree Construction, Nucleotide Divergence, Haplotype Network, and Gene Flow

Based on both newly generated sequences and those retrieved from GenBank, genetic divergence was calculated using the p-distance method, and a neighbor-joining (NJ) tree was constructed with MEGA 11 (Tamura et al., 2021).

To visualize genetic similarities and dissimilarities among samples, principal coordinate analysis (PCoA) was performed based on distance matrices using GenAlEx v6.503 (Peakall and Smouse, 2006).

The number of polymorphic sites, number of haplotypes, haplotype diversity, nucleotide diversity, Tajima’s D, Fu’s F_S, and fixation index based on pairwise differences (F_ST) were estimated using DnaSP v5.10.01 and Arlequin v3.5.2.2 (Librado and Rozas, 2009; Excoffier and Lischer, 2010). To determine the genetic relationships between paired species within each clade, haplotype networks were created using TCS and visualized with Hapsolutely v0.2.2 (Templeton et al., 1992; Vences et al., 2024).

The gene flow between closely related species was estimated as the number of migrants per generation (Nm) using MIGRATE-N 5.0.4 (Beerli, 2006; Beerli et al., 2019). Nm was calculated as Nm = θ × M, where N is the effective population size, m is the migration rate, θ is the mutation-scaled population size, and M is the mutation-scaled migration rate.

2.5. Mitogenome Sequence Comparison

The nucleotide and amino acid sequence similarities of mitochondrial genes were calculated using Geneious Prime v2023.0.1 (Biomatters). The ratio of nonsynonymous to synonymous substitutions (K_a/K_s) was measured using KaKs_Calculator v3.0 with the Yang–Nielsen model (Zhang, 2022).

3. Results

3.1. Datasets Prepared from Multi-Gene Sequences

We generated new sequences of the COI and 16S mitochondrial genes and H3 nuclear gene, as well as complete mitogenome sequences for Rimicaris variabilis and Rimicaris cf. variabilis (Tables S3).

Each gene was aligned individually using both the newly obtained sequences and those of Rimicaris chacei, Rimicaris hybisae, Rimicaris exoculata, and Rimicaris kairei retrieved from GenBank (Table S3). We were unable to include H3 sequences for R. exoculata and R. kairei, and mitogenome sequences for R. chacei and R. hybisae in our analyses because they were not available in public sequence databases.

3.2. Genetic Clusters of Rimicaris Species

Based on the NJ tree constructed from the partial gene datasets for six Rimicaris species, three distinct clades were identified (Figure 2). In each clade, the paired species originated from geographically distant oceanic regions or ridge systems: Clade I included R. chacei from the Mid-Atlantic Ridge (MAR) and R. hybisae from the Mid-Cayman Spreading Center (MCSC); Clade II included R. exoculata from MAR and R. kairei from the Carlsberg Ridge–Central Indian Ridge (CR-CIR); and Clade III included R. variabilis from SWP and R. cf. variabilis from nCIR.

The genetic divergence between paired species within each clade was 0.35–1.90% for COI and 0.04–0.30% for 16S (Table 1). By contrast, inter-clade divergences were substantially higher, at 6.95–8.74% for COI and 0.30–1.10% for 16S. The H3 marker lacked sufficient resolution to distinguish genetic divergence within or between clades, showing only a single nucleotide difference at the same position across all Rimicaris sequences.

PCoA based on distance matrices of COI revealed that the first principal coordinate (PCo1) accounted for 51.14% of the total genetic variance, and the second (PCo2) accounted for 30.15%. The resulting PCoA plot clearly separated the three clades, supporting the NJ tree topologies (Figure 3).

3.3. Genetic Connectivity Between Paired Rimicaris Species Within Each Clade

Based on COI sequence divergence, paired Rimicaris species within each clade fell within the range of intraspecific variation, as defined by species delimitation thresholds in DNA barcoding studies (Hebert et al., 2003b). To assess whether the paired species in each clade should be considered conspecific, we examined genetic clustering using haplotype networks, genetic structures, and gene flow estimates (Figure 4 and Figure 5, Table 2).

In Clade I, the haplotype network revealed four dominant haplotypes shared between R. chacei and R. hybisae, indicating overlapping genetic pools (Figure 4a). This result was supported by gene flow estimates, which showed moderate bidirectional exchange between the two species (Nm = 3.31 and 3.92; Figure 5a). Despite their shared haplotypes, the two species showed substantial genetic differentiation, with a pairwise F_ST of 0.47 (Table 2). Both species also had negative Tajima’s D and Fu’s F_S values, consistent with recent independent population expansion.

In Clade II, the haplotype network revealed no shared haplotypes between R. exoculata and R. kairei, indicating clear genetic separation (Figure 4b). This result was supported by gene flow estimates, which showed low bidirectional gene flow between the two species (Nm = 0.65 and 0.01; Figure 5b), and by genetic differentiation analyses, which yielded a high pairwise F_ST of 0.82 (Table 2). Both species also had negative Tajima’s D and Fu’s F_S values, suggesting recent independent population expansion. Notably, in each species, expansion appears to have occurred from a distinct ancestral haplotype.

In Clade III, the COI haplotype network had limited resolution in distinguishing R. variabilis and R. cf. variabilis, likely due to high haplotype diversity (~1.0) at the COI marker (Figure 4c, Table 2). Nevertheless, we detected a haplotype cluster shared by the two species. By contrast, the 16S haplotype network revealed a clearer structure, with a single dominant haplotype shared by both species and additional haplotypes restricted to R. variabilis (Figure 4d). The two species showed low genetic differentiation, with a pairwise F_ST of 0.10 (Table 2). Gene flow analysis revealed strong asymmetry, as migration from R. cf. variabilis (nCIR) to R. variabilis (SWP) was exceptionally high (Nm = 79.92) and the reverse flow was very low (Nm = 0.82; Figure 5c). Furthermore, R. variabilis in the SWP showed evidence of recent rapid population expansion, with θ = 0.097 (the highest value among the six Rimicaris species), negative Tajima’s D and Fu’s F_S values, and very high haplotype diversity (H_d = 0.96).

3.4. Mitogenomic Similarity Between Paired Rimicaris Species

Mitogenomic-level genetic similarities between the paired Rimicaris species within each clade were examined (Table 3). Within Clade II and Clade III, the amino acid sequences of all PCGs were nearly identical, with consistent gene lengths, start and stop codons, and overall sequence composition. Most nucleotide differences between paired species were synonymous substitutions, with a few notable non-synonymous changes observed in ND3 and ND6 between R. exoculata and R. kairei in Clade II.

4. Discussion

4.1. Clade-Specific Patterns of Genetic Similarity in Rimicaris

On longer geological timescales, tectonic processes, including the formation and subsequent breakup of Pangaea from the Paleozoic to Mesozoic eras, are thought to have profoundly influenced deep-sea migration, diversification, persistence, and extinction patterns in marine lineages (McClain and Hardy, 2010; Lins et al., 2012). The discontinuous distribution of deep-sea hydrothermal vents along mid-ocean ridges reinforces geographic isolation among vent communities (Bachraty et al., 2009; Vrijenhoek, 2010; Mullineaux et al., 2018). As a result, vent ecosystems exhibit remarkable regional endemism, with more than 85% of the approximately 700 known species considered endemic (McArthur and Tunnicliffe, 1998; Wolff, 2005). However, there are exceptions to these regional endemic patterns. For example, the vent mussel Bathymodiolus septemdierum and vent barnacle Leucolepas longa occur in widely separated vent fields (Watanabe et al., 2018; Mao et al., 2025). We also observed similar cross-regional connectivity in Rimicaris species. Each genetic clade of Rimicaris comprised two species originating from different oceanic regions or ridges, rather than from the same geographic area (Figure 2). Moreover, the three clades had distinct biological and ecological traits (Table 4).

In Clade I, R. chacei (MAR) and R. hybisae (MCSC) are morphologically similar in the adult stage, particularly in the final stomach volume and surface area of their mouthparts, suggesting that they represent a single species (Teixeira et al., 2013; Methou et al., 2024b). These two taxa share four haplotypes, two putatively originating from MAR and two from MCSC, and show bidirectional gene flow, indicating some genetic overlap (Figures 4a and 5a). A potential dispersal route connecting the MCSC and MAR is the Windward Passage of the Caribbean Sea (Connelly et al., 2012; Yearsley et al., 2020). However, the comprehensive genetic and ecological data analyzed in this study did not support uninterrupted bidirectional gene flow between the two regions, suggesting the presence of barriers to genetic exchange (Table 2 and Table 5). Thus, while the shared haplotypes suggest either historical connectivity or limited ongoing gene flow, the overall genetic structure implies partial reproductive isolation. At a finer local scale, the Snake Pit vent field on MAR appears to act as a keystone site, facilitating gene flow to other MAR sites and MCSC, but not in the reverse direction (Figure S2). Notably, the two MCSC-derived haplotypes were not included in specimens from the Snake Pit, whereas the two MAR-derived haplotypes were found in specimens from all sampled regions. This asymmetric pattern of gene flow and haplotype distribution resembles a ring species model, in which adjacent populations can interbreed, but the terminal populations become reproductively isolated due to geographic separation or asynchronous mating behaviors (Irwin, 2000; Irwin et al., 2005). However, further data on mating timing, larval dispersal, and reproductive compatibility are needed to confirm whether this system truly fits a ring species dynamic.

In Clade II, R. exoculata (MAR) and R. kairei (CR-CIR) are morphologically, genetically, ecologically, and biogeographically distinct (Figures 4b and 5b, Table 2 and Table 4). Morphologically, they differ in the conspicuousness of the surface setae on the carapace, and in the size of the pereopods and antennal flagellae (Watabe and Hashimoto, 2002). Genetically, they share no COI haplotypes, exhibit low bidirectional gene flow (Nm < 1), and show significant differentiation in pairwise F_ST values. However, while the COI and 16S markers indicate some divergence, this evidence is insufficient for species-level separation, and the mitogenome sequences of these two species are very similar, differing only by a few synonymous substitutions (Table 3). This limited sequence variation suggests that R. exoculata and R. kairei are likely in the early stages of speciation, i.e., incipient speciation. Divergence time estimates based on mitogenomes indicate that their common ancestor split into the two lineages approximately 5.38 Mya, whereas Alvinocarididae and Rimicaris originated much earlier, at 69.36 and 28.50 Mya, respectively (Sun et al., 2024). During the Miocene, the closure of the Tethys Ocean, which once connected the Mediterranean and Indian Oceans, has been implicated in driving vicariant speciation among marine taxa with Atlantic–Mediterranean–Indian distributions (Hrbek and Meyer, 2003; Hou and Li, 2018; Liu et al., 2018; Bialik et al., 2019; Agiadi et al., 2024). A plausible scenario for Clade II is that its common ancestor dispersed freely through the Tethys Sea, and its descendants became isolated in the Atlantic and Indian Oceans following the closure of the sea, initiating their subsequent divergence.

4.2. Adaptive Divergence, Eastward Dispersal, and Regional Barriers in Clade III

In Clade III, Rimicaris variabilis (SWP) and R. cf. variabilis (nCIR) exhibit a more complex evolutionary trajectory than those in Clades I and II. Notably, R. variabilis in the SWP utilizes a distinct energy source, reflected in its markedly low δ¹³C value, which is clearly separated from the other five Rimicaris species (Figure S4). This separation likely reflects differences in the chemosynthetic carbon fixation pathways of their symbionts and associated trophic interactions—specifically, Calvin–Benson–Bassham (CBB) cycle-based versus the reductive tricarboxylic acid (rTCA) cycle-based primary production (Suh et al., 2022b). These metabolic differences further suggest adaptive divergence in nutritional strategies driven by symbiotic associations: R. variabilis derives only partial nutrition from its cephalothoracic epibionts, whereas R. cf. variabilis depends on them entirely (Table 4; Lee et al., 2021).

However, we could not find any noticeable differences in rostrum or tail morphology between the Manus and other SWP populations of R. variabilis (unpublished data). In addition, we were unable to assess the morphological characteristics of R. cf. variabilis because the available specimens were damaged. Although sequence data from only nine R. cf. variabilis specimens were used in this study, 19 individuals were examined in total. Despite this limited sample size, our COI, 16S, and complete mitogenome analysis results consistently supported their classification as a single species (Table 1, Table 2 and Table 3). A notable finding is the highly asymmetric gene flow between regions; migration from nCIR to SWP was strong, whereas that in the reverse direction was almost negligible (Figure 5c). This asymmetry, combined with higher genetic diversity in SWP than in nCIR, suggests a historical or ongoing unidirectional dispersal route from nCIR into the SWP, followed by genetic expansion within the SWP population. Importantly, this eastward dispersal contrasts previous hypotheses proposing westward migration from the Pacific to the Indian Ocean, including the Pacific origins of alvinocaridid shrimps, neolepadid barnacles, and the vent mussel Bathymodiolus septemdierum (Moalic et al., 2012; Sun et al., 2018; Watanabe et al., 2018; Mao et al., 2025).

At a finer local scale, more complex genetic structures are evident within the SWP, particularly between the Manus Basin and other SWP sites such as the North Fiji Basin, Tonga Arc, and Futuna Arc. In the 16S haplotype network, haplotypes observed only in R. variabilis found predominantly in Manus Basin individuals (Figure 4d), and this population consistently showed the highest genetic diversity (Figures S1 and S3). Based on COI sequences, both the highest θ value and greatest intraspecific variation within R. variabilis were attributed to the Manus Basin population. Remarkably, only a single haplotype was shared between the Manus Basin and other SWP sites, a pattern that may reflect parallel mutations arising independently in each region rather than true connectivity. Similar trends have been reported in other vent-endemic taxa, including the gastropod Ifremeria nautilei, limpet Lepetodrilus aff. schrolli and Austinograea crabs (Thaler et al., 2011, 2014; Lee et al., 2019; Plouviez et al., 2019). These findings indicate the presence of a strong regional biogeographic barrier between the Manus Basin and other SWP vent systems, shaped by dispersal limitations imposed by the tectonic history and ocean circulation (Mitarai et al., 2016; Tunnicliffe et al., 2024).

Overall, these results support a plausible scenario for Clade III in which CIR was the source population, with dispersal proceeding along two routes: a strong direct pathway into the Manus Basin, and a weaker pathway through the Southern Ocean toward other SWP sites. These two lineages appear to have remained largely unconnected, with the Manus Basin lineage accumulating extensive genetic diversity, whereas populations along the latter route remained more interconnected across other SWP vent fields.

4.3. New Perspective on Vent Organism Conservation

Endemism in vent and seep fauna is thought to have occurred since the Paleozoic (Little et al., 1998; Campbell, 2006; Kiel, 2010). Much of the contemporary biogeographic endemism observed in vent and seep invertebrates is thought to have been shaped by Cenozoic tectonic events and changes in oceanic circulation (< 100 Mya; Van Dover et al., 2002). Within these groups, the alvinocaridid genus Rimicaris is among the most recently diversified lineages, with 16 described species originating since the Paleocene (Sun et al., 2024). Most of these species are distributed across hydrothermal vent fields in the Pacific and Indian Oceans, reflecting diversification over a relatively short tectonic history (Vereshchaka et al., 2015; Komai et al., 2016; Komai and Giguère, 2019; Methou et al., 2024a).

Understanding species diversity, distribution, genetic diversity, and connectivity provides a critical foundation in developing effective conservation strategies for hydrothermal vent ecosystems (Van Dover et al., 2002; Van Dover, 2012; Breusing et al., 2023; Tunnicliffe et al., 2024). Historically, migrations, distribution ranges, and the biogeographic structuring of vent fauna have been inferred largely from the faunal compositions of local communities (Mullineaux et al., 2018; Perez et al., 2021; Zhou et al., 2022). In the context of potential seabed mining, conservation plans have therefore emphasized the low connectivity among vent networks and the high degree of endemism, as species recruitment is unlikely to extend across vent system boundaries (Thomas et al., 2021; Tunnicliffe et al., 2024).

However, the present study reveals evidence of cross-regional connectivity among alvinocaridid shrimps spanning different oceanic regions and ridges. This finding challenges the prevailing view that vent species are strictly confined within provincial boundaries, and instead highlights that while some taxa are strongly constrained by dispersal barriers, others maintain connectivity on broader scales.

Our results suggest that effective conservation of vent ecosystems should be framed from a global perspective, rather than being restricted to single species or narrowly defined provinces. Treating vent ecosystems as a single homogeneous management unit risks overlooking their complex evolutionary and ecological dynamics. Effective conservation strategies should also recognize the importance of distinct biogeographic provinces, ensuring that management simultaneously addresses both global-scale connectivity and local-scale endemism.