EDP Sciences logo
Submit your paper
Sciengine logo
Search
Menu
All issues Volume 39 (2026) Aquat. Living Resour., 39 (2026) 21 Full HTML
Open Access
Issue
Aquat. Living Resour.
Volume 39, 2026
Article Number 21
Number of page(s) 20
DOI https://doi.org/10.1051/alr/2026013
Published online 22 July 2026

© M. Sekino et al., Published by EDP Sciences 2026

Licence Creative CommonsThis is an Open Access article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1 Introduction

While taxonomic classification is paramount across a broad range of biological disciplines (Padial and De la Riva, 2021), a substantial proportion of taxa remain affected by taxonomic uncertainty, often arising from the presence of cryptic species or lineages (Bickford et al., 2006; Struck et al., 2018). Such uncertainty poses a major challenge for wildlife and ecosystem conservation. When multiple cryptic species are lumped under a single taxonomic name, biodiversity is underestimated, potentially leading to the oversight of endangered species (Frankham et al., 2002; Bickford et al., 2006). Resolving this issue is also essential for accurately defining conservation units and prioritizing conservation efforts (Dimmick et al., 1999; Fraser and Bernatchez, 2001; Frankham et al., 2002). Furthermore, from the perspectives of comparative phylogeography and historical demography (Hewitt, 2004; Knowles, 2004), sympatric cryptic species provide valuable systems for assessing shared responses to past environmental oscillations.

True oysters (superfamily Ostreoidea Rafinesque, 1815) play a critical role as ecosystem engineers (Ruesink et al., 2005; Smith and Pruett, 2025). They provide key ecosystem services, including water filtration (Gottlieb and Schweighofer, 1996), nutrient cycling through benthic-pelagic coupling (Newell et al., 2005), and the formation of structurally complex biogenic habitats that support diverse communities (Beck et al., 2011; Richardson et al., 2025). Given these important ecological roles, accurately characterizing species composition within an ecosystem—and thus species diversity—is a fundamental prerequisite for effective ecosystem management (Pirot et al., 2000; Richardson et al., 2026). However, species identification in oysters is often complicated by high shell plasticity, which renders morphology-based diagnosis unreliable. Consequently, molecular markers have become the standard tools for taxonomic classification in oysters, or at least for species delimitation in taxonomically challenging groups (Guo et al., 2018).

The increasing availability of DNA sequence data has partially alleviated taxonomic uncertainty (Hebert et al., 2003, 2004; Rubinoff, 2006), although methodological limitations persist (Rubinoff, 2006; Collins and Cruickshank, 2013). The advent of molecular techniques has greatly refined the taxonomy of true oysters, leading to the discovery of new species (Lam and Morton, 2003; Wu et al., 2013; Xia et al., 2014; Cui et al., 2021), the resolution of synonymy between nominal species (Anderson and Adlard, 1994; Boudry et al., 1998; Reece et al., 2008; Polson et al., 2009), and the systematic revision within the genus Crassostrea Sacco, 1897 (Salvi et al., 2014; Salvi and Mariottini, 2017; for counterarguments, see Bayne et al., 2017, 2019; Guo et al., 2018). A recent collection of DNA data has also prompted a reassessment of the conventional morphology-based systematic framework for the genus Saccostrea Dollfus and Dautzenberg, 1920—a group primarily distributed in tropical and subtropical coastal regions.

Within the genus Saccostrea, 19 nominal species have been described, 15 of which occur in the Indo-West Pacific (IWP) (summarized in Tab. 2 of Sekino and Yamashita, 2016). A more recent study added another IWP species, Saccostrea mordoides Cui et al., 2021 (see also below). Owing to the notoriously variable shell morphology in Saccostrea, the identities of many nominal species remain ambiguous (Supplementary Material 1), giving rise to longstanding taxonomic controversy. For instance, an extreme early view (Harry, 1985) classified all IWP species into a single taxon, Saccostrea cuccullata (Born, 1778). In contrast, Inaba and Torigoe (2004) and Huber (2010) recognized 10 and nine IWP species, respectively, with only five species consistently regarded as valid across both studies (Sekino and Yamashita, 2016).

To address this taxonomic complexity, Lam and Morton (2006) laid the groundwork for DNA-based classification of IWP Saccostrea oysters. Recognizing the extreme difficulty of assigning established species names to mitochondrial DNA lineages, they defined nine lineages within a S. cuccullata species complex, comprising S. cuccullata lineages A–G, Saccostrea kegaki, and Saccostrea glomerata, where the lineages S. kegaki and S. glomerata align with the recognized species Saccostrea kegaki (Torigoe and Inaba, 1981) and Saccostrea glomerata (Gould, 1850), respectively. Their analysis of DNA sequences and shell characteristics also distinguished Saccostrea mordax (Gould, 1850), which Huber (2010) suggested is a synonym of Saccostrea scyphophilla (Péron and Lesueur, 1807), from the S. cuccullata complex. Lam and Morton further identified two lineages within S. mordax (A and B). A third S. mordax lineage, “C,” was subsequently reported (Sekino and Yamashita, 2013) and later described as S. mordoides (Cui et al., 2021). However, more recent work (Lukehurst et al., 2026) revised the taxonomy within the S. mordax group, proposing that S. mordax lineages A and B should be assigned to S. scyphophilla, whereas S. mordoides is a synonym of S. mordax (Supplementary Material 1). Given these ongoing discrepancies, this group is hereafter referred to as the conventionally defined “S. mordax” (lineages A–C).

Because S. cuccullata sensu stricto occurs outside the IWP (from the West Indies to West Africa; Inaba and Torigoe, 2004, or from West Africa to the Arabian Peninsula; Huber, 2010), Sekino and Yamashita (2016) tentatively used the term “non-mordax” instead of S. cuccullata to represent the IWP lineages (see Tan et al., 2025 for a recent redescription of “true” S. cuccullata). This was motivated by the fact that S. mordax can be distinguished genetically and morphologically from other IWP Saccostrea oysters. In addition to the Lam and Morton’s “non-mordax” lineages A–G, Sekino and Yamashita (2016) reported three additional lineages (H‒J) from Iriomote Island in the Ryukyu Island Arc, Japan (Fig. 1). However, they failed to detect the lineage previously referred to as Saccostrea malabonensis lineage 2, which was recorded from both Iriomote Island and the neighboring island of Ishigaki (Hamaguchi et al., 2014). As the basis for assigning the name S. malabonensis to this lineage is unclear, it is here provisionally designated as lineage “K” to avoid further taxonomic confusion. More recent DNA barcoding studies focusing on the Australian coast (Snow et al., 2023; McDougall et al., 2024) demonstrated that all sampled Saccostrea specimens fell within these predefined lineages.

While links between most alphabetical non-mordax lineages (A‒K) and recognized IWP species remain uncertain (Supplementary Material 1), lineage J has been confidently synonymized with Saccostrea spathulata (Lamarck, 1819) by Sekino and Yamashita (2016) and McDougall et al. (2024). This assignment is primarily based on its exceptionally large size for a Saccostrea oyster (up to 123 mm in shell height), as well as its thick, deeply cupped shells with conspicuous chomata (fine notches on the inner shell margin) and well-developed growth lamellae (squamae). McDougall et al. (2024) further suggested that the Australian black-lip rock oyster—previously identified as Saccostrea echinata (Quoy and Gaimard, 1835) —is in fact S. spathulata. Accordingly, we refer to non-mordax lineage J as S. spathulata. In summary, DNA-based classifications recognize three S. mordax lineages (A–C) and 13 non-mordax lineages (A–I, K, S. glomerata, S. kegaki, and S. spathulata) in IWP Saccostrea oysters.

Among true oysters, Saccostrea is a representative genus in subtropical coastal ecosystems along the Ryukyu Island Arc, which encompasses Amami-Oshima, Okinawa, Ishigaki, and Iriomote Islands (Fig. 1). The occurrence of all three S. mordax lineages and nine non-mordax lineages (A, C, F‒I, K, S. kegaki, and S. spathulata) has been confirmed in this region (Lam and Morton, 2006; Sekino and Yamashita, 2013, 2016; Hamaguchi et al., 2014). However, given the site-specific nature of molluscan fauna in the Japanese subtropical islands (Nawa, 2008), along with limited geographic coverage of sampling and small sample sizes, Sekino and Yamashita (2016) suggested that phylogenetic diversity in the Ryukyu Island Arc may be underestimated. For example, non-mordax lineages B, D, and E have been recorded from Taiwan (Lam and Morton, 2006), but not from Iriomote Island, despite its close geographic proximity (Fig. 1). Given that the two islands are only approximately 200 km apart, it is plausible that more extensive surveys would detect some or all of these unrecorded lineages on Iriomote Island. Furthermore, expanded sampling efforts may reveal additional, previously unrecognized lineages.

The present study had three main objectives. First, we re-evaluated the phylogenetic diversity of Saccostrea oysters on Ishigaki and Iriomote Islands, which harbor the highest marine biodiversity within the entire Ryukyu Island Arc (Reimer et al., 2019). To this end, we sequenced four mitochondrial genes for our specimens, consistent with sequence data obtained by Sekino and Yamashita (2016). This survey builds on previous studies by Lam and Morton (2006), Hamaguchi et al. (2014), and Sekino and Yamashita (2013, 2016).

Second, we examined the genetic (mitochondrial) diversity within each lineage. From a conservation perspective, such evaluation is critical (Frankham et al., 2002), but was not addressed in the previous studies cited above, primarily due to the insufficient sample sizes for robust estimation of genetic diversity indices. To improve sample representation for each lineage, we combined our sequence data with those from Sekino and Yamashita (2013, 2016).

Third, we inferred the demographic history of each lineage. In particular, we examined whether the Last Glacial Maximum (LGM; 30 to 19 ka, centered around 21 ka; Clark and Mix, 2002)—a period associated with habitat shifts and demographic changes in many taxa (Hewitt, 2000)—also influenced the demography of Saccostrea lineages. We hypothesized that comparisons of demographic histories among multiple sympatric lineages would enable a robust assessment of the potential impacts of the LGM.

In the present study, we focus on non-mordax lineages (see Zhang et al., 2025 for an assessment of genetic diversity in S. mordax lineages). For brevity, alphabetical non-mordax lineages are hereafter referred to as “Lineage” (e.g., Lineage A) following Snow et al. (2023), and S. mordax lineages, as “mordax” (e.g., mordax A). Additionally, our study islands (Fig. 1 and Tab. 1) are expressed without the term “Island” (e.g., Ishigaki for Ishigaki Island).

Thumbnail: Fig. 1 Refer to the following caption and surrounding text. Fig. 1

