Phylogeny and Diversification Patterns among Vesicomyid Bivalves

Vesicomyid bivalves are among the most abundant and diverse symbiotic taxa in chemosynthetic-based ecosystems: more than 100 different vesicomyid species have been described so far. In the present study, we investigated the phylogenetic positioning of recently described vesicomyid species from the Gulf of Guinea and their western Atlantic and Pacific counterparts using mitochondrial DNA sequence data. The maximum-likelihood (ML) tree provided limited support for the recent taxonomic revision of vesicomyids based on morphological criteria; nevertheless, most of the newly sequenced specimens did not cluster with their morphological conspecifics. Moreover, the observed lack of geographic clustering suggests the occurrence of independent radiations followed by worldwide dispersal. Ancestral character state reconstruction showed a significant correlation between the characters “depth” and “habitat” and the reconstructed ML phylogeny suggesting possible recurrent events of ‘stepwise speciation’ from shallow to deep waters in different ocean basins. This is consistent with genus or species bathymetric segregation observed from recent taxonomic studies. Altogether, our results highlight the need for ongoing re-evaluation of the morphological characters used to identify vesicomyid bivalves.


Introduction
Chemosynthetic ecosystems are found worldwide in the deep ocean and harbour specific communities that have a high level of endemism [1]. Some large families or genera, among which vesicomyid bivalves are one of the most remarkable examples, are distributed across a wide diversity of habitats including hydrothermal vents, cold seeps and whale carcasses. These inhabitants of the deep sea have developed intrinsic features that enable them to live in these extreme environments [2,3,4]. Similarly to many other taxa associated with chemosynthetic ecosystems, vesicomyids have chemoautotrophic bacteria in their gills that use sulphide components and provide their hosts with nutriments, thereby creating symbiotic relationships [5,6,7,8]. Since the late 1970s and the discovery of chemosynthetic ecosystems, vesicomyid species have been described on the basis of morphological characters, and a taxonomical framework has been established. Within the currently described 112 extant and 4 fossil species, 47 have been described from cold seeps, 67 from putative seep environments and 14 from hydrothermal vents [9,10], with some of these species sharing several environment types.
Recent revisions of the taxonomy of Vesicomyidae based on soft tissues and shell morphology include the description of several new species around the world [9,10,20,21,22,23]. These studies advocate the separation of the family into two subfamilies (Vesicomyinae and Pliocardiinae), re-description of the Calyptogena genus and the creation of new genera separating species into distinct genera. The subfamily Vesicomyinae comprises only the genus Vesicomya (small-sized bivalves without reduced gut and subfilamental tissue on the gills) and the subfamily Pliocardiinae (medium and large body size with subfilamental tissue on the gills), which includes 15 genera. Consequently, most vesicomyids from hydrothermal vent and cold seep ecosystems initially attributed to Calyptogena and Vesicomya genera have been assigned to new erected genera (e.g. Abyssogena) or old ones (e.g. Pliocardia) [23,24]. Nevertheless, the revision is under progress and some of these species are still assigned to Calyptogena. The Pliocardiinae are generally of large size, harbour sulphide-oxidizing bacteria and are highly specialized to sulphide-rich environments while symbiosis is not proved in small size Vesicomyinae not restricted to these environments [23,24]. Both Vesicomyinae and Pliocardinae are distributed worldwide with a broad vertical range [24]. However, the geographic distribution and depth range vary among genera. Most genera (11) have colonised at least two oceans and only five have a strictly regional distribution. However, the vertical distribution of each genus is narrow except for two genera (Vesicomya and Isorropodon) that have colonised zones ranging from the sublittoral to abyssal. Recent revisions of vesicomyid genera demonstrate that genera generally occur in only one bathymetric zone [23,24,25].
In the present study, we investigated the phylogenetic relationships among recently described vesicomyid species in the Gulf of Guinea and their western Atlantic and Pacific counterparts using mitochondrial DNA COI sequence data. The reconstructed phylogeny was used to: (1) compare the taxonomic status of recently described species based on morphological characters with molecular based taxonomy and discuss possible discrepancies between both approaches, and (2) place these species in their phylogenetic context and perform an ancestral character state reconstruction to test previous hypotheses suggesting a relationship between vesicomyid diversification and the habitat/depth that they occur.