Geographic map showing sampling locations. Sampling sites examined in the present study are indicated by italicized names in red. Other sites were surveyed by Sekino and Yamashita (2013, 2016). The warm, northward-flowing Kuroshio Current passes between Taiwan and Ishigaki/Iriomote and flows along the Pacific coast of Kyushu, Shikoku, and Honshu Islands, Japan. OT and YT denote the Okinawa Trough and the Yonaguni Depression, respectively. The map was produced using base maps processed with Ocean Data View version 5.6.5 (developed by Reiner Schlitzer; https://odv.awi.de/; last accessed on 6 August, 2025).
Table 1

Sampling sites of specimens analyzed in the present study.

2 Materials and methods

2.1 Sampling of oysters

Our sampling in Ishigaki and Iriomote avoided areas where the collection of marine organisms is prohibited by the Japanese Ministry of the Environment. The specific regulations and protected zones are provided on the official website for Iriomote-Ishigaki National Park (https://www.env.go.jp/en/nature/nps/park/iriomote/index.html; last accessed on 7 August, 2025).

We collected Saccostrea-like oysters by hand from eight sites (Fig. 1 and Tab. 1). The number of specimens per site ranged from 5 to 89 (total N = 250). A small number of individuals morphologically identified as S. mordax were included as an outgroup in phylogenetic analyses. Our study sites IRY (the mouth of the Yonata River, Iriomote) and ISN (Nagura Bay, Ishigaki), as well as site IRS (Iriomote) previously surveyed by Sekino and Yamashita (2016), were characterized by estuarine environments featuring mangrove forests. At all sites excluding IRY and ISN, we frequently encountered exceptionally large oysters. Based on their extremely large size and distinctive shell morphology, we hypothesized that those oysters belonged to S. spathulata. Accordingly, we preferentially sampled large oysters at the six sites. This biased sampling was beneficial for increasing the sample size of a particular lineage (presumably S. spathulata) for analyses of mitochondrial diversity and historical demography, while compromising the representativeness of the observed lineage frequencies at each site. Therefore, comparisons of lineage diversity among sites should be interpreted with caution.

Shell height (SH), defined as the maximum shell dimension, and wet body weight (BW), including both valves and the soft tissues, were recorded for all specimens except those severely damaged during sampling. Specimens were subsequently stored at −20°C. For DNA extraction, a small piece of frozen mantle or adductor muscle was excised with a clean scalpel blade. Each tissue sample was then incubated in 300 μl of TNES buffer containing 8 M urea and proteinase K (Asahida et al., 1996) at 37°C for approximately six days. Total genomic DNA was extracted from the lysate using Maxwell RSC Blood DNA Kit combined with Maxwell RSC Instrument (Promega).

2.2 DNA sequences of mitochondrial genes

Sekino and Yamashita (2016) investigated the phylogenetic diversity of their Saccostrea specimens based on the nucleotide sequences of four mitochondrial genes: the 3′-end region of split 16S ribosomal RNA (Milbury and Gaffney, 2005) (hereafter 16S), cytochrome c oxidase subunit I (COX1), cytochrome c oxidase subunit III (COX3), and NADH dehydrogenase subunit 3 (NAD3). For consistency with their dataset, we amplified the same genes from our specimens by polymerase chain reaction (PCR). PCR primer sets for amplification have been described elsewhere: partial 16S and COX3 in Sekino and Yamashita (2016), partial COX1 in Folmer et al. (1994), and full-length NAD3 in Sekino et al. (2012). However, the universal COX1 primers yielded poor amplification for most specimens from site IRY for unknown reasons (specifically, those with the prefix MA24; Tab. 1). For specimens with repeated COX1 amplification failures, we designed and used novel degenerate primers (SacCOXI-BF, 5′-CCTGTTYTRCTTTAYGGKAARCGRGAG-3′; and SacCOXI-BR, 5′-CACAAAACATGRGAAATTATYCCAAA-3′; K = G/T, R = A/G, and Y = C/T). We manually designed these primers based on alignments of full-length COX1 sequences derived from complete mitogenome sequences of seven nominal Saccostrea species available in public databases (DDBJ/EMBL/GenBank; Supplementary Material 2) as follows: S. cuccullata (KT992045), S. echinata (NC_036478), S. glomerata (NC_036483), S. kegaki (KT936590), S. malabonensis (ON649706; Mu, 2022), S. mytiloides (NC_036479), and S. mordax (FJ841968). However, species assignments to these mitogenomes are uncertain, except for S. mordax (see Results). Based on these mitogenomes, the expected PCR amplicon size for this primer set is 779 bp.

PCR amplification was performed using TaKaRa Ex Taq polymerase (Takara) or KOD FX Neo polymerase (Toyobo) in Life Touch or Life ECO thermal cyclers (Biore) (detailed in Supplementary Material 3). PCR amplicons were enzymatically purified by simultaneous treatment with Thermosensitive Alkaline Phosphatase (Promega) and Exonuclease I (New England Biolabs), following Sekino and Yamashita (2013). The purified amplicons were sequenced from both directions using BigDye Terminator v3.1 Cycle Sequencing Kit and a 3730xl DNA Analyzer (Applied Biosystems by Thermo Fisher Scientific). For poorly sequenced specimens, PCR products were subjected to agarose gel electrophoresis, and the corresponding amplicon bands were excised and purified prior to sequencing (MagExtractor–PCR & Gel Clean Up; Toyobo).

Forward and reverse sequences were trimmed and aligned using DNASIS Pro version 2.02 (Hitachi Software Engineering) to ensure consistency of sequence endpoints with those reported by Sekino and Yamashita (2016). Unlike COX1, COX3, and NAD3, the length of the 16S sequences varied among the specimens due to insertion/deletion (indel) polymorphisms. Therefore, multiple sequence alignment of our 16S sequences and 36 reference sequences (Supplementary Material 2) was performed using the E-INS-i algorithm implemented in MAFFT version 7.49 (Katoh and Standley, 2013). Ambiguously aligned regions, gap positions, and sites containing missing nucleotides in one or more sequences were removed with Gblocks version 0.91b (Castresana, 2000).

2.3 Phylogenetic tree reconstruction

We reconstructed phylogenetic trees to assign our specimens to predefined lineages using two approaches: neighbor-joining (NJ; Saitou and Nei, 1987) and Bayesian inference (BI; Huelsenbeck and Ronquist, 2001). Both methods were applied to individual genes as well as concatenated sequences of the four genes. Reference sequences available in public databases (Supplementary Material 2), including those from the complete mitogenomes described above, were incorporated into the tree reconstructions. For concatenated sequences, we re-performed multiple alignment on a reduced dataset for which fewer reference sequences were available (20 sequences) compared with the 16S-based analyses, and the 16S block was subsequently cleaned.

The NJ calculation with Kimura’s two-parameter model (K2P; Kimura, 1980) was conducted in MEGA version 6.06 (Tamura et al., 2013). Node reliability was evaluated by 1,000 bootstrap replicates. We used MrBayes version 3.2.7 (Ronquist et al., 2012) for BI analyses. Prior to BI computations, the best-fit evolutionary model for each gene was selected using jModelTest version 2.1.10 (Darriba et al., 2012) under default setting (11 substitution schemes, 88 models). The selected model was consistent across all four genes according to the Bayesian information criterion (see Results). For BI analyses, we conducted Markov chain Monte Carlo (MCMC) simulations (four chains per MCMC run in two parallel runs) until the average standard deviation of split frequencies fell below 0.01 (from 1.2 × 107 to 3.0 × 107 generations depending on the gene). Trees were sampled every 500 generations, and the first 25% of the samples were discarded as burn-in.

2.4 Inter- and intra-lineage mitochondrial diversity

In the statistical procedures described below, we used our sequences and those of known lineage identity from Sekino and Yamashita (2013, 2016) (Tab. 1; accession numbers are provided in Supplementary Material S1 of Sekino and Yamashita 2016). Unless otherwise stated, we analyzed concatenated sequences of the four genes after multiple alignment and 16S-trimming for this combined dataset.

As a measure of the nucleotide divergence between lineages, we estimated net nucleotide divergence between populations (Da; Nei, 1987) using DnaSP version 6 (Rozas et al., 2017). We performed the following analyses for each lineage to evaluate mitochondrial diversity within lineages. A statistical parsimony network among haplotypes (TCS method; Templeton et al., 1992) was constructed with PopART version 1.7 (Leigh and Bryant, 2015). Standard indices of mitochondrial diversity and two neutrality statistics—Tajima’s D (Tajima, 1989) and Fu’s FS (Fu, 1997)—were estimated using DnaSP. With the K2P distance model, pairwise FST (ΦST) between population samples was calculated in Arlequin version 3.5.2 (Excoffier et al., 2005), and significant departure of FST from zero was evaluated via 10,000 permutations of haplotypes. For FST estimation, we pooled specimens from multiple sampling sites within an island into an “island population sample” to increase sample size. An exception was made for the Amami population, which included neighboring Amami-Oshima and Kakeroma (Fig. 1). Specimens from Okinawa (sites A-Sa, G-Sa, and K-Sa; Fig. 1 and Tab. 1) were excluded due to small sample sizes (six in Lineage C and seven in Lineage F; Sekino and Yamashita, 2013, 2016). Accordingly, three population samples were defined: Ishigaki (ISF + ISI + ISN + ISU), Iriomote (IRF + IRF2 + IRS + IRSn + IRU + IRY), and Amami (AMS + AMK). Moreover, population samples with fewer than 10 specimens within a lineage were omitted from FST calculations. This resulted in 10 pairs of population samples available for comparisons across lineages and “sublineages” (see Results for the definition of sublineages).

2.5 Demographic history

Using BEAST version 2.6.6 (Bouckaert et al., 2014), we estimated historical changes in effective population size (female effective size; Nef) for each lineage based on Bayesian skyline plot (BSP) analysis (Drummond et al., 2005). For BSP computations, we set the best-fit evolutionary model identified using jModelTest. Ren et al. (2010) proposed that the Pacific oyster Crassostrea gigas (Thunberg, 1793) and the closely related Portuguese oyster Crassostrea angulata (Lamarck, 1819) diverged approximately 2.7 Mya based on nucleotide divergence between their mitogenomes. The two mitogenomes used in Ren et al. (gigas: EU672831; angulata: NC_012648) showed a nucleotide divergence of 2.26% in concatenated sequences of the four genes corresponding to those analyzed in the present study. Given these estimates, we applied a per-lineage molecular clock rate of 4.2 × 10−6 substitutions per site per kyr under a strict molecular clock model. MCMC simulations were run for 2.0 × 108 generations with an initial burn-in of 2.0 × 107. The resulting log files and trees were sampled every 3,000 steps.

3 Results

3.1 Lineage assignment

All sequences obtained in the present study have been deposited in public databases (see Supplementary Material 4 for accession numbers). The partial COX1, COX3, and full-length NAD3 sequences comprised 561, 765, and 354 bases, respectively, with no indel polymorphisms among our 250 specimens and reference sequences (COX1, 31 sequences; COX3, 21; and NAD3, 20; Supplementary Material 2). However, among the seven mitogenomes, the S. cuccullata mitogenome contained an anomalously short NAD3 gene (315 bases), most likely due to a sequencing or assembly error. Consequently, we excluded this mitogenome from phylogenetic analyses based on NAD3 and concatenated sequences.

The length of partial 16S sequences varied from 460 to 499 when all 36 reference sequences were included. Three of these reference sequences (Lineage K; Supplementary Material 2) were notably shorter, lacking more than 20 bases at the 3'-end. In contrast, the lengths of our sequences and the remaining reference sequences ranged from 485 to 499 bases. After multiple alignment and subsequent base-pruning, the final length of 16S sequences was reduced to 422 bases from an initial aligned length of 512 bases.

We used 20 reference sequences to analyze concatenated sequences (each derived from a single individual). Because fewer reference sequences were available, the length of the aligned and trimmed 16S block in the concatenated dataset increased to 463 bases, contributing to a total concatenated sequence length of 2,143 bases.

The Hasegawa-Kishino-Yano model (HKY; Hasegawa et al., 1985), incorporating a proportion of invariable sites (I; Churchill et al., 1992) and gamma-distributed rate heterogeneity (G; Yang, 1994), was consistently selected as the best-fit evolutionary model across all four genes. Accordingly, we uniformly applied the model HKY+I+G for BI tree reconstruction. The phylogenetic tree based on concatenated sequences (Fig. 2 for the BI tree; see Supplementary Material 5 for the NJ tree) revealed that all but one of our specimens were assigned to one of 11 predefined lineages (Lineages B‒D, F‒I, S. spathulata, mordax A–C; Tab. 2). A specimen from Iriomote (site IRY; ID: MA23-25) clustered near the Lineage H and S. kegaki clades but was distinct from both. We therefore designated this specimen as a novel lineage, “L”. Overall, clustering patterns were consistent across genes and tree reconstruction methods, although relationships among clades often varied (see Supplementary Materials 6 and 7 for single-gene BI and NJ trees, respectively).

The trees based on concatenated sequences (Fig. 2; Supplementary Material 5) also revealed several distinct clusters within Lineages C, D, and S. spathulata: two within C (C1 and C2), three within D (D1–D3), and two within S. spathulata (J1 and J2). These clusters were supported by high credibility values; however, their nucleotide divergence (Da) was very low compared with inter-lineage Da values (see below). Hence, these clusters are hereafter referred to as “sublineages.” The separation of these sublineages was not always evident in single-gene trees (Supplementary Materials 6 and 7). The occurrence of sublineages was not restricted to any particular island (Tab. 2; sublineages in specimens from Sekino and Yamashita, 2013, 2016, were assigned based on the following haplotype network construction).

Among our specimens, Lineage F was the most frequent (N = 74), followed by S. spathulata (N = 67) (Tab. 2). The high frequency of S. spathulata was directly attributable to our biased sampling strategy. All of nine observed lineages (excluding three mordax lineages) were identified among specimens from Iriomote, whereas Lineages B and L were absent from Ishigaki. Although Lineage K (i.e., S. malabonensis lineage 2) was previously reported from both Ishigaki and Iriomote (Hamaguchi et al., 2014), none of our specimens were placed in the Lineage K clade (16S-based trees; Supplementary Materials 6 and 7). Additionally, three sites—IRY and IRS (Iriomote) and ISN (Ishigaki)—characterized by mangrove estuaries exhibited relatively high lineage diversity (IRY, eight lineages; IRS, six; and ISN, seven; Tab. 2).

The complete mitogenomes other than those of S. mordax and S. malabonensis were assigned to Lineage B (Fig. 2; Supplementary Materials 5–7). This inconsistency likely reflects incorrect species identification of the source specimens used for these mitogenomes. Consequently, phylogenetic relationships among Saccostrea “species” inferred from these mitogenomes (Mu, 2022) should be re-evaluated with caution.

Thumbnail: Fig. 2 Refer to the following caption and surrounding text. Fig. 2

Bayesian phylogenetic tree of Saccostrea specimens based on concatenated sequences of four mitochondrial genes. The mordax group (mordax A, B, and C) was used as the outgroup. Lineages detected in the present study are indicated by bold, enlarged, and italicized names. Posterior probability values (≥ 0.80) are shown for nodes. Lineages were identified by cross-referencing single-gene trees (Supplementary Materials 6 and 7).
Table 2

Frequency of each lineage at each sampling site.

3.2 Shell measurement

Shells of at least one specimen from each observed lineage have been deposited in The University Museum, The University of Tokyo (Supplementary Material 4). Many specimens were not suitable for shell measurements, particularly for BW estimation, because their valves firmly attached to substrates or neighboring individuals were often severely damaged during sampling. Consequently, we measured SH and BW in 230 and 171 specimens, respectively. All specimens with large shell size (SH ≥ 100 mm) belonged exclusively to S. spathulata (Supplementary Material 8). The largest S. spathulata specimen measured 181.2 mm SH and 809.0 g BW. The distributions of both SH and BW largely overlapped among the remaining lineages. Among the identified lineages other than S. spathulata, the largest specimen was from Lineage F (79.0 mm SH). The specimen of the novel Lineage L was notably small (29.6 mm SH), smaller than the minimum SH observed in any other lineage. Some shells of the lineages are shown in Figure 3 (additional images of S. spathulata shells are provided in Supplementary Material 9). As discussed below, no definitive morphological keys were identified to distinguish among lineages, except for S. spathulata.

Thumbnail: Fig. 3 Refer to the following caption and surrounding text. Fig. 3

Shells of observed lineages. Each specimen is shown in up to four views and is identified by a unique label designated as follows: lineage or sublineage name, followed by the specimen ID and shell height (SH; mm, in parentheses). Abbreviations used in the figure are as follows: BV, both valves; OLV, outer left valve; ILV, inner left valve; ORV, outer right valve; IRV, inner right valve; GL, growth lamellae; HS, hyote spine.

3.3 Intra-lineage mitochondrial diversity

After excluding 11 specimens from the three mordax lineages (Tab. 2), we analyzed a combined dataset of our 239 sequences and 201 from Sekino and Yamashita (2013, 2016). In this dataset, concatenated sequences were 2,146 bases in length, with the cleaned 16S block contributing 466 bases.

We constructed a haplotype network for each lineage (Fig. 4), except for Lineages B (N = 2) and L (N = 1). The networks supported the validity of our sublineage definitions for Lineages C, D, and S. spathulata; these sublineages were separated by 10 or more mutational steps in all three lineages, when accounting for presumably unsampled or extirpated haplotypes. The haplotype networks exhibited a highly scattered pattern with numerous singleton haplotypes across lineages, consistent with a high level of mitochondrial diversity. Notably, haplotype sharing was absent or extremely limited within Lineages G, I, and S. kegaki, which consisted of divergent singletons separated by several mutational steps. The elevated mitochondrial diversity in the three lineages also became evident by estimating the observed number of haplotypes (NH) and unbiased haplotype diversity (H; Nei, 1987) (Tab. 3; see Supplementary Material 10 for gene-specific estimates). In these lineages, NH was equal to or very close to the respective sample size, resulting in saturated estimates of H (0.999–1.000).

Significantly negative values for either or both neutrality statistics (D and FS) were obtained across all lineages and sublineages (Tab. 3), indicating an excess of rare mutations. Neutrality statistic based on synonymous sites (DSyn) was largely consistent with D. Neutrality tests based on individual genes also yielded significantly negative values for D, FS, or both in most lineages and sublineages (Supplementary Material 10). In the 16S gene, slightly positive values were observed for D (0.65) and FS (0.87) in Lineage H, as well as for FS (0.31) in Lineage D (P > 0.1). Additionally, Sublineages C1 and D1 scored negative but non-significant D and FS values in the 16S gene (C1: D = −1.03, FS = −2.29; D1: D = −1.72, FS = −2.68; P > 0.05).

We subsequently estimated FST values between population samples within lineages and sublineages. When adhering to our definition of population samples and the sample-size criterion, 10 pairs of samples from five lineages (C, D, F, I, and S. spathulata) and three sublineages (C2, J1, and J2) were available for FST estimation (Tab. 4). Overall, FST values were low, ranging from <0.000 to 0.046. Only one pair yielded an FST value significantly greater than zero (Ishigaki vs. Iriomote in Lineage I; FST = 0.043), even after correction for multiple comparisons using the false discovery rate (Benjamini and Hochberg, 1995; q = 0.040)

Thumbnail: Fig. 4 Refer to the following caption and surrounding text. Fig. 4

Statistical parsimony (TCS) networks of haplotypes in each Saccostrea lineage based on concatenated sequences of four mitochondrial genes. The number of tick marks on each edge represents the number of mutational steps between nodes (no tick mark is placed for single step mutations). Node size is proportional to haplotype frequency. Lineages B and L with small sample sizes (N = 2 and 1, respectively) are omitted.
Table 3

Mitochondrial diversity within lineages and sublineages based on concatenated sequences of the four mitochondrial genes (2,146 bases).

Table 4

FST (ΦST) estimated between population samples within each lineage and sublineage.

3.4 Inter-lineage nucleotide divergence

We estimated Da values among lineages and sublineages (the Da matrix is provided in Supplementary Material 11). Inter-lineage values ranged from 4.3% (Lineage H vs. S. kegaki) to 16.4% (Lineage G vs. S. spathulata). The nucleotide divergence of the novel Lineage L from the closely related Lineage H and S. kegaki was 5.2% and 6.0%, respectively. Regarding Lineage C and S. spathulata, comparisons between the respective sublineages yielded Da values of less than 1.0% (0.5% for C1 vs. C2; 0.6% for J1 vs. J2). These values were comparable to the Da value between mordax A and B (0.7%; not shown in the Da matrix). Inter-sublineage Da values within Lineage D were slightly higher (0.9% for D1 vs. D2; 1.3% for D1 vs. D3; and 1.1% for D2 vs. D3). All of these five inter-sublineage values fell within the reported range of mean intraspecific nucleotide diversity in oysters (0.0–1.4% based on COX1 sequences; Liu et al., 2011).

Based on the Da values combined with the per-lineage molecular clock rate defined above, we estimated divergence times (Td) among all identified lineages and sublineages (Supplementary Material 11). Inter-lineage Td values ranged from 5.1 Myr (Lineage H vs. S. kegaki) to 19.5 Myr (Lineage G vs. S. spathulata). Most lineages were estimated to have diverged more than 10 Mya, except for Lineages H, L, and S. kegaki, which likely diverged more recently (Td: 5.1–7.1 Myr). Inter-sublineage divergence times were much lower, with Td values ranging from 0.6 Myr (C1 vs. C2) to 1.5 Myr (D1 vs. D3).

3.5 Demographic history

Given generally low estimates of FST, we considered population structuring within lineages to be negligible. Accordingly, we pooled specimens from multiple islands within each lineage for BSP analyses. Using the HKY+I+G model, we estimated BSPs for eight lineages with sample sizes ≥ 16 (C, D, F, G, H, I, S. kegaki, and S. spathulata). Additionally, we separately conducted BSP analyses for sublineages within Lineages C (C1 and C2), S. spathulata (J1 and J2), and D (D1), but not for Sublineages D2 and D3 (N ≤ 6).

BSP analyses provided little evidence for drastic demographic changes around or during the LGM (Fig. 5). Lineage D exhibited a rapid increase in Nef at approximately 150 ka, although the BSP for its sublineage (D1) remained relatively constant. A similar increase in Nef was observed in Lineage F around the same time period. In Lineages G, I, and S. kegaki, demographic expansion was inferred to have occurred earlier (>350 ka). Both Lineage C and S. spathulata exhibited a stepwise increase in Nef, a pattern not evident in their sublineages (C1, C2, J1, and J2).

Thumbnail: Fig. 5 Refer to the following caption and surrounding text. Fig. 5

Bayesian skyline plots (BSPs) estimated for Saccostrea lineages and sublineages. Shaded areas represent the 95% highest posterior density interval. BSP analyses were not performed for the following lineages and sublineages owing to small sample sizes: Lineages B (N = 2) and L (N = 1), and two sublineages within Lineage D (D2, N = 6; and D3, N = 3).

4 Discussion

4.1 Identified lineages and sublineages

Among our specimens from Ishigaki and Iriomote, we identified nine lineages: Lineages B‒D, F‒I, L, and S. spathulata. The present study represents the first record of Lineages B and D from Japanese waters, while Lineage L is described here for the first time. The consistent separation of Lineage L from the other lineages, along with its close yet clearly differentiated relationships with Lineage H and S. kegaki across multiple phylogenetic trees, strongly indicates that this finding is attributable to neither sequencing error nor the presence of nuclear mitochondrial pseudogenes (Song et al., 2008; Buhay, 2009).

Additionally, sublineages were identified within Lineages C, D, and S. spathulata. We defined these sublineages as being analogous to shallow mitochondrial lineages frequently observed within species. Although the mechanisms driving the shallow inter-sublineage divergence remain unclear, geographical isolation during recent paleoceanographic transitions is considered the most plausible explanation (discussed later). These sublineages subsequently experienced secondary contact, as evidenced by their co-occurrence across all islands where the parental lineages were found. While these sublineages were identified based on concatenated sequences, their separation was not consistently recovered in single-gene analyses. The present study therefore highlights the advantage of analyzing longer mitochondrial sequences for detecting weak phylogenetic signals. In this context, a fine-scale phylogenetic portrait of Saccostrea oysters may be further elucidated through analyses of mitogenomes (Keis et al., 2013) and/or multiple loci (Dupuis et al., 2012).

4.2 A hotspot for phylogenetic diversity

In the present study, three IWP lineages previously reported from the Ryukyu Island Arc were not detected: Lineage A (Okinawa; Lam and Morton, 2006), Lineage K (Ishigaki and Iriomote; Hamaguchi et al., 2014), and S. kegaki (Iriomote and Kakeroma; Sekino and Yamashita (2016)). Combining our data with these previous records, we conclude that at least 12 lineages—namely A‒D, F‒I, K, L, S. kegaki, and S. spathulata—inhabit the Ryukyu Island Arc.

Among the 14 explicitly defined IWP lineages (A‒I, K, L, S. glomerata, S. kegaki, and S. spathulata), S. glomerata is generally considered restricted to Oceanian coasts, including Australia and New Zealand (Lam and Morton, 2006; Huber, 2010). Except for S. glomerata, Lineage E (Lam and Morton, 2006) is the only lineage not yet reported from the Ryukyu Island Arc. Lineage E occurs along the coasts of northern Taiwan and northern Australia, as inferred from phylogenetic relationships (Fig. 4 of Lam and Morton, 2006). Given this broad distribution and the geographical proximity between Taiwan and the southern Ryukyu Islands (Ishigaki and Iriomote), the absence of records for Lineage E from the Ryukyu Island Arc to date is more likely attributable to insufficient geographical sampling or limited sample sizes, notwithstanding our extended survey, rather than its true absence from the region.

The phylogenetic diversity observed in the Ryukyu Island Arc is exceptionally high compared with other IWP regions. For example, Lam and Morton (2006) identified five lineages among their specimens from northern Taiwan (B, D‒F, and G), three from southern China (B, D, and F), and two from southern Malaysia (B and F). According to Wu et al. (2019), S. kegaki also occurs along the Taiwanese coast, increasing the total to at least six lineages confirmed in Taiwanese waters. In Australia, eight lineages (A, B, E‒G, I, S. glomerata, and S. spathulata) have been recognized (Lam and Morton, 2006; Snow et al., 2023, McDougall et al., 2024). Although more extensive field surveys will surely uncover additional lineages in these regions, the identification of 12 lineages from the Ryukyu Island Arc thus far underscores its significance as a global hotspot for phylogenetic diversity. Notably, as many as 11 lineages (B‒D, F‒I, K, L, S. kegaki, and S. spathulata)—representing nearly 80% of all known IWP lineages—have been recorded from Iriomote alone, despite its small size (area, 289.6 km2, and coastline, 129 km; Statistics Bureau of Japan: https://www.stat.go.jp/english/data/nenkan/74nenkan/index.html; last accessed on 24 April, 2025). The natural environment of Iriomote has been protected by strict environmental regulations since its designation as a national park by the Japanese Government in 1972 (https://www.env.go.jp/en/nature/nps/park/iriomote/index.html; last accessed on 13 June, 2025). Furthermore, Iriomote has been registered as a Natural World Heritage Site under the UNESCO World Heritage Convention since 2021, owing to its diverse ecosystems and unique biota (https://whc.unesco.org/en/natural-world-heritage; last accessed on 13 June, 2025). It is therefore plausible that the pronounced phylogenetic diversity observed on Iriomote has been preserved, at least in part, by limited anthropogenic habitat degradation. Finally, samples collected from mangrove estuaries consistently harbor high lineage diversity, suggesting that well-preserved mangrove ecosystems play a critical role in maintaining this rich phylogenetic diversity.

4.3 Elusive shell morphology

Sekino and Yamashita (2016) and McDougall et al. (2024) referred to Lineage J as S. spathulata, distinguishing it from the other lineages based on its exceptionally large size and the distinct shell morphology as described in Introduction. Our S. spathulata specimens exhibited shell characteristics consistent with their observations. Notably, the largest specimen in the present study attained 181 mm SH, substantially exceeding previously reported maximum sizes, including that of the holotype (145 mm SH; McDougall et al., 2024). While well-grown S. spathulata individuals can be readily identified morphologically, McDougall et al. (2024) noted that the shell features of smaller individuals overlap with those of Lineages B, G, and S. glomerata, highlighting the challenge of accurately identifying S. spathulata juveniles without molecular data.

Apart from S. spathulata, few lineage-diagnostic morphological keys were available for the remaining lineages, as previously discussed (Lam and Morton, 2006; Sekino and Yamashita 2016; McDougall et al., 2024). Below, we provide brief morphological notes for several observed lineages with reference to Figure 3.

Lineage F: This lineage is the most common in the Ryukyu Island Arc and has a wide geographical distribution in the West Pacific, including, but very likely not limited to, Malaysia, Taiwan, Hainan Island (southern China), Honshu and Shikoku Islands (Japan), Ryukyu Island Arc, and Australia (Lam and Morton, 2006; Sekino and Yamashita, 2013, 2016; Hamaguchi et al., 2014; McDougall et al., 2024; this study). The highly variable shell forms—ribbed, deeply cupped, elongated, or flat—make morphological identification of this lineage particularly challenging.

Lineage G: Some specimens possessed an extended, paper-thin shell margin (G_MA24-13), resembling Lineage I (see below). However, other specimens exhibited a wavy or saw-like shell margin and ribbed left valve (G_MA23-07 and G_MA23-08), posing difficulties in differentiating this lineage from Lineage C (C1_MA23-09 and C2_MA23-17).

Lineage H: Sekino and Yamashita (2016) noted that small individuals of this lineage often exhibited hyote (tubular) spines, an important taxonomic key for morphologically identifying S. kegaki. Their observation was supported by the present study; of our two Lineage H specimens, one possessed noticeable spines while the other did not (H_SH23-02 and H_Nm24-37, respectively). Spines have also been observed in some specimens of Lineages A and B (Lam and Morton, 2006; McDougall et al., 2024). The phylogenetically close relationships among Lineages A, B, H, and S. kegaki hint that the trait is underlain by genetic factors; however, the mechanisms underlying individual variation in spine expression, even under apparently similar environmental conditions (Sekino and Yamashita 2016), remain unclear.

Lineage I: This lineage was frequently found on mangrove trees, consistent with the observations by Sekino and Yamashita (2016) and McDougall et al. (2024). Sekino and Yamashita (2016) reported a remarkable extension of growth lamellae that appear to mimic tree bark. McDougall et al. (2014) characterized this lineage by broad, lobe-like squamae. The shells of our specimens generally align with these descriptions (I_MA23-21), but the presence of such an extended shell margin was variable (I_Nm24-22). Interestingly, Lineage I exhibited a very small but significantly positive FST value between the Iriomote and Ishigaki population samples, whereas FST estimates were effectively zero in all pairwise comparisons within the broadly distributed Lineage F. This suggests that Lineage I may experience impeded inter-island connectivity as a habitat specialist favoring patchy mangrove environments.

Lineage L: Our single specimen from this lineage (L_MA23-25) was the smallest among the specimens examined. The left valve was moderately cupped, and the right valve featured a serrated edge. Weak but recognizable chomata was present, as is typical for Saccostrea oysters (but see Hamaguchi et al., 2014). Despite its close relationship with Lineage H and S. kegaki, no spines were observed, although the presence/absence of spines may be individually variable as with Lineage H.

4.4 Did the evolution of the Ryukyu Island Arc contribute to lineage diversification?

The Ryukyu Island Arc is often referred to as “the Galapagos of the Orient,” as its constituent islands harbor many endemic terrestrial animals that evolved from ancestors originating on the Eurasian (East Asian) Continent (Watanabe et al., 2023). Rifting of the Ryukyu Island Arc from the eastern continental margin commenced during the Late Miocene (10‒6 Ma). Subsequent back-arc spreading of the Okinawa Trough (Fig. 1) led to isolation of the arc from the East China Sea (ECS) Shelf during the Early Pleistocene (2.6‒0.8 Ma) (Kinoshita et al., 2025). An alternative hypothesis suggests a more recent separation at approximately 1.5 Ma (Osozawa et al., 2012); however, this timeline is inconsistent with the high level of terrestrial endemism in the Ryukyu Island Arc, which is largely attributed to pre-Pleistocene vicariant events (Kinoshita et al., 2025).

Geographical isolation of the Ryukyu Island Arc is a well-established driver of endemic terrestrial diversification in the region (Pfingstl et al., 2019; Watanabe et al., 2023). Whether this model also explains the extensive phylogenetic diversity of Saccostrea oysters in the region remains an open question. Given the estimated divergence times among the phylogenetically close Lineages H, L, and S. kegaki (5.1–7.1 Myr), together with the observations of Lineages H and L only from the Ryukyu Island Arc (H in Iriomote and Ishigaki; L in Iriomote), it is possible that the divergence among these three took place in concert with the geographical evolution of the Ryukyu Island Arc. Furthermore, the divergence between sublineages within Lineages C, D, and S. spathulata dates to the Early to Middle Pleistocene (1.5–0.6 Mya), suggesting that the evolution of the Ryukyu Island Arc may have contributed, at least in part, to this diversification. However, we cannot exclude the alternative possibility that these lineages and sublineages colonized the Ryukyu Island Arc after diverging in remote regions. Unlike terrestrial animals, oysters can disperse over long distances via planctonic larvae transported by major current systems such as the Kuroshio Current (Fig. 1). Testing these hypotheses will require broader geographical surveys beyond the Ryukyu Island Arc to evaluate the degree of endemism of these lineages and sublineages.

4.5 Limited influence of the Last Glacial Maximum on historical demography

Glacial-interglacial cycles over the past 800 kyr have caused repeated sea-level fluctuations (Lambeck et al., 2002). During the LGM, global sea levels experienced a pronounced regression (Clark and Mix, 2002; Lambeck et al., 2002). Multiple lines of palaeoceanographic evidence indicate that, during the LGM, the sea level in the ECS fell by 120‒160 m relative to present-day levels (Fig. 1). This regression resulted in subaerial exposure of the vast majority of the ECS Shelf (Emery et al., 1971; Wang and Wang, 1980; Saito, 1998; Xiao et al., 2004; Matsuzaki et al., 2019).

The population dynamics of marine organisms in the West Pacific, including the ECS, have been strongly influenced by palaeoceanographic oscillations stemming from eustatic climatic changes (Liu et al., 2006; Song et al., 2010; Ni et al., 2014). Such species typically exhibit genetic signatures of pre-LGM demographic expansion (reviewed in Ni et al., 2014), although some taxa underwent population expansion during the LGM (Okinawa-lineage of Eriocheir crab; Xu et al., 2009) or after the LGM (the Portuguese oyster in Taiwan; Hsiao et al., 2016). The non-star-like haplotype networks and high haplotype diversity observed in most lineages are indicative of long-term population stability. Consistent with this pattern, our BSP analyses showed that the onset of demographic expansion in the tested lineages—excluding Lineage H with a flat BSP—predated the LGM, a result comparable to patterns reported for various West Pacific marine species (Ni et al., 2014). Population sizes also appeared to remain largely stable throughout and after the LGM. While Lineage C and S. spathulata, each of which included two sublineages (C1 and C2 within C; J1 and J2 within S. spathulata), exhibited stepped BSP profiles with a relatively recent demographic expansion (the second expansion around 100 ka), this finding may be an analytical artifact. Step-formed BSP curves can arise if genetically divergent populations are pooled for reconstruction, potentially producing spurious demographic signals (Fagundes et al., 2008). Based on FST estimates, it is very unlikely that our population samples of Lineage C and S. spathulata contained specimens from divergent populations. Nevertheless, the observed stepped BSPs could be a spurious result of pooling specimens from sublineages that evolved in allopatry and were recently admixed, complicating interpretation of their demographic histories. Similarly, the apparent rapid expansion inferred for Lineage D around 150 ka should be interpreted with caution, as the dataset comprised three sublineages (Grant, 2015). Indeed, these punctuated BSP patterns disappeared when BSPs were reconstructed separately for individual sublineages, although the nearly flat BSPs in Sublineages C1, J1, and D1, as well as in Lineage H, warrant careful consideration, given the limited power to detect population expansion associated with small sample sizes (Grant, 2015).

Our BSP results, albeit with some uncertainty, along with the non-star-like haplotype networks and extensive haplotype diversity observed across lineages, suggest that the LGM had only a limited impact on the population abundance of Saccostrea oysters in the Ryukyu Island Arc, despite the extensive subaerial exposure of most of the ECS Shelf during that period. This interpretation is not entirely concordant with the significantly negative neutrality statistics detected across lineages, which are often interpreted as signatures of recent and sudden population expansion (Ramos-Onsins and Rozas, 2002). However, negative neutrality statistics can also arise from processes other than demographic change that violate the assumptions of the standard neutral model, including natural selection. In particular, background selection, whereby purifying selection removes mildly deleterious mutations, can reduce diversity at linked neutral sites and generate an excess of rare variants, thereby producing negative neutrality statistics (Charlesworth et al., 1993; Charlesworth, 2012). Given the BSP results, the presence of numerous low-frequency haplotypes that do not form star-like networks, and the consistently negative neutrality statistics—including DSyn—across most gene-lineage combinations, our neutrality test results are more plausibly explained by background selection acting on the nonrecombining mitochondrial genome in large populations with continual input of mildly deleterious variants, rather than by recent demographic expansion.

The geographical position of the Ryukyu Island Arc, which extends along the deep Okinawa Trough (Fig. 1), likely limited the impact of the LGM on demographic trajectories. With a maximum depth of approximately 2,200 m (Letouzey and Kimura, 1985), the Okinawa Trough remained a deep-water marine environment even during the LGM. This ensured that surrounding regions were never subaerially exposed, despite drastic sea-level regressions (Xu and Oda, 1999). Some hypotheses propose that the warm Kuroshio Current, which enters the ECS northwards through the Yonaguni Depression between Taiwan and the southern Ryukyu Islands (Fig. 1), was deflected east of the Ryukyu Island Arc during the LGM owing to the presence of a shallow land bridge at the depression (Xu and Oda, 1999; Kao et al., 2006; Diekmann et al., 2008). However, a growing body of evidence indicates a persistent, albeit weakened, inflow into the paleo-ECS through this passage (Kawahata and Ohshima, 2004; Ijiri et al., 2005; Lee et al., 2013; Vincent-Vogt and Mitarai, 2020). The thermal and trophic transport provided by the Kuroshio Current likely sustained coral reef ecosystems along the Ryukyu Island Arc in the paleo-ECS even during the LGM (Vincent-Vogt and Mitarai, 2020). Given this paleoceanographic context, we argue that the paleo-ECS surrounding the Okinawa Trough served as a crucial refugium for marine life benefiting from the continued inflow of warm current. We further propose that the coastal waters around the paleo-Ryukyu Island Arc provided suitable habitats for tropical and subtropical Saccostrea oysters to persist without requiring large-scale range shifts or colonization of new habitats. Under such equable environmental conditions, both population sizes and mitochondrial diversity would have remained relatively stable.

4.6 Concluding remarks

The pronounced phylogenetic diversity of IWP lineages is likely driven by the accumulation of genetic differences through mutation and stochastic drift, potentially coupled with local adaptation triggered by palaeoceanographic fluctuations (Gavrilets, 2003). Natural hybridization may have further facilitated this diversification process (Abbott et al., 2013). Although mitochondrial data alone cannot fully disentangle the complex evolutionary histories underlying this elevated diversity, future population genomic approaches—encompassing both inter- and intra-lineage comparisons—will provide deeper insights into the mechanisms shaping diversification in this group.

The present study documented a remarkably high level of phylogenetic diversity among Saccostrea oysters in the Ryukyu Island Arc. Nevertheless, the detection of only a single specimen belonging to the novel Lineage L, despite our expanded sampling efforts, suggests that the true diversity of IWP Saccostrea may still be underestimated. More extensive surveys, including regions beyond the Ryukyu Island Arc, will be essential for constructing a more comprehensive picture of IWP Saccostrea diversity. We particularly recommend prioritizing field surveys in mangrove estuaries, which are likely to harbor substantial lineage variation.

A critical caveat in DNA-based specimen identification—especially for Saccostrea—is the unreliability of species names in public DNA databases (Collins and Cruickshank, 2013). Many entries contain errors arising from misidentification by original submitters, a problem further highlighted in the present study (see also Hamaguchi et al., 2014; Sekino and Yamashita 2016). The indiscriminate use of registered species names and their associated nucleotide sequences as references therefore risks perpetuating taxonomic confusion within this genus.

Acknowledgements

This study was inspired by, and is respectfully dedicated to, the late Hiroyoshi Yamashita, whose profound knowledge of shellfish taxonomy and deep insights into oyster shell morphology laid the foundation for our research. We are grateful to Masaei Kanematsu, Miki Osoko, and Tomomi Kagiyama of Yaeyama Field Station, Fisheries Technology Institute, The Japan Fisheries Research and Education Agency, for their assistance in collecting oyster specimens on Ishigaki Island. Dr. Takenori Sasaki of The University Museum, The University of Tokyo (UMUT), supported M. Sekino with shell photography and the deposition of shell materials in UMUT. The manuscript benefited greatly from constructive feedback provided by the handling editors and three anonymous reviewers. This study was conducted with support from The Japan Fisheries Research and Education Agency through the project “The search for oysters suitable for aquaculture in tropical islands and the technical development of their aquaculture (2AE210)”.

Funding

The present study was financially supported by The Japan Fisheries Research and Education Agency (2AE210).

Conflicts of interest

The authors have no competing interests to disclose.

Data availability statement

All mitochondrial DNA sequences obtained in the present study have been deposited in the DDBJ/EMBL/GenBank databases (accession numbers can be found in Supplementary Material 4).

Supplementary Material

Supplementary Material 1. Summary of synonymy and mitochondrial lineages of Indo-West Pacific (IWP) Saccostrea oysters.

Supplementary Material 2. List of reference sequences used for phylogenetic analyses.

Supplementary Material 3. PCR amplification of mitochondrial genes.

Supplementary Material 4. Accession numbers of DNA sequences in DDBJ/EMBL/GenBank databases.

Supplementary Material 5. NJ tree based on concatenated sequences of four mitochondrial genes.

Supplementary Material 6. Single-gene BI trees.

Supplementary Material 7. Single-gene NJ trees.

Supplementary Material 8. Distributions of shell height (SH) and wet body weight (BW) across identified Saccostrea lineages.

Supplementary Material 9. Shells of S. spathulata.

Supplementary Material 10. Statistics of mitochondrial polymorphisms within lineages and sublineages estimated for each gene.

Supplementary Material 11. Nucleotide divergence and divergence time between lineages.

Access Supplementary Material

References

  • Abbott R, Albach D, Ansell S, Arntzen JW, Baird SJ, Bierne N, Boughman J, Brelsford A, Buerkle CA, Buggs R, Butlin RK, Dieckmann U, Eroukhmanoff F, Grill A, Cahan SH, Hermansen JS, Hewitt G, Hudson AG, Jiggins C, Jones J, Keller B, Marczewski T, Mallet J, Martinez-Rodriguez P, Möst M, Mullen S, Nichols R, Nolte AW, Parisod C, Pfennig K, Rice AM, Ritchie MG, Seifert B, Smadja CM, Stelkens R, Szymura JM, Väinölä R, Wolf JB, Zinner D. 2013. Hybridization and speciation. J Evol Bio 26: 229–246. [Google Scholar]
  • Anderson TJ, Adlard RD. 1994. Nucleotide sequence of a rDNA internal transcribed spacer supports synonymy of Saccostrea commercialis and S. glomerata. J Molluscan Stud 60: 196–197. [Google Scholar]
  • Asahida T, Kobayashi T, Saitoh K, Nakayama I. 1996. Tissue preservation and total DNA extraction from fish stored at ambient temperature using buffers containing high concentration of urea. Fish Sci 62: 727–730. [Google Scholar]
  • Bayne B, Ahrens M, Allen SK, Anglès d'auriac M, Backeljau T, Beninger P, Bohn R, Boudry P, Davis J, Green T, Guo X, Hedgecock D, Ibarra A, Kingsley-Smith P, Krause M, Langdon C, Lapegue S, Li C, Manahan D, Mann R, Perez-Paralle L, Powell EN, Rawson PD, Speiser D, Sanchez JL, Shumway S, Wang H. 2017. The proposed dropping of the genus Crassostrea for all Pacific cupped oysters and its replacement by a new genus Magallana: a dissenting view. J Shellfish Res 36: 545–547. [Google Scholar]
  • Bayne B, Anglès d'Auriac M, Backeljau T, Beninger P, Boudry P, Carnegie R, Davis J, Guo X, Hedgecock D, Krause M, Langdon C, Lapègue S, Manahan D, Mann R, Powell E, Shumway S. 2019. A scientific name for Pacific oysters. Aquaculture 499: 373. [Google Scholar]
  • Beck MW, Brumbaugh RD, Airoldi L, Carranza A, Coen LD, Crawford C, Defeo O, Edgar GJ, Hancock B, Kay MC, Lenihan HS, Luckenbach MW, Toropova CL, Zhang G, Guo X. 2011. Oyster reefs at risk and recommendations for conservation, restoration, and management. BioScience 61: 107–116. [CrossRef] [Google Scholar]
  • Benjamini Y, Hochberg Y. 1995. Controlling the false discovery rate: a practical and powerful approach to multiple testing. J R Stat Soc Ser B 57: 289–300. [Google Scholar]
  • Bickford D, Lohman DJ, Sodhi NS, Ng OKL, Meier R, Winker K, Krista K. Ingram KK, Das I. 2006. Cryptic species as a window on diversity and conservation. Trends Ecol Evol 22: 148–155. [Google Scholar]
  • Bouckaert R, Heled J, Kühnert D, Vaughan T, Wu CH, Xie D, Suchard MA, Rambaut A, Drummond AJ. 2014. BEAST 2: A software platform for Bayesian evolutionary analysis. PLoS Comput Biol 10: e1003537. [Google Scholar]
  • Boudry P, Heurtebise S, Collet B, Cornette F, Gérard A. 1998. Differentiation between populations of the Portuguese oyster, Crassostrea angulata (Lamark) and the Pacific oyster, Crassostrea gigas (Thunberg), revealed by mtDNA RFLP analysis. J Exp Mar Biol Ecol 226: 279–291. [CrossRef] [Google Scholar]
  • Buhay JE. 2009. “COI-like” sequences are becoming problematic in molecular systematic and DNA barcoding studies. J Crustac Biol 29: 96–110. [Google Scholar]
  • Castresana J. 2000. Selection of conserved blocks from multiple alignments for their use in phylogenetic analysis. Mol Biol Evol 17: 540–552. [Google Scholar]
  • Charlesworth B, Morgan MT, Charlesworth D. 1993. The effect of deleterious mutations on neutral molecular variation. Genetics 134: 1289–1303. [Google Scholar]
  • Charlesworth B. 2012. The effects of deleterious mutations on evolution at linked sites. Genetics 190: 5–22. [Google Scholar]
  • Churchill GA, von Haeseler A, Navidi WC. 1992. Sample size for a phylogenetic inference. Mol Biol Evol 9: 753–769. [Google Scholar]
  • Clark PU, Mix AC. 2002. Ice sheets and sea level of the Last Glacial Maximum. Quat Sci Rev 21: 1–7. [Google Scholar]
  • Collins RA, Cruickshank RH. 2013. The seven deadly sins of DNA barcoding. Mol Ecol Resour 13: 969–975. [Google Scholar]
  • Cui Z, Hu L, Li C, Zhang Z, Guo X, Wang H. 2021. Identification of Saccostrea mordax and a new species Saccostrea mordoides sp. nov. (Bivalvia: Ostreidae) from China. J Shellfish Res 40: 63–75. [Google Scholar]
  • Darriba D, Taboada GL, Doallo R, Posada D. 2012. jModelTest 2: more models, new heuristics and parallel computing. Nat Methods 9: 772. [Google Scholar]
  • Diekmann B, Hofmann J, Henrich R, Fütterer DK, Röhl U, Wei KY. 2008. Detrital sediment supply in the southern Okinawa Trough and its relation to sea-level and Kuroshio dynamics during the late Quaternary. Mar Geol 255: 83–95. [Google Scholar]
  • Dimmick WW, Ghedotti MJ, Grose MJ, Maglia AM, Meinhardt DJ, Pennock DS. 1999. The importance of systematic biology in defining units of conservation. Conserv Biol 13: 653–660. [Google Scholar]
  • Drummond AJ, Rambaut A, Shapiro B, Pybus OG. 2005. Bayesian coalescent inference of past population dynamics from molecular sequences. Mol Biol Evol 22: 1185–1192. [Google Scholar]
  • Dupuis JR, Roe AD, Sperling FAH. 2012. Multi-locus species delimitation in closely related animals and fungi: one marker is not enough. Mol Ecol 21: 4422–4436. [Google Scholar]
  • Emery KO, Niino H, Sullivan B. Post-Pleistocene levels of East China Sea, in: K.K. Trekian (Ed.), Late Cenozoic Glacial Ages, Yale University Press, New Heaven, USA, 1971, pp. 381–390. [Google Scholar]
  • Excoffier L, Laval G, Schneider S. 2005. Arlequin ver. 3.0: An integrated software package for population genetics data analysis. Evol Bioinform Online 1: 47–50. PMID: 19325852; PMCID: PMC2658868. [Google Scholar]
  • Fagundes NJ, Kanitz R, Bonatto SL. 2008. A reevaluation of the Native American mtDNA genome diversity and its bearing on the models of early colonization of Beringia. PLoS One 17: e3157. [Google Scholar]
  • Folmer O, Black M, Hoeh W, Lutz R, Vrijenhoek R. 1994. DNA primers for amplification of mitochondrial cytochrome c oxidase subunit I from diverse metazoan invertebrates. Mol Mar Biol Biotechnol 3:294–299. PMID: 7881515 [PubMed] [Google Scholar]
  • Frankham R, Ballou JD, Briscoe DA. Introduction to Conservation Genetics, Cambridge University Press, Cambridge, UK, 2002. [Google Scholar]
  • Fraser DJ, Bernatchez L. 2001. Adaptive evolutionary conservation: towards a unified concept for defining conservation units. Mol Ecol 10: 2741–2752. [Google Scholar]
  • Fu YX. 1997. Statistical tests of neutrality of mutations against population growth, hitchhiking and background selection. Genetics 147: 915–925. [CrossRef] [PubMed] [Google Scholar]
  • Gavrilets S. 2003. Perspective: Models of speciation: What have we learned in 40 years? Evolution 57: 2197–2215. [Google Scholar]
  • Gottlieb SJ, Schweighofer ME. 1996. Oysters and the Chesapeake Bay ecosystem: A case for exotic species introduction to improve environmental quality? Estuaries 19: 639–650. [Google Scholar]
  • Grant WS. 2015. Problems and cautions with sequence mismatch analysis and Bayesian skyline plots to infer historical demography. J Hered 106:333–346. [Google Scholar]
  • Guo X, Li C, Wang H, Xu Z. 2018. Diversity and evolution of living oysters. J Shellfish Res 37: 755–771. [Google Scholar]
  • Hamaguchi M, Shimabukuro H, Usuki H, Hori M. 2014. Occurrences of the Indo-West Pacific rock oyster Saccostrea cucullata in mainland Japan. Mar Biodivers Rec 7: e84. [Google Scholar]
  • Harry HW. 1985. Synopsis of the supraspecific classification of living oysters (Bivalvia: Gryphaeidae and Ostreidae). Veliger 28: 121–158. [Google Scholar]
  • Hasegawa M, Kishino H, Yano T. 1985. Dating of the human-ape splitting by a molecular clock of mitochondrial DNA. J Mol Evol 22: 166–174. [Google Scholar]
  • Hebert PDN, Cywinska A, Ball SL, deWaard JR. 2003. Biological identifications through DNA barcodes. Proc R Soc Lond Ser B, Biol Sci 270: 313–321. [Google Scholar]
  • Hebert PDN, Stoeckle MY, Zemlak TS, Francis CM. 2004. Identification of birds through DNA barcodes. PLoS Biol 2: e312. [Google Scholar]
  • Hewitt G. 2000. The genetic legacy of the Quaternary ice ages. Nature 405: 907–913. [Google Scholar]
  • Hewitt GM. 2004. Genetic consequences of climatic oscillations in the Quaternary. Philos Trans R Soc Ser B 359: 183–195. [Google Scholar]
  • Hsiao ST, Chuang SC, Chen KS, Ho PH, Wu CL, Chen CA. 2016. DNA barcoding reveals that the common cupped oyster in Taiwan is the Portuguese oyster Crassostrea angulata (Ostreoida; Ostreidae), not C. gigas. Sci Rep 6: 34057. [Google Scholar]
  • Huber M. Compendium of Bivalves, ConchBooks, Hackenheim, Germany, 2010. [Google Scholar]
  • Huelsenbeck JP, Ronquist F. 2001. MRBAYES: Bayesian inference of phylogenetic trees. Bioinformatics 17: 754–755. [CrossRef] [PubMed] [Google Scholar]
  • Ijiri A, Wang L, Oba T, Kawahata H, Huang Chen Y, Huang Chi Y. 2005. Paleoenvironmental changes in the northern area of the East China Sea during the past 42,000 years. Palaeogeogr Palaeoclimatol Palaeoecol 219: 239–261. [Google Scholar]
  • Inaba A, Torigoe K. 2004. Oysters in the world, part 2: systematic description of the recent oysters. Bulletin of the Nishinomiya Shell Museum 3: 1–63 (in Japanese). [Google Scholar]
  • Kao SJ, Wu CR, Hsin YC, Dai M. 2006. Effects of sea level change on the upstream Kuroshio Current through the Okinawa Trough. Geophys Res Lett 33: L16604. [Google Scholar]
  • Katoh K, Standley DM. 2013. MAFFT multiple sequence alignment software version 7: improvements in performance and usability. Mol Biol Evol 30: 772–780. [CrossRef] [PubMed] [Google Scholar]
  • Kawahata H, Ohshima H. 2004. Vegetation and environmental record in the northern East China Sea during the late Pleistocene. Global Planet Change 41: 251–273. [Google Scholar]
  • Keis M, Remm J, Ho SYW, Davison J, Tammeleht E, Tumanov IL, Saveljev AP, Männil P, Kojola I, Abramov AV, Margus T, Saarma U. 2013. Complete mitochondrial genomes and a novel spatial genetic method reveal cryptic phylogeographical structure and migration patterns among brown bears in north-western Eurasia. J Biogeogr 40: 945–920. [Google Scholar]
  • Kimura M. 1980. A simple method for estimating evolutionary rates of base substitutions through comparative studies of nucleotide sequences. J Mol Evol 16: 111–120. [CrossRef] [PubMed] [Google Scholar]
  • Kinoshita G, Suda K, Aoki D, Ohnishi N, Jogahara T, Sato JJ, Yamada F, Suzuki H. 2025. Island population dynamics since the Late Miocene: comparative phylogeography of mammalian species in three genera (Pentalagus, Diplothrix, and Tokudaia) endemic to the Central Ryukyu Islands. Mamm Study 50: 143–150. [Google Scholar]
  • Knowles LL. 2004. The burgeoning field of statistical phylogeography. J Evol Biol 17: 1–10. [Google Scholar]
  • Lam K, Morton B. 2003. Mitochondrial DNA and morphological identification of a new species of Crassostrea (Bivalvia: Ostreidae) cultured for centuries in the Pearl River Delta, Hong Kong, China. Aquaculture 228: 1–13. [Google Scholar]
  • Lam K, Morton B. 2006. Morphological and mitochondrial-DNA analysis of the Indo-West Pacific rock oysters (Ostreidae: Saccostrea species). J Molluscan Stud 72: 235–245. [Google Scholar]
  • Lambeck K, Esat TM, Potter EK. 2002. Links between climate and sea levels for the past three million years. Nature 419: 199–206. [Google Scholar]
  • Lee KE, Lee HJ, Park JH, Chang YP, Ikehara K, Itaki T, Kwon HK. 2013. Stability of the Kuroshio path with respect to glacial sea level lowering. Geophys Res Lett 40: 392–396. [Google Scholar]
  • Leigh JW, Bryant D. 2015. PopART: Full-feature software for haplotype network construction. Methods Ecol Evol 6: 1110–1116. [Google Scholar]
  • Letouzey J, Kimura M. 1985. Okinawa Trough genesis: structure and evolution of a backarc basin developed in a continent. Mar Pet Geol 2: 111–130. [Google Scholar]
  • Liu JX, Gao TX, Zhuang ZM, Jin XS, Yokogawa K, Zhang YP. 2006. Late Pleistocene divergence and subsequent population expansion of two closely related fish species, Japanese anchovy (Engraulis japonicus) and Australian anchovy (Engraulis australis). Mol Phylogenet Evol 40: 712–723. [Google Scholar]
  • Liu J, Li Q, Kong L, Yu H, Zheng X. 2011. Identifying the true oysters (Bivalvia: Ostreidae) with mitochondrial phylogeny and distance-based DNA barcoding. Mol Ecol Resour 11: 820–830. [Google Scholar]
  • Lukehurst SS, Tan SK, Tan KS, Dehm J, Hewavitharane C, Wells FE. 2026. Clarifying the identities of Saccostrea scyphophilla (Péron & Lesueur, 1807) and S. mordax (Gould, 1850) (Bivalvia: Ostreidae) from the Indo-Pacific. Zootaxa 5799: 1–31. [Google Scholar]
  • Matsuzaki KM, Itaki T, Tada R. 2019. Paleoceanographic changes in the Northern East China Sea during the last 400 kyr as inferred from radiolarian assemblages (IODP Site U1429). Prog Earth Planet Sci 6: 22. [Google Scholar]
  • McDougall C, Nenadic N, Richardson M, Healy JM. 2024. Molecular identification of intertidal rock oyster species in north-eastern Australia reveals new candidates for aquaculture. Aquaculture 587: 740838 [Google Scholar]
  • Milbury CA, Gaffney PM. 2005. Complete mitochondrial DNA sequence of the Eastern oyster Crassostrea virginica. Mar Biotechnol 7: 697–712. [Google Scholar]
  • Mu W. 2022. The complete mitochondrial genome of Saccostrea malabonensis (Ostreida: Ostreidae): characterization and phylogenetic position. Mitochondrial DNA B Resour 7: 1945–1947. [Google Scholar]
  • Nawa J. 2008. Molluscan fauna of tidal flats in the Ryukyu Islands, part 2: Okinawa, Miyako and Yaeyama Islands. Bulletin of the Nishinomiya Shell Museum 5: 1–81 (in Japanese with English abstract). [Google Scholar]
  • Nei M. Molecular Evolutionary Genetics, Columbia University Press, New York, USA, 1987. [Google Scholar]
  • Newell RIE, Fisher TR, Holyoke RR, Cornwell JC. Influence of eastern oysters on nitrogen and phosphorus regeneration in Chesapeake Bay, USA, in: R. Dame, S. Olenin (Eds.), The Comparative Roles of Suspension Feeders in Ecosystems, NATO Science Series: IV–Earth and Environmental Sciences 47, Springer, Dordrecht, Netherlands, 2005, pp. 93–120. [Google Scholar]
  • Ni G, Li Q, Kong L, Yu H. 2014. Comparative phylogeography in marginal seas of the northwestern Pacific. Mol Ecol 23: 534–548. [Google Scholar]
  • Osozawa S, Shinjo R, Armid A, Watanabe Y, Horiguchi T, Wakabayashi J. 2012. Palaeogeographic reconstruction of the 1.55 Ma synchronous isolation of the Ryukyu Islands, Japan, and Taiwan and inflow of the Kuroshio warm current. Int Geol Rev 54: 1369–1388. [Google Scholar]
  • Padial JM, De la Riva I. 2021. A paradigm shift in our view of species drives current trends in biological classification. Biol Rev 96: 731–751. [Google Scholar]
  • Pfingstl T, Wagner M, Hiruta SF, Koblmüller S, Hagino W, Shimano S. 2019. Phylogeographic patterns of intertidal arthropods (Acari, Oribatida) from southern Japanese islands reflect paleoclimatic events. Sci Rep 9: 19042. [Google Scholar]
  • Pirot J-Y, Meynell PJ, Elder D. Ecosystem Management: Lessons from Around the World. A Guide for Development and Conservation Practitioners, International Union for Conservation of Nature and Natural Resources (IUCN), Gland, Switzerland and Cambridge, UK, 2000. [Google Scholar]
  • Polson MP, Hewson W, Eernisse D, Baker PK, Zacherl DC. 2009. You say conchaphila, I say lurida: Molecular evidence for restricting the Olympia oyster (Ostrea lurida Carpenter 1864) to temperate western North America. J Shellfish Res 28: 11–21. [Google Scholar]
  • Ramos-Onsins SE, Rozas J. 2002. Statistical properties of new neutrality tests against population growth. Mol Biol Evol 19: 2092–2100. [Google Scholar]
  • Reece KS, Cordes JF, Stubbs JB, Hudson KL, Francis EA. 2008. Molecular phylogenies help resolve taxonomic confusion with Asian Crassostrea oyster species. Mar Biol 153: 709–721. [Google Scholar]
  • Reimer JD, Biondi P, Lau YW, Masucci GD, Nguyen XH, Santos MEA, Wee HB. 2019. Marine biodiversity research in the Ryukyu Islands, Japan: Current status and trends. PeerJ 7: e6532 [Google Scholar]
  • Ren J, Liu X, Guo X, Liu B. 2010. Unusual conservation of mitochondrial gene order in Crassostrea oysters: Evidence for recent speciation in Asia. BMC Evol Biol 10: 394. [Google Scholar]
  • Richardson MA, Healy J, Nenadic N, Boyle A, Connolly RM, McDougall C. 2026. Beneath the surface: the overlooked subtidalreef-building oysters of the east coast of Australia. Restor Ecol: e70382 [Google Scholar]
  • Richardson MA, Buelow C, Connolly RM, Gillies CL, Nenadic N, Porter R, Traurig M, McDougall C. 2025. Characterising tropical oyster reefs: invertebrate-environment associations and a newly documented reef building species. Mar Environ Res 208: 107136. [Google Scholar]
  • Ronquist F, Teslenko M, van der Mark P, Ayres DL, Darling A, Höhna S, Larget B, Liu L, Suchard MA, Huelsenbeck JP. 2012. MrBayes 3.2: efficient Bayesian phylogenetic inference and model choice across a large model space. Syst Biol 61: 539–542. [Google Scholar]
  • Rozas J, Ferrer-Mata A, Sánchez-DelBarrio JC, Guirao-Rico S, Librado P, Ramos-Onsins SE, Sánchez-Gracia A. 2017. DnaSP 6: DNA sequence polymorphism analysis of large datasets. Mol Biol Evol 34: 3299–3302. [Google Scholar]
  • Ruesink JL, Lenihan HS, Trimble AC, Heiman KW, Micheli F, Byers JE, Kay MC. 2005. Introduction of non-native oysters: Ecosystem effects and restoration implications. Annu Rev Ecol Evol Syst 36: 643–89. [Google Scholar]
  • Rubinoff D. 2006. Utility of mitochondrial DNA barcodes in species conservation. Conserv Biol 20: 1026–1033. [Google Scholar]
  • Saito Y. 1998. Sea levels of the last glacial in the East China Sea continental shelf. Quat Res 37: 235–242 (in Japanese with English abstract). [Google Scholar]
  • Saitou N, Nei M. 1987. The neighbor-joining method: a new method for reconstructing phylogenetic trees. Mol Biol Evol 4: 406–425. [Google Scholar]
  • Salvi D, Macali A, Mariottini P. 2014. Molecular phylogenetics and systematics of the bivalve family Ostreidae based on rRNA sequence-structure models and multilocus species tree. PLoS One 9: e108696. [Google Scholar]
  • Salvi D, Mariottini P. 2017. Molecular taxonomy in 2D: a novel ITS2 rRNA sequence-structure approach guides the description of the oysters' subfamily Saccostreinae and the genus Magallana (Bivalvia: Ostreidae). Zool J Linn Soc 179:263–276. [Google Scholar]
  • Sekino M, Sato S, Hong JS, Li Q. 2012. Contrasting pattern of mitochondrial population diversity between an estuarine bivalve, the Kumamoto oyster Crassostrea sikamea, and the closely related Pacific oyster C. gigas. Mar Biol 159: 2757–2776. [Google Scholar]
  • Sekino M, Yamashita H. 2013. Mitochondrial DNA barcoding for Okinawan oysters: a cryptic population of the Portuguese oyster Crassostrea angulata in Japanese waters. Fish Sci 79: 61–76. [Google Scholar]
  • Sekino M, Yamashita H. 2016. Mitochondrial and nuclear DNA analyses of Saccostrea oysters in Japan highlight the confused taxonomy of the genus. J Molluscan Stud 82: 492–506. [Google Scholar]
  • Smith RS, Pruett JL. 2025. Oyster restoration to recover ecosystem services. Annu Rev Mar Sci 17: 83–113. [Google Scholar]
  • Snow M, Fotedar S, Wilson NG, Kirkendale LA. 2023. Clarifying the natural distribution of Saccostrea Dollfus and Dautzenberg, 1920 (edible rock oyster) species in Western Australia to guide development of a fledgling aquaculture industry. Aquaculture 566: 739202. [Google Scholar]
  • Song H, Buhay JE, Whiting MF, Crandall KA. 2008. Many species in one: DNA barcoding overestimates the number of species when nuclear mitochondrial pseudogenes are coamplified. Proc Natl Acad Sci USA 105: 13486–13491. [Google Scholar]
  • Song N, Zhang XM, Sun XF, Yanagimoto T, Gao TX. 2010. Population genetic structure and larval dispersal potential of spotted tail goby Synechogobius ommaturus in the north-west Pacific. J Fish Biol 77: 388–402. [Google Scholar]
  • Struck TH, Feder JL, Bendiksby M, Birkeland S, Cerca J, Gusarov VI, Kistenich S, Larsson KH, Liow LH, Nowak MD, Stedje B, Bachmann L, Dimitrov D. 2018. Finding evolutionary processes hidden in cryptic species. Trends Ecol Evol 33: 153‒163. [Google Scholar]
  • Tajima F. 1989. Statistical method for testing the neutral mutation hypothesis by DNA polymorphism. Genetics 123: 585‒595. [CrossRef] [PubMed] [Google Scholar]
  • Tamura K, Stecher G, Peterson D, Filipski A, Kumar S. 2013. MEGA6: Molecular Evolutionary Genetics Analysis v6.0. Mol Biol Evol 30: 2725–2729. [Google Scholar]
  • Tan SK, Wells FE, Tan KS, Lukehurst SS, Morgan M, Fotedar S. 2025. Identity of the enigmatic oyster Saccostrea cuccullata (Bivalvia: Ostreidae). J Molluscan Stud 91: eyaf007. [Google Scholar]
  • Templeton AR, Crandall KA, Sing CF. 1992. A cladistic analysis of phenotypic associations with haplotypes inferred from restriction endonuclease mapping and DNA sequence data. III. Cladogram Estimation. Genetics 132: 619–633. [Google Scholar]
  • Vincent-Vogt NS, Mitarai S. 2020. A persistent Kuroshio in the glacial East China Sea and implications for coral paleobiogeography. Paleoceanogr Paleoclimatol 35: e2020PA003902. [Google Scholar]
  • Wang J, Wang P. 1980. Relationship between sea-level changes and climatic fluctuations in East China Sea since late Pleistocene. Acta Geogr Sin 35: 299–312 (in Chinese with English abstract). [Google Scholar]
  • Watanabe N, Arai K, Otsubo M, Toda M, Tominaga A, Chiyonobu S, Sato T, Ikeda T, Takahashi A, Ota H, Iryu Y. 2023. Geological history of the land area between Okinawa Jima and Miyako Jima of the Ryukyu Islands, Japan, and its phylogeographical significance for the terrestrial organisms of these and adjacent islands. Prog Earth Planet Sci 10: 40. [Google Scholar]
  • Wu JH, Ju YM, Hsiao ST. 2019. The complete mitochondrial genome of Saccostrea kegaki (Pterioida, Ostreidae). Mitochondrial DNA B Resour 4: 642–643. [Google Scholar]
  • Wu X, Xiao S, Yu Z. 2013. Mitochondrial DNA and morphological identification of Crassostrea zhanjiangensis sp. nov. (Bivalvia: Ostreidae): a new species in Zhanjiang, China. Aquat Living Resour 26: 273–280. [Google Scholar]
  • Xia J, Wu X, Xiao S, Yu Z. 2014. Mitochondrial DNA and morphological identification of a new cupped oyster species Crassostrea dianbaiensis (Bivalvia: Ostreidae) in the South China Sea. Aquat Living Resour 27: 41–48. [Google Scholar]
  • Xiao S, Li A, Jiang F, Li T, Wan S, Huang P. 2004. The history of the Yangtze River entering sea since the Last Glacial Maximum: A review and look forward. J Coast Res 20: 599–604. [Google Scholar]
  • Xu J, Chan TY, Tsang LM, Chu KH. 2009. Phylogeography of the mitten crab Eriocheir sensu stricto in East Asia: Pleistocene isolation, population expansion and secondary contact. Mol Phylogenet Evol 52: 45–56. [Google Scholar]
  • Xu X, Oda M. 1999. Surface-water evolution of the eastern East China Sea during the last 36,000 years. Mar Geol 156: 285–304. [Google Scholar]
  • Yang Z. 1994. Maximum likelihood phylogenetic estimation from DNA sequences with variable rates over sites: Approximate methods. J Mol Evol 39: 306–314. [Google Scholar]
  • Zhang Y, Ma H, Han C, Zhao Z, Wei Z, Duan X, Li J, Qin Y, Li J, Yu Z, Pan Y, Zhang Y. 2025. Comparison of genetic diversity and population structure of Saccostrea mordax lineages A, B, and C across South China Sea and Celebes Sea. Front Mar Sci 12: 1593382. [Google Scholar]

Cite this article as: Sekino M, Oguchi M, Sanda T, Iwasaki T. 2026. Phylogenetic diversity and comparative historical demography within a species complex of Indo-West Pacific Saccostrea oysters. Aquat. Living Resour. 39: 21. https://doi.org/10.1051/alr/2026013

All Tables

Table 1

Sampling sites of specimens analyzed in the present study.

In the text
Table 2

Frequency of each lineage at each sampling site.

In the text
Table 3

Mitochondrial diversity within lineages and sublineages based on concatenated sequences of the four mitochondrial genes (2,146 bases).

In the text
Table 4

FST (ΦST) estimated between population samples within each lineage and sublineage.

In the text

All Figures

Thumbnail: Fig. 1 Refer to the following caption and surrounding text. Fig. 1

Geographic map showing sampling locations. Sampling sites examined in the present study are indicated by italicized names in red. Other sites were surveyed by Sekino and Yamashita (2013, 2016). The warm, northward-flowing Kuroshio Current passes between Taiwan and Ishigaki/Iriomote and flows along the Pacific coast of Kyushu, Shikoku, and Honshu Islands, Japan. OT and YT denote the Okinawa Trough and the Yonaguni Depression, respectively. The map was produced using base maps processed with Ocean Data View version 5.6.5 (developed by Reiner Schlitzer; https://odv.awi.de/; last accessed on 6 August, 2025).
In the text
Thumbnail: Fig. 2 Refer to the following caption and surrounding text. Fig. 2

Bayesian phylogenetic tree of Saccostrea specimens based on concatenated sequences of four mitochondrial genes. The mordax group (mordax A, B, and C) was used as the outgroup. Lineages detected in the present study are indicated by bold, enlarged, and italicized names. Posterior probability values (≥ 0.80) are shown for nodes. Lineages were identified by cross-referencing single-gene trees (Supplementary Materials 6 and 7).
In the text
Thumbnail: Fig. 3 Refer to the following caption and surrounding text. Fig. 3

Shells of observed lineages. Each specimen is shown in up to four views and is identified by a unique label designated as follows: lineage or sublineage name, followed by the specimen ID and shell height (SH; mm, in parentheses). Abbreviations used in the figure are as follows: BV, both valves; OLV, outer left valve; ILV, inner left valve; ORV, outer right valve; IRV, inner right valve; GL, growth lamellae; HS, hyote spine.
In the text
Thumbnail: Fig. 4 Refer to the following caption and surrounding text. Fig. 4

Statistical parsimony (TCS) networks of haplotypes in each Saccostrea lineage based on concatenated sequences of four mitochondrial genes. The number of tick marks on each edge represents the number of mutational steps between nodes (no tick mark is placed for single step mutations). Node size is proportional to haplotype frequency. Lineages B and L with small sample sizes (N = 2 and 1, respectively) are omitted.
In the text
Thumbnail: Fig. 5 Refer to the following caption and surrounding text. Fig. 5

Bayesian skyline plots (BSPs) estimated for Saccostrea lineages and sublineages. Shaded areas represent the 95% highest posterior density interval. BSP analyses were not performed for the following lineages and sublineages owing to small sample sizes: Lineages B (N = 2) and L (N = 1), and two sublineages within Lineage D (D2, N = 6; and D3, N = 3).
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.