Specimens
Several specimens of nine putative vesicomyid species were collected worldwide in cold seeps, including five in the Gulf of Guinea, and identified morphologically to the species or genus level and recently assigned to existing or new genus [9,10,25] Figure 1). Sample collection was performed using deep submergence vehicles (Nautile during the DIAPISUB and NAUTI-MATE cruises, the ROV Victor 6000 during the BIOZAIRE 2 cruise, and Quest 4000 (MARUM) during the RV Meteor expedition M70-3 [28] and during the expedition M76-3b GUINECO (MARUM) as well as by surface dredges (BIOZAIRE 3 cruise) ( Table 1). In most cases, specimens were dissected after retrieval or frozen whole at 280uC. Species identification was performed by a taxonomist (R. Von Cosel) [10] or by the authors for new specimens of already described species from the Gulf of Guinea (KO, CD). Specimens unidentified to the species level and collected in the Middle America Trench off Mexico in the Pacific Ocean were named C. n sp (CI) and C. n sp (CII) (Tables 1 and S1).

DNA extraction, amplification and sequencing
Total genomic DNA from the adductor muscles of 29 Ch. regab, three L. chuni, two I. perplexum, and one individual for each of the six remaining species, was isolated using the CTAB extraction method (Doyle and Doyle 1987). Purified nucleic acids were resuspended in DNase-free water and stored at 220uC. Universal (LCOI1490+HCOI2198 [29]) and vesicomyid-specific (VesL-CO+VesHCO [12]) primers were used to amplify a partial fragment of the mitochondrial COI gene. PCR amplifications were carried out in 50 ml reactions with the following final concentrations: 16 Taq DNA polymerase buffer (Promega), 2.5 mM MgCl2, 0.20 mM dNTPs, 0.6 mM each primer, 0.5 U Taq DNA polymerase (Promega) and 5 ml of purified DNA. The following profile was used: incubation at 94uC for 90 s, 35 cycles of incubation at 94uC for 30 s, 45uC for 60 s, 72uC for 75 s and a final extension at 72uC for 10 min. PCR products (40) were then sent to GATC Biotech for sequencing in both directions. Accession numbers are given in Table S1.

Molecular taxonomy and phylogenetic analyses
Mitochondrial COI gene sequences from all vesicomyid bivalves available in the GenBank database were included in all analyses. As Genbank contains data obtained prior to the new taxonomic revision, or even data from specimens either unnamed or later re-identified with a different species name, we renamed unidentified or misidentified specimens of sequences retrieved from GenBank according to their identification in the most recent taxonomic revisions by expert taxonomists [20,23,24]. For example Calyptogena sp. ''Clam #2 [12] became Abyssogena southwardae following [23]. These new names given in table S1 are those used in the phylogenetic tree.
To assess the phylogenetic relationships of vesicomyids, we used 87 COI nucleotide sequences (including 12 sequences obtained in this study). Given their close relationship with vesicomyids, species belonging to families Arcticidae (Arctica islandica) Corbiculidae (Corbicula fluminaris and Corbicula japonica), Glauconomidae (Glauconome chinensis), Veneridae (Venus antique, Anomalocardia brasiliana and Mercenaria mercenaria) were chosen as outgroups [30]. Alignment of nucleotide sequences was conducted using Geneious Pro version 4.6 [31], and verified by eye to maximise positional homology. The Akaike's information criterion (AIC) [32] implemented in Modeltest v. 3.7 [33] was used to determine the evolutionary model that best fit the data set. PhyML v. 2.4.4 [34] was used to estimate the maximum likelihood (ML) tree using the model supported by Modeltest (GTR+I+G) model [35], and to test the robustness of the inferred trees using non-parametric bootstrap proportions (BPs) with 1000 pseudoreplicates. ML analyses were carried out through the freely available web-based portal Bioportal (http://www.bioportal.uio.no). Specimens provisionally assigned to existing genera or whose description is in progress following [24], are underlined by different colours.
Ancestral character state reconstruction for the evolution of depth limits and habitat during vesicomyid evolution (53 species) was performed using MacClade v. 4.08 [36] to test whether the evolutionary pattern of habitat and depth was significantly correlated with the molecular phylogeny. The analyses based on the maximum depth limits of each taxon were performed using 13 depth intervals (see Figure 2 for further details) whereas ancestral reconstruction of the character ''habitat'' included four different categories (Seeps, Vents, Seeps and Vents, and Coastal). Each depth limit and different habitat was mapped onto the reconstructed ML topology. To obtain the null distribution of the number of transitions of the character ''depth'' and ''habitat'' 1000 random trees were generated. Tree lengths were compared with the mean number of steps of the null distribution of each reconstruction to test whether the observed number of transitions was significantly more correlated with the phylogeny than expected at random.

Molecular versus morphological taxonomy
The ML analysis based on a 630 bp fragment of the mitochondrial COI gene yielded the tree (2lnL = 5653.83) shown in Figure 2. Despite the rather low support of the deepest nodes, more recently diversified clusters of genera and species showed reasonable bootstrap support (.75%) on the ML tree.
All the species newly analysed for this study (Ch. regab, L. chuni, C. valdiviae, E. guiness, W. guineensis, I. perplexum, L. myriamae, C. n sp.I and C.n sp. II), identified by expert taxonomists (R. Von Cosel and K. Olu), did actually group with their conspecifics from GenBank. C. valdiviae grouped with the species complex C. pacifica/V. lepta. The species identified as L. chuni was included within a clade of several species from the Japan Trench (BP = 92%). A relationship between an unidentified species from the Middle America Trench (C. n. sp I) and E. extenta (AF114388) is depicted in the ML tree. Two other unidentified specimens, one from the same area (C. n sp CII) and one from the Costa Rica accretionary wedge (AF114390) grouped together ( Figure 2) with high bootstrap support (100%).
The molecular taxonomical reconstruction reported here provided limited support even for the recent taxonomic revision of vesicomyids based on morphological criteria. Upon 14 morphologically determined genera present in the data set only three were recovered as monophyletic: Christineconcha, Abyssogena and Pliocardia. The four haplotypes belonging to the genus Christineconcha clustered together with a maximum divergence of 0.6%, supporting the occurrence of a single species. Species belonging to the newly described genus Abyssogena [23] were monophyletic as well as species belonging to the genus Pliocardia.   The monophyly of the genus Calyptogena defined according to morphological criteria is broken due to the grouping of C. fausta with V. lepta. This clustering is quite unexpected but was also identified by Kojima et al. [16] as one of the pairs of trans-pacific species supported by bootstrap values greater than 90% (here 95%, Figure 2). This species may nevertheless deserve to be reexamined together with V. lepta specimens during the on-going genera revision. Species provisionally assigned to other genera by Krylova and Sahling [24] are not included in this cluster (e.g. C. magnifica). Similarly, one of the species provisionally assigned to the genus Phreagena (C. nankaiensis) is not included in the cluster gathering the other three species of this genus, C. soyoae, C. kilmeri and C. okutani [16].
Despite recent revisions [20,21,24], the reconstructed phylogeny presented here did not correspond to the current morphologybased taxonomy for the remaining genera, including Wareniconcha, Laubiericoncha and Isorropodon. Indeed, the genus Wareniconcha whose type species is W. guineensis from the Gulf of Guinea did not form a monophyletic cluster with the other provisionally assigned species V. lepta. Similarly, species of the genus Phreagena did not group together. Moreover, both species L. myriamae and L. chuni of the recently described genus Laubiericoncha do not form a monophyletic group, and I. perplexum and C. fossajaponica belonged to Isorropodon. Nevertheless these tree genera are represented by only two species in this study and additional taxa would be needed to confirm this observation.
In contrast, this phylogenetic reconstruction revealed wellsupported groups that have yet to be defined, but confirmed some suggestions from previous studies. Indeed, Kojima et al. [16] identified 6 pairs of closely related (90 to 100% of similarity) transpacific species that were also observed here (e.g. C. soyoae and C. kilmeri, V. magnocultellus and V. gigas). The newly sequenced L. chuni from West Africa is closely related (divergence ,2%) to an unidentified species sampled off Sanriku, Japan (BP = 100%) that suggest also Atlantic-Pacific similarities. Morphological identification of this unidentified specimen and phylogeographic analysis with more samples from both geographic areas would help determining whether these specimens belong to a single species exhibiting recent or contemporary exchanges, or sub-species with on-going divergence and lineage sorting. Moreover, both specimens were also close to C. laubieri and V. gigas from western and eastern Pacific Ocean, although they do not belong to the morphologically defined genera thus-far identified in this part of the world. Nevertheless, L. chuni and C. laubieri share many characters according to current study by Krylova (pers.com.) Further study on the systematic status of these specimens is needed. Other studied species were not related with wellsupported nodes to previously known specimens. Barcoding and phylogenetic positioning of other species, particularly those from the Atlantic Ocean, may be necessary for a better understanding of the relationships among species.
As pointed out by Kojima et al. (2006) [19], C. solidissima and C. kawamurai appear as closely related and the species identified as C. soldissima is probably a junior synonym of C. kawamurai. Moreover, C. fossajaponica and both unidentified species from the Kurile Trench reported by Kojima et al. 2004 [16] and Okutani et al. 2009 [37] were very closely related. These three specimens may be conspecific as previously suggested by Okutani et al. 2009 [37].
Two cryptic species complexes have been highlighted by previous phylogenetic studies [12,38]. The first complex corresponds to the species V. lepta and C. pacifica phylogenetically analysed by Peek et al. [15] who attributed the divergence to their respective distribution in vents versus seeps. In a later analysis including other lineages, Goffredi et al. [13] identified three more lineages within the same complex and suggested their bathymetric segregation. The divergence of the five lineages they identified in total is confirmed in the present study, with an obvious cluster gathering C. pacifica specimens (together with C. storobogatovi from hydrothermal vents) and another cluster of 4 lineages (V. lepta, mtII, mtIII, mtV) supported by high bootstrap values ( Figure 2). Additionally, the species C. fausta also clustered in this group with the mtII lineage as its closest relative. C. kilmeri, initially included in the C. kilmeri/V. gigas complex of species, is clearly separated from V. gigas but undistinguishable from its Western Pacific counterpart C. soyoae, and closely related to another Japanese margin species (C. okutanii). In contrast, V. gigas are spread throughout three paraphyletic lineages supported by high bootstrap values, some also gathering Western Pacific sister species (C. laubieri, V. magnocultellus) and the new Western Africa Laubiericoncha chuni. (Figure 2).
To clarify the taxonomy of the two ''new'' species cryptic complexes ''V'' gigas and ''V.'' lepta, further morphologic and genetic studies are needed.
Both genera Archivesica and Laubiericoncha have gills with two pairs of demibranchs, a character also shared with ''C.''kilmeri and E. extenta (in fact the whole cluster was supported by a 99% bootstrap value, except C. okutani, C. kawamurai and C. soldissima). On the contrary, the genus Calyptogena, the ''V''. lepta complex, and the genera Abyssogena, Christineconcha and Isorropodon have a single pair of demibranchs. The presence of two demibranchs is considered as an ancestral character [20]. Besides this group of species, the species''C.''magnifica (new genus [24]) and ''C'' packardana branching in distinct clusters also shared this putative ancestral character. Ambiguous generic diagnoses may result from morphological convergences. The common occurrence of reduction in the hinge structures that would have evolved independently several times in the family, have been suggested to explain the lack of monophyly of generic names [24]. Soft part characters, such as specialised gill related to symbiosis (reduced outer gill demibranch, presence of wide interlamella septa), absence of pallial sinus or number of pallial apertures, may deserve careful consideration when describing new genera/species and it would also likely be interesting to compare their distribution with molecular phylogeny.
The discrepancies observed here between the most recent morphological taxonomy and our reconstructed molecular phylogeny highlights the necessity to carry on a re-evaluation of the morphological characters used in the identification of vesicomyid bivalves by performing as much as possible a systematic parallel morphological and molecular characterization. Furthermore, the molecular characterization of species assigned to some genera such as Vesicomya, Isorropodon and Callogonia, although not yet available, would increase the resolution of the phylogenetic tree.

Phylogenetic perspective
The family Vesicomyidae constitutes a monophyletic group of relatively recent origin (50-80 million years) [12,39]. However, results obtained here support a scenario of multiple migrations between eastern and western Pacific Ocean, as suggested by Kojima et al. (2004) [16]. Accordingly, several species, such as V. sp mtII and C. fausta and C. kilmeri and C. soyoae, cluster together despite their distinct distributional ranges across the Pacific Ocean.
The currently known biogeographic provinces and pathways, such as the Atlantic Equatorial Belt, that have been identified and are shared among various chemosynthetic taxa [40,41,42,43], do not seem to apply to the set of vesicomyid species analysed in this study. Several species or species complexes of Bathymodiolus bivalves, gastropods and crustaceans exhibit an amphi-Atlantic distribution [40], whereas the vesicomyids collected from the eastern to western Atlantic (Gulf of Mexico, Barbados, Mid-Atlantic Ridge, Gulf of Guinea) are not closely related. Instead, they are scattered across the phylogenetic tree, except one specimen from the West Florida Escarpment (AF008280) that resolved as the sister lineage of a clade including A. southwardae from the Barbados Accretionary Prism and V. sp from the Mid-Atlantic Ridge. Further, none of the five vesicomyids sampled in the Gulf of Guinea across a large depth range (700-4000 m) seem to be related to any Atlantic seep or vent counterparts. This may be the result of either incomplete sampling or limited exchanges between the western and eastern parts of the Atlantic. Moreover, other species collected from the Gulf of Guinea did not form a monophyletic group either, although some of them occur in sympatry. For instance, Ch. regab and L. chuni were sampled on the Regab pockmark in the same species aggregate but are quite distant in the phylogenetic tree. Similarly, C. valdiviae and E. guiness occur on the same pockmark in the Guiness area, but were also distantly related in both phylogenetic and morphological [10] analyses. W. guineensis did not cluster with other species from the Gulf of Guinea. The lack of clustering of these sympatric species is similar to patterns observed for species that inhabit distinct geographic regions and oceans. This result suggests the occurrence of multiple series of speciation events followed by large-scale dispersal.
The phylogenetic reconstruction performed in the present study supports the findings of previous studies, suggesting that the speciation of the C. pacifica/V. lepta complex seems to be both related to the habitat (seep/vent) [15] and depth [13]. So far, C. pacifica seems restricted to seeps whereas V. lepta (mtIV) restricted to vents. The recent sampling of this complex at both seeps and vents in the basin would complete this analysis (Arnaud-Haond unpubl.). On the contrary, the lineages mtIII and mtV were only sampled at seeps while mtII is distributed in both ecosystems from Monterey Bay to Juan de Fuca along a relatively large depth range [13].
Interestingly both complexes of apparently cryptic species V.lepta and V. gigas have a wide bathymetric range and likely diverged in several lineages submitted to distinct environmental conditions. Both species complexes occur both at seeps and vents. In contrast C. pacifica and C. kilmeri have narrow bathymetric ranges and have been observed only in seeps thus far [20]).

Preliminary attempt of character tracing
The influence of environmental parameters on the diversification and speciation in Vasicomyids was supported by the ancestral character reconstruction using the habitat and depth ranges of vesicomyid species (Figure 3). Both analyses showed a significant correlation between patterns of diversification and habitat/ maximum depth where the species occur regardless the lack of resolution of some of the basal branches.
Vesicomyid bivalves exhibit a worldwide distribution and have colonized different chemosynthetic environments: hydrothermal vents, cold seeps, and whale falls [2,3,4]. The results reported here suggest that other factors rather than geographic boundaries have contributed to shape speciation patterns in vesicomyids. Recurrent events of 'stepwise speciation' seem to have occurred from shallow to deep waters in different ocean basins [13], [44]. Ancestral character state reconstructions performed here support this general trend (Fig. 2). In a recent review, Krylova and Sahling [24], reported the geographic and vertical distributions of vesicomyid genera and the species of the Gulf of Guinea seem to be depth-segregated [10]. Most genera occur in a single bathymetric zone; some have a bathyal (e.g. Callogonia, Ectenagena, Elenaconcha, etc.) distribution, whereas others are restricted to the abyssal zone (Abyssogena and Christineconcha) [24]. The deep Abyssogena genus clusters in the ML phylogenetic tree with V. kuroshimana and V. crenulomarginata (Figure 3), which are shallower species. Other species such as E. guiness, C. elongata and C. kawamurai show a more basal position in the phylogenetic tree although assigned to genera with upper bathyal distributions. These results would be in agreement with a putative pattern of diversification from the upper bathyal to abyssal depths. In other chemiosynthetic bivalve families such as Bathymodiolinae, phylogenetic analyses indicate that a seep species ancestor ''gave rise to deep-sea hydrothermal vent mytilids'' in an evolutionary progression from shallow-to deep-water habitats [45]. More recently, the ''wooden steps to deep-sea vent'' hypothesis has been supported by the basal position of sunken-wood species in the phylogeny [46,47]. None of the whale-fall specimens (e.g., V. gigas) studied here and also inhabiting seeps and vents, showed a basal position in the phylogenetic tree but instead clustered with other conspecific specimens (Table S1). W. guineensis sampled in the vicinity of the Congo River channel could be associated to terrestrial organic matter accumulations. However, only vesicomyids larvae have been collected to date on sunken wood, during an in situ colonisation experiment with wood cubes deployed in the Nile deep-sea fan [48]. The absence of vesicomyid adults on sunken wood may either reflect their real status, incomplete sampling or extinction of the shallower species. The ''wooden step'' hypothesis in vesicomyids remains therefore an open question. The low resolution of the deepest nodes would render any conclusive interpretation of this character analysis highly speculative. The trend observed here on more recent and better-supported nodes opens however interesting perspectives. Also, adding further markers to the analysis could improve the resolution of the phylogeny and our understanding about the origin and evolutionary history of those widely distributed chemosynthetic taxa.