Skip to main content

A hybrid zone between Bathymodiolusmussel lineages from eastern Pacific hydrothermal vents



The inhabitants of deep-sea hydrothermal vents occupy ephemeral island-like habitats distributed sporadically along tectonic spreading-centers, back-arc basins, and volcanically active seamounts. The majority of vent taxa undergo a pelagic larval phase, and thus varying degrees of geographical subdivision, ranging from no impedance of dispersal to complete isolation, often exist among taxa that span common geomorphological boundaries. Two lineages of Bathymodiolus mussels segregate on either side of the Easter Microplate, a boundary that separates the East Pacific Rise from spreading centers connected to the Pacific-Antarctic Ridge.


A recent sample from the northwest flank of the Easter Microplate contained an admixture of northern and southern mitochondrial haplotypes and corresponding alleles at five nuclear gene loci. Genotypic frequencies in this sample did not fit random mating expectation. Significant heterozygote deficiencies at nuclear loci and gametic disequilibria between loci suggested that this transitional region might be a ‘Tension Zone’ maintained by immigration of parental types and possibly hybrid unfitness. An analysis of recombination history in the nuclear genes suggests a prolonged history of parapatric contact between the two mussel lineages. We hereby elevate the southern lineage to species status as Bathymodiolus antarcticus n. sp. and restrict the use of Bathymodiolus thermophilus to the northern lineage.


Because B. thermophilus s.s. exhibits no evidence for subdivision or isolation-by-distance across its 4000 km range along the EPR axis and Galápagos Rift, partial isolation of B. antarcticus n. sp. requires explanation. The time needed to produce the observed degree of mitochondrial differentiation is consistent with the age of the Easter Microplate (2.5 to 5.3 million years). The complex geomorphology of the Easter Microplate region forces strong cross-axis currents that might disrupt self-recruitment of mussels by removing planktotrophic larvae from the ridge axis. Furthermore, frequent local extinction events in this tectonically dynamic region might produce a demographic sink rather than a source for dispersing mussel larvae. Historical changes in tectonic rates and current patterns appear to permit intermittent contact and introgression between the two species.


Hybrid zones, regions of intergradation between genetically distinct populations or species, have been documented for a number of terrestrial and aquatic organisms. The locations, shapes and ages of these zones, and the unique dispersal, behavioral and fitness characteristics of the participating organisms have motivated a number of models that invoke endogenous (i.e. genetic) or exogenous (i.e. environmental) forces reviewed by [15]. Intergradation can result from primary or secondary processes that are difficult to distinguish [6]. Primary intergradation, which occurs at the contact boundaries between distinct habitats or along steep environmental gradients, is associated with parapatric speciation. Allelic and phenotypic frequencies often exhibit step-clines shaped by disruptive selection across the contact zone. Secondary intergradation, resulting from contact between previously allopatric populations, can produce similar step-clines. Different histories of dispersal and selection can produce similar patterns of concordance in allelic frequency clines at multiple gene loci, though subtle and discriminating footprints may remain in the genetic architecture of intergrading lineages [1].

Most examples of hybrid zones come from terrestrial and freshwater environments that have experienced past climatic changes — for example, at suture zones where formerly isolated biotic assemblages attain secondary contact [7, 8]. Hybrid zones are uncommon in the open ocean, because marine environments tend to be more homogeneous temporally and spatially, and many marine species have high dispersal capabilities [9]. Nonetheless, marine hybrid zones are found in coastal environments that experience steep environmental gradients or were impacted by past climatic changes [10]. Hybridization is prevalent among subtidal mussels, Mytilus edulis, M. galloprovincialis, and M. trossulus, members of the “blue mussel complex” [11, 12]. These closely related species are known to readily hybridize wherever two blue mussels are sympatric [1315]. A mid-ocean hybrid zone was first discovered in deep-sea hydrothermal vent mussels (Mollusca: Bivalvia: Mytilidae) living along the Mid-Atlantic Ridge [16]. Hydrothermal vents differ from other open oceanic habitats because they occur as discrete island-like habitats distributed along the global mid-ocean ridge system (Figure 1), in disjunct back-arc spreading centers, and on active seamounts [17]. Because vent communities rely almost entirely on reduced volcanic gases (primarily H2S and CH4,) and chemosynthetic bacteria for primary productivity, habitat quality can vary greatly in time and space [18]. The larvae of most vent animals are dispersed by hydrographic circulation that can be disrupted by geomorphological structures, such as ridge offsets (e.g., transform faults), intervening seamounts, and other bathymetric features [19, 20]. Consequently, genetic studies have often revealed geographic subdivision and species boundaries associated with such features [21, 22].

Figure 1

Sample localities for eastern Pacific Bathymodiolus are marked by a color filled circle in accordance with the following regions: NEPR denotes the northern East Pacific Rise; SEPR = southern East Pacific Rise; GAR = Galápagos Rift; and PAR = Pacific-Antarctic Ridge. Numbers indicate spreading rates of ridge segments in mm/year.

The Mid-Atlantic Ridge (MAR) mussels, Bathymodiolus azoricus and B. puteoserpentis, mussels segregate latitudinally and hybridize at an intermediate location that is relatively inhospitable for mussels [16, 23]. The hybrid population exhibits an excess of the parental cytonuclear genotypes, leading Won et al.[24] to hypothesize that the contact region is Tension Zone [sensu 1] maintained primarily by immigration from the parental zones and possibly hybrid unfitness. A coalescence analysis of the nuclear GF1B locus identified recombination events between branch tips belonging to the parental lineages [25]. The degree of mitochondrial divergence and the apparent absence of earlier inter-lineage recombination events at GF1B were interpreted as evidence for a period of allopatric differentiation beginning about 0.760 million years (Myr) ago, followed by prolonged secondary contact and introgression.

Discovery of the MAR hybrid zone stimulated our interests in other potential contact zones along the global mid-ocean ridge system. Our target was the Easter Microplate boundary (EMB, Figure 1), where vent fauna associated with the 7000-km long East Pacific Rise (EPR) contact a related fauna associated with a NE segment of the Pacific-Antarctic Ridge (PAR) [22]. The EMB separates geminate species of bythograeid crabs [26, 27], lepetodrilid limpets [28], and Bathymodiolus mussels [29]. The EMB does not contribute to significant geographical partitioning of all lepetodrilid species [28] and annelid species [30, 31], however. Instead, this region acts as a “variable dispersal filter” with affects that depend on the life history characteristics and dispersal modes of individual taxa [22]. To better define this contact zone, we explored and sampled vents that were previously reported to exist in the region between 23°S on the NW flank of the Easter Microplate and 38°S latitude along the PAR [32, 33].

Herein, we report the discovery of a zone of intergradation between two mussel lineages that contact one another at the EMB. Bathymodiolus thermophilus Kenk & Wilson [34] is distributed on the Galápagos Rift (GAR) and East Pacific Rise (EPR) between 13°N and 21°S latitude. A related lineage, Bathymodiolus aff. thermophilus[29], was found at 31–32°S latitude on segments of the EPR axis connecting the Easter and the Juan de Fernandez microplates. The mussel lineages are reciprocally monophyletic for distinct mitochondrial COI haplotypes that are 4.6% divergent, and they exhibit significant shifts in allozyme frequencies at five polymorphic loci. These concordant differences at multiple, independent, gene loci led Won et al.[29] to suggest that the southern lineage might warrant consideration as a distinct species, but additional samples and analyses of independent nuclear gene loci were needed to assess this hypothesis. With these goals, we examined mitochondrial COI and portions of five nuclear genes from mussels sampled along a 5600-km range from 13°N to 38°S latitude (Figure 1). Multilocus genotypes constructed from the DNA sequences were used to assess geographical structure and subdivision, to identify putative hybrids, and to test for gametic disequilibrium and recombination in the contact zone.



Expeditions conducted between 1990 and 2005 explored segments of the GAR, EPR and PAR (Figure 1) to obtain specimens for multi-species phylogeographic studies reviewed in [22]. Exploration in the Easter Microplate region was guided by shipboard coordinates and towed-video surveys reported in a PhD thesis [35]. Mussel samples were obtained from 14 localities, subsequently labeled with abbreviations that include the approximate latitude, e.g., N13, N11, … S31, S38 (Table 1). Samples collected with the human-occupied vehicle (HOV) Alvin were placed in an insulated ‘biobox’ containing ambient seawater at 2–4°C. Adult mussels (up to ~15 cm in length) were sampled with scoops, nets or directly with Alvin’s robotic manipulators. The S23 sample of juveniles (≤25 mm) was attached by byssal threads to a rock sampled for geological studies. Upon recovery of the submersible, all samples were stored temporarily in refrigerated seawater prior to dissection of gill and adductor muscle tissues that were frozen at −80°C or preserved in 95% ethanol.

Table 1 Sample coordinates for eastern Pacific Bathymodiolus mussels and sample sizes examined for six gene loci

DNA methods

DNA extraction and purification, general PCR conditions, amplicon purification and DNA sequencing used methods that were previously reported for mussels [37]. Six primer sets were used to amplify DNA targets from mussels (Table 2). All PCR products were diluted in 40 μl of sterile water and purified with a Multiscreen HTS PCR 96 vacuum manifold system (Millipore Corp. Billerica, MA). PCR products were sequenced bi-directionally on an ABI3130 sequencer with BigDye Terminator v3.1 (Life Technologies Corp., Carlsbad, CA) chemistry and primers used in PCR.

Table 2 PCR primers and methods for six gene loci in Bathymodiolus mussels

Molecular statistics

Computer programs used to format the data, obtain estimates of population genetic parameters, and conduct statistical tests are listed along with procedural references (Table 3). We used Phase v. 2.1.1 to ascertain alleles in individuals heterozygous for nuclear genes. Characterization of alleles in individuals heterozygous for insertions and deletions (indels) required cloning and sequencing of at least five clones per individual methods in [39]. Analyses of genic diversity, Hardy-Weinberg equilibrium, gametic disequilibrium, Isolation-with-Migration (IMa2), and Isolation-by-Distance (IBD) were limited to the exon portions of four nuclear genes that contained introns. None of the exon regions in these sequences exhibited stop-codons. The states of a 9-bp indel in the EF1α intron were treated as two-alleles. Tests for intra-genic recombination were conducted with Phase v. 2.1.1. Phylogenetic analyses of historical recombination within nuclear genes were conducted with Beagle on complete sequences, including introns and exons.

Table 3 Programs used to estimate genetic parameters and conduct statistical tests and procedures

Gene networks for the exon portions of each locus were constructed with Network Publisher v. Inferred haplotypes were re-labeled as discrete alleles (a1, a2, , an) and concatenated manually into a multilocus genotype (MLG) file (a1/a3, b2/b2, , g1/g1). Uninformative singleton alleles were lumped with their most closely related haplotype. The MLG file was used as an input file for Genepop v. 4.0.10. Exact tests for Hardy-Weinberg equilibrium (HWE) were conducted according to the method of Raymond and Rousset [49]. The MLG file was also used as an input file for Arlequin v. to conduct pairwise likelihood ratio tests for gametic disequilibrium and Mantel tests for Isolation-by-Distance (IBD). Sequential Bonferroni corrections were used to adjust α-levels [58].

The MLG file was also used to conduct assignment tests. Structure v. 2.3.3 analyses were conducted with and without prior information about population samples and admixture, and with correlated and uncorrelated allelic frequencies. The first 103 generations (burn-in) were discarded and 107 generations were iterated for 3 to10 times per estimate of K, the number of genotypic clusters. Assuming a uniform prior on K, we estimated a Bayes factor for each value of K after averaging the ln Pr(X/K) values between simulations. NewHybrids was used to assign the multilocus genotypes to putative parental, F1, F2, or backcross categories. We used default genotypic classes with no prior information on allelic frequencies and included uniform and Jefferey’s priors [59] for θ and π. Five separate analyses were conducted with 105 sweeps of burn-in and 107 sweeps of data collection. Threshold values (Tq) of 0.9 [56] and 0.75 [60] were assigned separately for the genotypic categories if q ≥ Tq, and they were left unassigned if q < Tq[61].

We used the Isolation-with-Migration method, IMa2[50, 62, 63] to estimate the following demographic parameters: (τ) the time of population splitting; (θA, θN, θS) the sizes of ancestral and descendant populations; and (mN, mS) immigration into the descendant populations. An HKY substitution model with back mutation was assumed. An inheritance scalar to adjust the relative expected effective population size was 0.25 for mitochondrial COI and 1.0 for nuclear loci. Analyses were conducted under a two-population model based on Structure results. The two-population model lumped the N13–S23 and GAR samples into a northern group, and the S31–S38 samples into a southern group. The S23 sample was treated as northern due to the preponderance of northern alleles. A second set of analyses was also run without the S23 population. Analyses involved at least 108 steps with the first 104 steps discarded as burn-in (−b) [64], with 50 attempts at chain swapping per step (−k), 50–80 chains (−n) with geometric heating (−f), and g1 and g2 values of 0.99 and 0.3, respectively. Analyses were conducted multiple times and the number of steps between saving genealogies (−d) was increased to 1000. Initially we set upper bounds for each parameter, as recommended in the IMa2 manual, then performed multiple runs using random number seeds with more appropriate bounds based on the initial output. We assessed the convergence of MCMC samplings of population parameters in multiple ways by allowing analyses to run ensuring effective sample sizes were a minimum of 100 for all parameters, examining output graphs for unimodality, ensuring genealogies were updating, and checking that multiple runs gave similar estimates. Multiple runs based on the same priors were then combined in L-mode and log-likelihood ratio tests were performed to test for isolation.


Adult mussels were abundant at most of the active vent fields between 13°N and 38°S latitude (Figure 1; Table 1). The S23 locality on the NW flank of the Easter Microplate was an exception. During three dive-days (bottom-time = 14.1 hours), we observed a small vent field littered with the shells of dead bivalves. Post-expedition inspection of the dive videos identified one living adult that was not sampled; however, we inadvertently collected 21 juvenile mussels (≤ 25 mm length) that were attached to the underside of a sulfide rock sampled for geological studies. The S23 site, located in a deep basaltic crevice, was recently more active. Sulfide rubble from dead chimneys was scattered about, but we found only one active vent chimney with a high temperature of 296°C. The dominant animals were typical of the fauna found at the periphery of active vents — e.g., serpulid annelids, caridean shrimp, anemones, enteropneusts, crinoids, brisingid sea stars, holothurians, sponges, and siphonophores [65]. Typical vent fauna, like alvinellid polychaetes were observed but not abundant. Two living Riftia pachyptila tubeworms were observed and one was collected [31]. Remnant tubes of Tevnia jerichonana occurred on sulfide rubble covering an area of diffuse low-temperature (≤ 8.5°C) venting, but no living worms were collected. Overall, the S23 locality exhibited hallmarks of a senescent vent field.

Geographical subdivision and admixture

Analyses were conducted with sequence polymorphisms in mitochondrial COI and the exon portions of four nuclear genes, SAHH, ANT, Cat and Col-1 (Table 4). A 9-bp indel in EF1α intron was treated as two discrete alleles. Haplotype networks for all six loci revealed differentiation between northern (N) and southern (S) mussel lineages (Figure 2). Distinct N and S haplotypes for COI were reciprocally monophyletic and separated by at least 20 substitutions. Though less divergent, predominantly N and S alleles segregated at all five nuclear loci (Figure 2). For each locus, the S23 sample (black pie slices, Figure 2) included N and S alleles. Latitudinal step-clines for five loci had common inflection points centered on the S23 locality (Figure 3a). The SAHH cline was more complex (Figure 3b). Frequencies of the northern alleles, SAHH*1 and *2, declined to zero at S17, counterbalanced by increases of SAHH*7 and *8. The SAHH*4 allele only occurred in the intermediate region, S17 and S23.

Table 4 Gene characteristics and GenBank accession numbers for six gene loci in eastern Pacific Bathymodiolus mussels
Figure 2

Gene networks for six loci in eastern Pacific Bathymodiolus. Numbers represent alleles used in the statistical analyses and black dots represent missing haplotypes. Colored pie slices denote regions, as indicated in Figure 1: NEPR = blue; GAR = green; SEPR = yellow; S23 = black S23; and PAR = red.

Figure 3

Latitudinal clines. (a) Step-clines for the predominant northern allele (haplotype: *1) at six loci toward southern direction from left to right (b) Allelic frequencies of SAHH alleles.

Indices of genic diversity differed among the genes (Table 5). Diversity within localities (h) was greatest for COI and smallest for ANT, which was defined by one SNP, and Col-1, which was not polymorphic at the NEPR localities. Among the nuclear loci, SAHH locus had the highest diversity and richness. Expected heterozygosity and allelic richness peaked in the intergrade region, between S17 and S31, for each of the nuclear genes (Figure 4a and b). Heterozygote deficiencies (F IS ) also peaked in the intergrade zone (Figure 4c). All statistically significant values for Fu’s Fs were negative (Table 5). They were most pronounced for COI, where nine of 13 values were significantly negative. Negative Fs values reflect significant excesses of rare alleles in the simple star-like nodes for COI (Figure 2).

Table 5 Genic diversity in eastern Pacific Bathymodiolus mussels
Figure 4

Statistical estimates of (a) gene diversity, (b) allelic richness, and (c) heterozygote deficiencies ( F IS ) for nuclear loci for mussel populations at each locality.

Except for S23, the other samples did not exhibit evidence for non-random mating. Sporadically elevated values for F IS in these samples were not significant after sequential Bonferroni correction of the α-levels. Similarly, sporadic gametic disequilibria were not significant following correction for multiple tests. In contrast, violations of Hardy-Weinberg equilibrium were pronounced in the S23 sample. Three of the five nuclear loci exhibited significant heterozygote deficiencies (Table 6). Cytonuclear disequilibrium (CND) was tested in the S23 sample after recoding the mitochondrial haplotypes (Figure 2) as N versus S “mitotypes”. CND was significant for all five pairwise gametic combinations. Two of the tests remained significant following sequential Bonferroni correction (Table 6). Significant gametic disequilibrium was evident in 10 pairwise comparisons of nuclear loci. Two contrasts remained significant following sequential Bonferroni correction (Table 6).

Table 6 Tests for random mating expectations in the S23 sample of eastern Pacific Bathymodiolus mussels

Pairwise sequence divergence within (d w ) and among (d a ) sample localities was estimated separately for each locus (Additional file 1: Table S1). The mean d w values were small for each locus: COI = 0.47%; SAHH = 0.46%; ANT = 0.01%; Cat = 0.38%; and Col-1 = 0.25%. The mean d a values were comparably low among localities within the northern (N13–S17) and southern (S31–S38) regions. Mean d a values between samples from the different regions were typically much greater than the d w values. The admixed S23 sample varied among loci in its affinities to the northern and southern regions.

These regional differences were also revealed by multilocus genotypic assignment methods. The most probable partitioning of the total sample resulted in K = 3 (Additional file 1: Table S2). The northern and southern partitions were self-evident from the previous analyses. An intermediate partition composed of the S17 and S23 was influenced by the SAHH locus, perhaps due to recombination or incomplete lineage sorting. Exclusion of SAHH from the analysis resulted in an estimate of K = 2, with S17 and S23 nested in the northern partition. The NewHybrids analysis (Figure 5) identified four of the S23 individuals with high posterior probabilities (PP > 0.95) as hybrids, and other potential hybrids (PP < 0.75) also were noted. Due to the limited number of loci, insufficient power existed to unequivocally assign the hybrid individuals to first-generation (F1), second-generation (F2) or specific backcross categories [60]. Several putative hybrids also occurred in the S32 sample. The S17 sample had a single individual that was assigned to the southern partition (PP > 0.95).

Figure 5

NewHybrids analyses of multi-locus nuclear genotypes in eastern Pacific Bathymodiolus . Sample localities are ordered along the Y-axis. Vertical lines mark 0.75 and 0.90 values of posterior probabilities for assignments to the following genotypic categories: blue = northern parental; red = southern parental; black = F1 hybrid; white = F2 hybrid; light blue = northern backcross hybrid; and pink = southern backcross hybrid.

Isolation-with-migration (IMa2) analyses were conducted with exon regions of genes that exhibited no evidence for recombination. Analyses were conducted under a two-population model, defined by the previous Structure and NewHybrid analyses, a northern (N) population pooled the N13–S23 and GAR samples, and a southern (S) population pooled the S31–S38 samples. Maximum likelihood estimates (MLEs) for the effective population size scaled by mutation rate (θ = 2Nμ) were four-times larger for the northern population, as defined (Table 7, Additional file 2: Figure S1). Pooling would limit our interpretation of these estimates, however, if significant structuring occurred within a partition. To assess this potential problem, we excluded the S23 sample from the IMa2 analysis. This exclusion reduced the estimates of effective population size by one-half (Table 7), but θN still remained four-times greater than θS. Decisions about pooling greatly influenced estimates of the immigration vectors (M i  = m i /μ). Immigration into the northern population was greater than the reverse (MN > MS) if the S23 sample was included in the northern population. Log likelihood ratio tests on migration rates showed exclusion of S23 from the analysis resulted in no significant bidirectional migration. Estimates of ancestral population sizes (θA) and population splitting time (τ = ) were unresolved, regardless of pooling. A sharp peak near zero in the highest posterior density (HPD) curve for τ is indicative of recent divergence; however, an escalating tail of the curve prevented an estimate for population splitting time from the multi-locus data.

Table 7 Maximum likelihood estimates (MLE) of population parameters from IM a 2 analyses from B. thermophilus cx. Lower (L) and upper (U) 95% HPD (highest posterior density) are given


We reconstructed recombination histories from the complete sequence data (with introns) for four nuclear loci. Recombination could not be tested for ANT because it exhibited only one polymorphic site. Coalescence analyses conducted with Beagle[45] identified large numbers of recombination events in SAHH, Cat, Col-1 and EF1a (indicated with stars in Additional file 3: Figure S2). The northern and southern lineages were color-coded, as in Figure 1; however, homoplasy among the lineages obscured the identification of within- versus between-lineage recombination events. Nonetheless, the temporal distribution of recombination events suggests that the parental lineages may have had a long history of parapatric contact.


Previous studies reported that B. thermophilus populations distributed between 13°N and 17°S exhibited evidence for isolation-by-distance [29, 36]. We re-tested this hypothesis with the present mitochondrial and nuclear polymorphisms. The combined pairwise F ST values for nuclear loci were not correlated with geographic distance (r = 0.008; P < 0.567, 1000 permutations) (Additional file 4: Figure S3). In contrast, ϕ ST values estimated from mitochondrial sequences exhibited a significant correlation (r = 0.447; P < 0.019, 1000 permutations). Nonetheless, the ϕ ST correlation resulted from contrasts involving the S17 sample. Excluding S17 from the analysis eliminated the correlation. Based on the Structure and IMa2 analyses, the S17 sample appears to be influenced by introgressed southern alleles.

Taxonomic implications

Based on allozyme and mitochondrial DNA differences, Won et al.[29] suggested that the northern and southern B. thermophilus lineages might warrant recognition as distinct species. To assess this suggestion, we compared them with Bathymodiolus azoricus Cosel et al.[66] and B. puteoserpentis Cosel et al.[67], named geminate-species that hybridize along the Mid Atlantic Ridge. Although individuals of the MAR species can be discriminated with 95% confidence based on shell dimensions alone [24, 68], the B. thermophilus lineages exhibit no comparable differences [29]. Also, sequence divergence for mitochondrial genes (COI and ND4) is greater between two MAR lineages, but differentiation at the nuclear loci is consistently less (Table 8). The SAHH, Col-1 and Ef1α sequences examined in both species complexes are less between the MAR species (Table 8). Different numbers of allozyme loci were examined in the two complexes; so, we could only compare mean F ST values for the subsets of polymorphic genes. The proportion of allozyme variance that exists between the MAR species is half that of the B. thermophilus lineages. To summarize, shell dimensions and mitochondrial sequences indicate that the MAR species are more divergent than B. thermophilus lineages, but the opposite is true for nuclear genes. Because adaptively neutral divergence in protein-coding nuclear genes is generally believed to accumulate much more slowly than in mitochondrial genes, the consensus of evidence suggests that the EPR lineages split earlier than the MAR species.

Table 8 Genetic divergence between northern and southern lineages of B. thermophilus and B. antarcticus (N/S), and between B. azoricus and B puteoserpentis ( az/pu ) from the Mid-Atlantic Ridge

Assessing the species-status of allopatric and parapatric populations is an old problem in biology. Here we adopt the “metapopulation species concept” which equates species with separately evolving segments of the broader metapopulation [70]. In this sense, the widely distributed metapopulation known as B. thermophilus sensu latu is clearly partitioned into separately evolving metapopulation segments. Mixed nuclear genotypes in the hybrid zone reveal that northern and southern lineages can still hybridize. Yet, the capability to hybridize does not prove they are conspecific, as many recognized plant and animal species hybridize at contact zones or in disturbed environments [5]. Consequently, we hereby elevate the southern lineage to species status, and restrict the application of B. thermophilus sensu stricto to populations occurring to the north of the Easter Microplate hybrid zone.


Family MYTILIDAE Rafinesque, 1815

Subfamily Bathymodiolinae Kenk & Wilson, 1985

Type genus Bathymodiolus Kenk & Wilson, 1985

Bathymodiolus antarcticus Johnson & Vrijenhoek new species

(Figure 6)

Figure 6

Bathymodiolus antarcticus n. sp. paratype (USNM: 1196520). A. lateral view of right valve. B. Posterior view of hinge, left and right valves.

Bathymodiolus aff. thermophilus[29, 37, 7173]Bathymodiolus cf. thermophilus[74]

Type material. — Holotype: sampled at 37.78°S/110.90°W on the Pacific-Antarctic Ridge, 2200 meters depth during Alvin dive A4088 (Table 1). Preserved in 10% sea-water buffered formalin, then transferred to ethanol (USNM 1196518). Paratypes: sampled at 31°51S/112°W, on the East Pacific Rise, 2336 meters depth during Alvin dive A4092 (Table 1). Preserved in ethanol (6 USNM 1196519; 6 SMNH 8414–8415). Species registered with Zoobank:

Diagnosis. — This species differs from B. thermophilus s.s. by the following combination of character states of its mitochondrial COI barcode: 126C, 160A, 217C, 276 T, 282 G/A, 330 T, 390A, 403 T, 504 T, 519 T, 547C, 558C, 570A, and 603C. It also differs diagnostically for exon regions the SAHH gene. Morphological differences are unknown.

Distribution. — Known from hydrothermal vent localities between 31°09S on the East Pacific Rise to 37°47S on the Pacific-Antarctic Ridge.

Etymology. — Latin adjective antarcticus (ant’arc.ti.cus), meaning opposite the North Star, and referring to its distribution on the Pacific-Antarctic Ridge.

Bathymodiolus thermophilus sensu stricto, Kenk & Wilson, 1985, as described

Type material. — Holotype: USNM 803661, preserved in 70% ethanol. Collected 20 January 1979 by Ellis and Ballard on Alvin dive 879 from Mussel Bed (0 47.89’N; 86 9.21’W) at 2495 m depth on the Galápagos Rift.

Revised diagnosis. — This species differs from B. antarcticus n. sp. by the following combination of character states of its mitochondrial COI barcode: 126 T, 160 G, 217 T, 276C, 282 T, 330C, 390C, 403C, 504 G/A, 519C/A, 547 T, 558 T, 570 G, and Type material. — Holotype: United States National Museum (USNM 803661), preserved in 70% ethanol. Collected 20 January 1979 by Ellis and Ballard on Alvin dive 879.

Distribution. — Known from hydrothermal vents on the Galápagos Rift between 0.80°N/86.15 and 0.80°N/86.122 W, and from the East Pacific Rise (EPR) between 12.80°N/103.93°W and 21.57°S/114.4 W.


This is the second report of a mid-ocean hybrid zone involving hydrothermal vent mussels [16, 29]. As in the previous case involving the Mid-Atlantic Ridge (MAR) mussels, the present example from the eastern Pacific Ocean also exhibits the properties of a Tension Zone: (1) coincident clines at multiple loci; (2) heterozygote deficiencies in the contact zone; and (3) significant gametic disequilibrium among nuclear and cytoplasmic genes (Figures 2, 3, 4, Table 6). Tension Zones are maintained by hybrid unfitness (i.e. endogenous selection) and recurrent immigration from the parental regions [1, 75, 76]. The abrupt step-clines at six independent gene loci suggest that introgression has been limited to adjacent localities.

We could not ascertain from the present genetic data whether the excess proportions of parental genotypes in the hybrid zone were products of selection against hybrids, or consequences of admixture (i.e. Wahlund effects). Comparison of juvenile and adult cohorts might help to identify ongoing selection against hybrids e.g., [77], but we are not aware of future efforts to explore this remote region with submersible vehicles. Simple admixture of immigrant recruits from the adjacent parental regions would generate the observed heterozygote deficiencies, gametic phase disequilibrium [78], and putative hybrid genotypes in the S23 sample (Figure 5). Although the S23 vent field appeared to be in a senescent stage of the hydrothermal cycle a single active vent chimney supported small populations of chemosynthetically dependent animals. We observed one living adult mussel (unsampled) while reviewing the dive videos from 14 hrs of Alvin ‘bottom-time’. Shell debris from recently deceased bivalves littered the floor of the crevice that contained the active chimney. Other vents hosting more robust mussel populations might exist nearby, but they remain undiscovered.

“Sweepstakes Reproductive Success” (SRS) as observed in some near-shore bivalves and other organisms reviewed in [42], might be invoked to explain the heterozygote deficiencies, but the hybrid zone juveniles were not the progeny of few females. The sample of 21 juveniles contained 11 distinct mitochondrial haplotypes (Table 5). Haplotypic diversity (h = 0.819) was not different from the mean diversity (h = 0.807 ± 0.077) for all other samples, which exhibited no evidence for deviations from random-mating expectations. The planktotrophic larvae of Bathymodiolus mussels appear to mix thoroughly before they settle at vents. Their capacity for long-distance dispersal is consistent with a lack of geographical subdivision among EPR (N13–S17) and GAR populations spanning thousands of kilometers [79] and this study. A previous report of subdivision across the Equator is an artifact of sampling gaps cf. [21, 36]. Nonetheless, broad-scale homogeneity could coincide with small-scale patchiness in genotypic frequencies under the SRS model [42]. Hypervariable genetic markers might reveal localized patchiness due to episodic recruitment of cohorts from a limited number of parents in different source populations. For example, AFLP variation provides preliminary evidence for small-scale patchiness in the eastern Pacific vent annelid Riftia pachyptila[80], despite broad-scale homogeneity of nuclear and mitochondrial DNA markers [31]. In contrast, microsatellite DNAs provided no evidence for small-scale patchiness in the western Pacific vent snail, Ifremeria nautileii[81].

Divergence and intergradation

Reciprocally monophyletic differences for mitochondrial haplotypes and step-clines in allozyme frequencies suggest a relatively long history of isolation between B. antarcticus n. sp. and B. thermophilus s.s. Based on borrowed molecular clock rates for COI (1–2% divergence per million years), we estimated that the northern and southern mitochondrial lineages split 2.1– 4.3 Myr ago, which temporally overlaps a span of dates encompassing orogeny of the Easter Microplate, 2.5–5.3 Myr ago [82, 83]. An ongoing study that uses fossil data to calibrate major splitting events in the subfamily Bathymodiolinae estimates a similar time for the most recent common ancestor of the mitochondrial lineages (J. Lorien, personal communication). Nonetheless, mitochondrial ‘gene trees’ do not necessarily equal ‘species trees’ [84]. Different mutation rates among nuclear and mitochondrial genes, the amounts of ancestral polymorphism in various genes, and differing intensities of selection on independent characters can greatly complicate relationships between gene and organismic trees. Disagreements between these trees are further complicated by introgressive hybridization and reticulation [85]. Thus, it is not surprising that IMa2 analysis failed to resolve a splitting time (τ) from the multi-locus data for these mussels.

Whether intergradation between the northern and southern lineages results from primary or secondary processes remains uncertain. A coalescent analysis of the MAR mussel species provided evidence for “an initial period of allopatric differentiation during which recombination was blocked between lineages” followed by “a long history of gene flow” [25]. In contrast, our analysis of past recombination in the eastern Pacific mussels is consistent with a long history of parapatric contact. However, these contacts probably are episodic in this dynamic region of the global ridge system. Chaotic local extinctions and colonizations of vents as tectonic activities move up and down a ridge axis can create transient contact zones that allow gene flow [86]. The observed recombination history of these mussel lineages is consistent with episodic isolation and re-contact as tectonic conditions vary temporally.

Isolating mechanisms

Density troughs (regions of low population size) can function as powerful traps to stabilize hybrid zones and limit introgression [8789]. The MAR hybrid zone coincides with a density trough at ‘Broken Spur’ [16, 24], a geographically intermediate vent field that is relatively inhospitable for mussels due to precipitation of toxic polymetalic sulfides believed to interfere with respiration [23]. Mussel population densities were also low in the Easter microplate region during our 1999 and 2005 expeditions, but the causes are different. First, Won et al.[29] hypothesized that strong cross-axis currents will remove planktotrophic mussel larvae from the ridge axis, limiting local self-recruitment. These superfast-spreading centers have highly inflated axial calderas and essentially no lateral walls to constrain axial circulation [90]. The predominantly cross-axis currents in this region drive buoyant hydrothermal plumes to the west [91] presumably carrying larvae that rise above the axial calderas. Secondly, Coykendall et al.[31] hypothesized that vent habitats along the southern EPR are subject to high rates of disturbance. The habitats are frequently “repaved” with lava flows that drive local extinctions and create new vents [74]. We explored a number of hydrothermally active vent fields between 18.40°and 21.73°S latitude during our January 1999 expedition (Table 1), but none of them supported robust mussel populations. Soon after, Van Dover [92] visited the ‘Animal Farm’ vent field at 18.60°S and found it “clearly in a waning stage of the hydrothermal cycle’ (p. 142). Venting water temperatures were only slightly elevated above ambient and most of the mussels were dead. The living individuals were smaller on average compared to dense mussel populations just to the north (Table 1). The BIOSPEEDO expedition in 2004 sampled mussels at 18.40°and 21.42°S, but population densities were not reported [36]. Demographic fluctuations result in variance-effective population sizes that might be very small for some taxa, leading to losses of genetic variation [22]. For example, southern EPR populations of Riftia pachyptila are nearly devoid of nuclear and mitochondrial sequence variation [31], but other annelids are not similarly affected [30, 36]. Together, cross-axis currents and metapopulation processes might reduce the genetically effective population size of mussel populations in this region, impeding along-axis dispersal and limiting introgression across the Easter Microplate.


Although B. thermophilus s.s. and B. antarcticus n. sp. cannot be distinguished with shell measurements, they differ in having reciprocally monophyletic mitochondrial lineages. Polymorphisms at nuclear loci (allozymes and DNA sequences) are not completely diagnostic but the two species comprise discrete genotypic clusters with predominantly non-overlapping geographical distributions. Some shared alleles might be remnants of incomplete lineage sorting between the metapopulation lineages, but these alleles are not widely distributed geographically, as expected for ancestral polymorphisms. Instead, they only penetrate localities that flank the Easter Microplate, a pattern expected for introgressive hybridization or disruptive selection along an environmental gradient [1].

Divergence between the two species probably initiated with orogeny of the Easter Microplate, 2.5–5.3 Myr ago. Superfast tectonic spreading rates in this dynamic part of the global mid-ocean ridge system result in frequent habitat destruction, leaving episodic gaps in the distribution of many vent-restricted species. Prolonged periods of partial to complete isolation would facilitate the accumulation of genetic differences in nuclear and mitochondrial genes exposed to genetic drift and potentially contrasting selection regimes. We suggest two hypotheses about possible adaptions that might result in disruptive selection across the Easter Microplate boundary. First, we hypothesize that novel life history traits emerged to accommodate different current regimes in the northern and southern regions. Northern B. thermophilus s.s. populations are subjected to a circulation pattern dominated by mesoscale eddies and along-axis currents [20, 9395]. There, larval buoyancy and duration probably evolved to favor the retention of the planktotrophic prodissoconch stages. South of the Easter Microplate, B. antarcticus n. sp. is subjected to strong cross-axis vectors forced by the Antarctic Circumpolar Current [29]. Larval buoyancies and durations that function in the northern current regime would result in larval removal in this contrasting regime [96]. Though we currently have no developmental data on the larval life histories of these mussels, this hypothesis is testable. Early settling post-larvae retain their developmental history in prodissoconch sizes and morphologies that can provide significant clues regarding larval life history [97].

Secondly, we hypothesize that coadaptation between the mussel hosts and associated thiotrophic endosymbionts might result in disruptive selection. B. thermophilus s.s. and B. antarcticus n. sp. differ with respect to the symbiont strains they host in specialized cells (bacteriocytes) contained mostly in gill tissues [98]. Though we cannot at this time exclude some degree of vertical transmission, the endosymbionts are mostly acquired locally from the environment in which mussel larval settle [99]. An ongoing multi-gene analysis of symbiont diversity has only found the B. thermophilus s.s. symbiont strain in the hybrid zone mussels (Y-J Won, unpublished data), even in host individuals with predominantly southern genotypes. Subtle interactions involving these intracellular microbes might have affected evolution in mitochondrial and nuclear genes of the hosts. According to the “Red King” model, infectious horizontal symbionts and their hosts are expected to exhibit decelerated evolutionary rates, because the participants are subject to constraints that limit any form of change [100, 101]. Episodic range expansions into the contact zone will result in secondary contacts that re-establish opportunities for gene flow and introgression in the hosts, and possible symbiont transfers. Yet, conflicts between divergent larval strategies, symbiont coadaptation, and other unknown processes would determine the strength of selection against intergrades. Disruptive processes involving the mussel hosts and their obligate endosymbionts might contribute to partial isolation and incipient speciation of these mussels.

Author’s contributions

SBJ conducted most of the DNA sequencing and statistical analyses with the assistance of YJW and JBH. RCV led oceanographic expeditions that with the help of SBJ and YJW collected mussel specimens. SBJ, RCV and YJW drafted the manuscript. All authors read and approved the final manuscript.


  1. 1.

    Barton NH, Hewitt GM: Analysis of hybrid zones. Annu Rev Ecol Syst. 1985, 16: 113-148. 10.1146/

    Google Scholar 

  2. 2.

    Harrison RG: Hybrid Zones and the Evolutionary Process. 1993, Oxford: Oxford Press

    Google Scholar 

  3. 3.

    Hewitt GM: Hybrid zones—natural laboratories for evolutionary studies. Trends Ecol Evol. 1988, 3: 158-167. 10.1016/0169-5347(88)90033-X.

    CAS  PubMed  Google Scholar 

  4. 4.

    Moore WS, Price TJ: Nature of selection in the northern flicker hybrid zone and its implications for speciation theory. Hybrid Zones and the Evolutionary Process. Edited by: Harrison RG. 1993, Oxford: Oxford Press, 196-225.

    Google Scholar 

  5. 5.

    Arnold ML: Natural Hybridization and Evolution. 1997, Oxford, UK: Oxford University Press

    Google Scholar 

  6. 6.

    Endler J: Geographic Variation, Speciation and Clines. 1977, Princeton (NJ): Princeton University Press

    Google Scholar 

  7. 7.

    Remington CL: Suture-zones of hybrid interaction between recently joined biotas. Evolutionary Biology. Edited by: Dobzhansky T, Hecht MK, Steere WC. 1968, New York: Plenum Press, 321-428.

    Google Scholar 

  8. 8.

    Hewitt G: The genetic legacy of the Quaternary ice ages. Nature. 2000, 405 (6789): 907-913. 10.1038/35016000.

    CAS  PubMed  Google Scholar 

  9. 9.

    Palumbi SR: Genetic divergence, reproductive isolation, and marine speciation. Annu Rev Ecol Syst. 1994, 25: 547-572. 10.1146/

    Google Scholar 

  10. 10.

    Gardner JPA: Hybridization in the sea. Adv Mar Biol. 1997, 31: 1-78.

    Google Scholar 

  11. 11.

    Koehn RK: The genetics and taxonomy of species in the genus Mytilus. Aquaculture. 1991, 94: 125-145. 10.1016/0044-8486(91)90114-M.

    Google Scholar 

  12. 12.

    Rawson PD, Hilbish TJ: Evolutionary relationships among male and female mitochondrial DNA lineages in the Mylilus edulis species complex. Mol Biol Evol. 1995, 12: 893-901.

    CAS  PubMed  Google Scholar 

  13. 13.

    McDonald JH, Seed R, Koehn RK: Allozymes and morphometric characters of three species of Mytilus in the Northern and Southern Hemispheres. Mar Biol. 1991, 111: 323-333. 10.1007/BF01319403.

    Google Scholar 

  14. 14.

    Rawson PD, Agrawal V, Hilbish TJ: Hybridization between the blue mussels Mytilus galloprovincialis and M. trossulus along the pacific coast of North America: evidence for limited introgression. Mar Biol. 1999, 134: 201-211. 10.1007/s002270050538.

    Google Scholar 

  15. 15.

    Gosling EM: The genetics of Mytilus. The mussel, Mytilus: ecology, physiology, genetics and culture, developments in aquaculture and fisheries science. Edited by: Gosling EM. 1992, Amsterdam: Elsevier, 309-382.

    Google Scholar 

  16. 16.

    O’Mullan GD, Maas PAY, Lutz RA, Vrijenhoek RC: A hybrid zone between hydrothermal vent mussels (Bivalvia: Mytilidae) from the Mid-Atlantic Ridge. Mol Ecol. 2001, 10: 2819-2831.

    PubMed  Google Scholar 

  17. 17.

    Van Dover CL, German CR, Speer KG, Parson LM, Vrijenhoek RC: Evolution and biogeography of deep-sea vent and seep invertebrates. Science. 2002, 295 (5558): 1253-1257. 10.1126/science.1067361.

    CAS  PubMed  Google Scholar 

  18. 18.

    Van Dover CL: The ecology of Deep-Sea Hydrothermal Vents. 2000, Princeton: Princeton University Press

    Google Scholar 

  19. 19.

    Thomson RE, Mihály SF, Rabinovich AB, McDuff RE, Veirs SR, Stahr FR: Constrained circulation at Endeavour ridge facilitates colonization by vent larvae. Nature. 2003, 24: 545-549.

    Google Scholar 

  20. 20.

    Mullineaux LS, Adams DK, Mills SW, Beaulieu SE: Larvae from afar colonize deep-sea hydrothermal vents after a catastrophic eruption. Proc Natl Acad Sci. 2010, 107 (17): 7829-7834. 10.1073/pnas.0913187107.

    CAS  PubMed Central  PubMed  Google Scholar 

  21. 21.

    Audzijonyte A, Vrijenhoek R: When gaps really are gaps: Statistical phylogeography of hydrothermal vent invertebrates. Evolution. 2010, 64 (8): 2369-2384.

    PubMed  Google Scholar 

  22. 22.

    Vrijenhoek RC: Genetic diversity and connectivity of deep-sea hydrothermal vent metapopulations. Mol Ecol. 2010, 19: 4391-4411. 10.1111/j.1365-294X.2010.04789.x.

    PubMed  Google Scholar 

  23. 23.

    Desbruyères D, Almeida A, Comtet T, Khripounoff A, Le Bris N, Sarradin P-M, Segonzac M: A review of the distribution of hydrothermal vent communities along the northern Mid-Atlantic ridge: Dispersal vs. environmental controls. Hydrobiologia. 2000, 440 (1/3): 201-216. 10.1023/A:1004175211848.

    Google Scholar 

  24. 24.

    Won Y-J, Hallam SJ, O’Mullan GD, Vrijenhoek RC: Cytonuclear disequilibrium in a hybrid zone involving deep-sea hydrothermal vent mussels of the genus Bathymodiolus. Mol Ecol. 2003, 12 (11): 3185-3190. 10.1046/j.1365-294X.2003.01974.x.

    CAS  PubMed  Google Scholar 

  25. 25.

    Faure B, Jollivet D, Tanguy A, Bonhomme F, Bierne N: Speciation in the deep sea: Multi-locus analysis of divergence and gene flow between two hybridizing species of hydrothermal vent mussels. PLoS One. 2009, 4 (8): e6485-10.1371/journal.pone.0006485.

    PubMed Central  PubMed  Google Scholar 

  26. 26.

    Guinot D, Hurtado LA: Two new species of hydrothermal vent crabs of the genus Bythograea from the southern east pacific rise and from the galapagos rift (Crustacea Decapoda Brachyura Bythograeidae). Comptes Rendus Biologies. 2003, 326: 423-439. 10.1016/S1631-0691(03)00126-4.

    PubMed  Google Scholar 

  27. 27.

    Guinot D, Hurtado LA, Vrijenhoek RC: New genus and species of brachyuran crab from the southern east pacific rise (Crustacea Decapoda Brachyura Bythograeidae). Comptes Rendus Biologies. 2002, 325 (11): 1119-1128. 10.1016/S1631-0691(02)01529-9.

    Google Scholar 

  28. 28.

    Johnson S, Warén A, Vrijenhoek RC: DNA barcoding of Lepetodrilus limpets reveals cryptic species. J Shellfish Res. 2008, 27 (1): 43-51. 10.2983/0730-8000(2008)27[43:DBOLLR]2.0.CO;2.

    Google Scholar 

  29. 29.

    Won Y-J, Young CR, Lutz RA, Vrijenhoek RC: Dispersal barriers and isolation among deep-sea mussel populations (Mytilidae: Bathymodiolus) from eastern Pacific hydrothermal vents. Mol Ecol. 2003, 12 (1): 169-184.

    CAS  PubMed  Google Scholar 

  30. 30.

    Hurtado LA, Lutz RA, Vrijenhoek RC: Distinct patterns of genetic differentiation among annelids of eastern Pacific hydrothermal vents. Mol Ecol. 2004, 13 (9): 2603-2615. 10.1111/j.1365-294X.2004.02287.x.

    CAS  PubMed  Google Scholar 

  31. 31.

    Coykendall DK, Johnson SB, Karl SA, Lutz RA, Vrijenhoek RC: Genetic diversity and demographic instability in Riftia pachyptila tubeworms from eastern Pacific hydrothermal vents. BMC Evol Biol. 2011, 11: 96-10.1186/1471-2148-11-96.

    PubMed Central  PubMed  Google Scholar 

  32. 32.

    Stoffers P, Worthington T, Petersen S, Hannington M, Tükay M, Ackermand D, Borowski C, Danert S, Fretdorff S, Haase K, et al: Cruise Report SONNE 157, FOUNDATION 3: Magmatic and Hydrothermal Processes at a Spreading Axis influenced by a Hotspot: the Pacific-Antarctic Ridge and Off-Axis Seamounts near 37°S, Berichte - Reports Institut für Geowissenschafen Nr 17. 2001, Kiel, Germany: Christrian-Albrechts Universität zu Kiel

    Google Scholar 

  33. 33.

    Stoffers P, Worthington T, Hekinian R, Petersen S, Hannington M, Tükay M, the SIO 157 Shipboard Scientific Party: Silicic volcanism and hydrothermal activity documented at Pacific-Antarctic Ridge. EOS, Transactions of the American Geophysical Union. 2002, 83 (28): 32-38.

    Google Scholar 

  34. 34.

    Kenk VC, Wilson BR: A new mussel (Bivalvia, Mytilidae) from hydrothermal vents in the Galapagos Rift Zone. Malacologia. 1985, 26 (3180): 253-271.

    Google Scholar 

  35. 35.

    Holler G: Ph.D. Dissertation. Erschnelnungsformer hydrothermler Activität am Ostpazifischen Rücken zwischen 6° und 30° Süd. 1993, Göttingen: Georg August Universität

    Google Scholar 

  36. 36.

    Plouviez S, Shank TM, Faure B, Daguin-Thiebaut C, Viard F, Lallier FH, Jollivet D: Comparative phylogeography among hydrothermal vent species along the East Pacific Rise reveals vicariant processes and population expansion in the South. Mol Ecol. 2009, 18 (18): 3903-3917. 10.1111/j.1365-294X.2009.04325.x.

    CAS  PubMed  Google Scholar 

  37. 37.

    Génio L, Johnson SB, Vrijenhoek RC, Cunha MR, Tyler PA, Kiel S, Little CTS: New record of “Bathymodiolus” mauritanicus Cosel,2002 from the Gulf of Cadiz (NE Atlantic) mud volcanoes. J Shellfish Res. 2008, 27 (1): 53-61. 10.2983/0730-8000(2008)27[53:NROBMC]2.0.CO;2.

    Google Scholar 

  38. 38.

    Audzijonyte A, Vrijenhoek RC: Three nuclear genes for phylogenetic, SNP and population genetic studies of molluscs and other invertebrates. Mol Ecol Resour. 2010, 10: 200-204. 10.1111/j.1755-0998.2009.02737.x.

    CAS  PubMed  Google Scholar 

  39. 39.

    Johnson SB, Young CR, Jones WJ, Warén A, Vrijenhoek RC: Migration, isolation, and speciation of hydrothermal vent limpets (Gastropoda; Lepetodrilidae) across the Blanco Transform Fault. Biol Bull. 2006, 210 (2): 140-157. 10.2307/4134603.

    PubMed  Google Scholar 

  40. 40.

    Excoffier L, Lischer HEL: Arlequin suite ver. 3.5: A new series of programs to perform population genetics analyses under Linux and Windows. Mol Ecol Resour. 2010, 10: 564-567. 10.1111/j.1755-0998.2010.02847.x.

    PubMed  Google Scholar 

  41. 41.

    Fu YX: Statistical tests of neutrality of mutations against population growth, hitchhiking and background selection. Genetics. 1997, 147 (2): 915-925.

    CAS  PubMed Central  PubMed  Google Scholar 

  42. 42.

    Slatkin M: Linkage disequilibrium in growing and stable populations. Genetics. 1994, 137 (1): 331-

    CAS  PubMed Central  PubMed  Google Scholar 

  43. 43.

    Slatkin M, Barton NH: A comparison of three indirect methods for estimating average levels of gene flow. Evolution. 1989, 43 (7): 1349-1368. 10.2307/2409452.

    Google Scholar 

  44. 44.

    Smouse PE, Long JC: Matrix correlation analysis in anthropology and genetics. Yearbook of Physical Anthropology. 1992, 35: 187-213. 10.1002/ajpa.1330350608.

    Google Scholar 

  45. 45.

    Lyngso RB, Song YS, Hein J: Workshop on Algorithms in Bioinformatics. Minimum recombination histories by branch and bound. 2005, 239-250.

    Google Scholar 

  46. 46.

    Drummond A, Rambaut A: BEAST. v. 1.0.3. 2003, Oxford, UK: University of Oxford

    Google Scholar 

  47. 47.

    Asmussen MA, Basten CJ: Sampling theory for cytonuclear disequilibria. Genetics. 1994, 138: 1351-1363.

    CAS  PubMed Central  PubMed  Google Scholar 

  48. 48.

    Rousset F: Genepop’007: a complete reimplementation of the Genepop software for Windows and Linux. Mol Ecol Resour. 2008, 8: 103-106. 10.1111/j.1471-8286.2007.01931.x.

    PubMed  Google Scholar 

  49. 49.

    Raymond M, Rousset F: GENEPOP (Ver. 1.2): population genetics software for exact tests and ecumenicism. J Hered. 1995, 86: 248-249.

    Google Scholar 

  50. 50.

    Hey J: Isolation with migration models for more than two populations. Mol Biol Evol. 2010, 27: 905-920. 10.1093/molbev/msp296.

    CAS  PubMed Central  PubMed  Google Scholar 

  51. 51.

    Drummond A, Ashton B, Buxton S, Cheung M, Cooper A, Duran C, Field M, Heled J, Kearse M, Markowitz S, et al: Geneious v5.4. 2011, Available from,

    Google Scholar 

  52. 52.

    Bandelt H-J, Forster P, Röhl A: Median-joining networks for inferring intraspecific phylogenies. Mol Biol Evol. 1999, 16: 37-48. 10.1093/oxfordjournals.molbev.a026036.

    CAS  PubMed  Google Scholar 

  53. 53.

    Anderson EC, Thompson EA: A model-based method for identifying species hybrids using multilocus genetic data. Genetics. 2002, 160: 1217-1229.

    CAS  PubMed Central  PubMed  Google Scholar 

  54. 54.

    Stephens M, Smith NJ, Donnelly P: A new statistical method for haplotype reconstruction from population data. Am J Hum Genet. 2001, 68 (4): 978-989. 10.1086/319501.

    CAS  PubMed Central  PubMed  Google Scholar 

  55. 55.

    Stephens M, Donnelly P: A comparison of Bayesian methods for haplotype reconstruction from population genotype data. Am J Hum Genet. 2003, 73: 1162-1169. 10.1086/379378.

    CAS  PubMed Central  PubMed  Google Scholar 

  56. 56.

    Pritchard JK, Stephens M, Donnelly P: Inference of population structure using multilocus genotype data. Genetics. 2000, 155: 945-959.

    CAS  PubMed Central  PubMed  Google Scholar 

  57. 57.

    Hubisz MJ, Falush D, Stephens M, Pritchard JK: Inferring weak population structure with the assistance of sample group information. Mol Ecol Resour. 2009, 9: 1322-1332. 10.1111/j.1755-0998.2009.02591.x.

    PubMed Central  PubMed  Google Scholar 

  58. 58.

    Rice WR: Analyzing tables of statistical tests. Evolution. 1989, 43: 223-225. 10.2307/2409177.

    Google Scholar 

  59. 59.

    Gelman A: Bayesian model building by pure thought. Some principals and examples. Statista Sinica. 1996, 6: 215-232.

    Google Scholar 

  60. 60.

    Vähä J, Primmer CR: Efficiency of model-based Bayesian methods for detecting hybrid individuals under different hybridization scenarios and with different numbers of loci. Mol Ecol. 2006, 15: 63-72.

    PubMed  Google Scholar 

  61. 61.

    Olivera R, Godinho R, Randi E, Ferrand N, Alves PC: Molecular analysis of hybridisation between wild and domestic cats (Felis silvestris) in Portugal: implications for conservation. Conserv Genet. 2008, 9: 1-11.

    Google Scholar 

  62. 62.

    Hey J, Nielsen R: Multilocus methods for estimating population sizes, migration rates and divergence time, with applications to the divergence of Drosophila pseudoobscura and D. persimilis. Genetics. 2004, 167 (2): 747-760. 10.1534/genetics.103.024182.

    CAS  PubMed Central  PubMed  Google Scholar 

  63. 63.

    Hey J, Nielsen R: Integration within the Felsenstein equation for improved Markov chain Monte Carlo methods in population genetics. Proc Natl Acad Sci USA. 2007, 104 (8): 2785-2790. 10.1073/pnas.0611164104.

    CAS  PubMed Central  PubMed  Google Scholar 

  64. 64.

    Nielsen R, Wakeley J: Distinguishing migration from isolation: a Markov chain Monte Carlo approach. Genetics. 2001, 158: 885-896.

    CAS  PubMed Central  PubMed  Google Scholar 

  65. 65.

    Desbruyères D, Segonzac M, Bright M: Handbook of Deep-Sea Hydrothermal Vent Fauna. Brest. 2006, France: Zugleich Kataoge der Oberöstereichischen Lanmuseen, N.S. 43-

    Google Scholar 

  66. 66.

    von Cosel R, Comtet T, Krylova E: Bathymodiolus (Bivalvia: Mytilidae) from hydrothermal vents on the Azores triple junction and the Logatchev hydrothermal field, Mid-Atlantic ridge. Veliger. 1999, 42 (3): 218-248.

    Google Scholar 

  67. 67.

    von Cosel R, Métivier B, Hashimoto J: Three new species of Bathymodiolus (Bivalvia: Mytilidae) from hydrothermal vents in the Lau basin and the Fiji basin, western Pacific, and the Snake Pit area Mid-Atlantic ridge. Veliger. 1994, 37: 374-392.

    Google Scholar 

  68. 68.

    Maas PAY, O’Mullan GD, Lutz RA, Vrijenhoek RC: Genetic and morphometric characterization of mussels (Bivalvia: Mytilidae) from Mid-Atlantic hydrothermal vents. Biol Bull. 1999, 196: 265-272. 10.2307/1542951.

    CAS  PubMed  Google Scholar 

  69. 69.

    Jones WJ, Won YJ, Maas PAY, Smith PJ, Lutz RA, Vrijenhoek RC: Evolution of habitat use by deep-sea mussels. Mar Biol. 2006, 148 (4): 841-851. 10.1007/s00227-005-0115-1.

    Google Scholar 

  70. 70.

    de Queiroz K: Species concepts and species delimitation. Syst Biol. 2007, 56 (6): 879-886. 10.1080/10635150701701083.

    PubMed  Google Scholar 

  71. 71.

    Duperron S: Topics in Geobiology. vol. 33. The Diversity of Deep-Sea Mussels and Their Bacterial Symbioses. 2010, Dordrecht: Springer, Netherlands, 137-167.

    Google Scholar 

  72. 72.

    Petersen JM, Zielinski FU, Pape T, Seifert R, Moraru C, Amann R, Hourdez S, Girguis PR, Wankel SD, Barbe V, et al: Hydrogen is an energy source for hydrothermal vent symbioses. Nature. 2011, 476 (7359): 176-180. 10.1038/nature10325.

    CAS  PubMed  Google Scholar 

  73. 73.

    van der Heijden K, Petersen JM, Dubilier N, Borowski C: Genetic connectivity between north and south Mid-Atlantic ridge chemosynthetic bivalves and their symbionts. PLoS One. 2012, 7 (7): e39994-10.1371/journal.pone.0039994.

    CAS  PubMed Central  PubMed  Google Scholar 

  74. 74.

    Hey R, Massoth G, Vrijenhoek R, Rona P, Lupton J, Butterfield D: Hydrothermal vent geology and biology at Earth’s fastest spreading rates. Mar Geophys Res. 2006, 27: 137-153. 10.1007/s11001-005-1887-x.

    Google Scholar 

  75. 75.

    Key KHL: The concept of stasipatric speciation. Syst Zool. 1968, 17: 14-22. 10.2307/2412391.

    Google Scholar 

  76. 76.

    Bazykin AD: Hypothetical mechanism of speciation. Evolution. 1969, 23: 685-687. 10.2307/2406862.

    Google Scholar 

  77. 77.

    Dowling TE, Moore WS: Evidence for selection against hybrids in the family Cyprinidae (genus Notropis). Evolution. 1985, 39: 152-158. 10.2307/2408524.

    Google Scholar 

  78. 78.

    Nei M, Li WH: Linkage disequilibrium in subdivided populations. Genetics. 1973, 75: 213-219.

    CAS  PubMed Central  PubMed  Google Scholar 

  79. 79.

    Craddock C, Hoeh WR, Lutz RA, Vrijenhoek RC: Extensive gene flow in the deep-sea hydrothermal vent mytilid Bathymodiolus thermophilus. Mar Biol. 1995, 124: 137-146. 10.1007/BF00349155.

    Google Scholar 

  80. 80.

    Shank TM, Halanych KM: Toward a mechanistic understanding of larval dispersal: insights from genomic fingerprinting of the deep-sea hydrothermal vent tubeworm Riftia pachyptila. Marine Ecology. 2007, 28 (1): 25-35. 10.1111/j.1439-0485.2007.00146.x.

    CAS  Google Scholar 

  81. 81.

    Thaler AD, Zelnio K, Saleu W, Schultz TF, Carlsson J, Cunningham C, Vrijenhoek RC, Dover CLV: The spatial scale of genetic subdivision in populations of Ifremeria nautilei, a hydrothermal-vent gastropod from the southwest Pacific. BMC Evol Biol. 2012, 11 (11): 372-

    Google Scholar 

  82. 82.

    Rusby RI, Searle RC: A history of the Easter Microplate, 5.25 Ma to present. J Geophys Res. 1995, 100 (B7): 12,617-12,640.

    Google Scholar 

  83. 83.

    Naar DF, Hey RN: Tectonic evolution of the Easter Microplate. J Geophys Res. 1991, 96: 7961-7993. 10.1029/90JB02398.

    Google Scholar 

  84. 84.

    Pamilo P, Nei M: Relationships between gene trees and species trees. Evolution. 1988, 5: 568-583.

    CAS  Google Scholar 

  85. 85.

    Machado CA, Hey J: The causes of phylogenetic conflict in a classic Drosophila species group. Proceedings of the Royal Society B: Biological Sciences. 2003, 270: 1193-1202. 10.1098/rspb.2003.2333.

    CAS  PubMed Central  PubMed  Google Scholar 

  86. 86.

    Jollivet D, Chevaldonné P, Planque B: Hydrothermal-vent alvinellid polychaete dispersal in the Eastern Pacific. 2. A metapopulation model based on habitat shifts. Evolution. 1999, 53: 1143-1156. 10.2307/2640818.

    Google Scholar 

  87. 87.

    Barton NH: The dynamics of hybrid zones. Heredity. 1979, 43: 341-359. 10.1038/hdy.1979.87.

    Google Scholar 

  88. 88.

    Hewitt G, Barton N: The structure and maintenance of hybrid zones as exemplified by Podisma pedestris. Insect cytogenetics. Edited by: Blackman RL, Hewitt GM, Ashburner M. 1980, Oxford: Blackwell Scientific, 149-169.

    Google Scholar 

  89. 89.

    Harrison R: Hybrid zones: windows on evolutionary process. Oxf Surv Evol Biol. 1990, 7: 69-128.

    Google Scholar 

  90. 90.

    MacDonald KC: Linkages between faulting, volcanism, hydrothermal activity and segmentation on fast spreading centers. Faulting and Magmatism at Mid-Ocean Ridges. Edited by: Buck W, Delaney P, Karson J, Lagabrielle Y. 1998, Washington, DC: American Geophysical Union, 27-58.

    Google Scholar 

  91. 91.

    Lupton JE: Hydrothermal helium plumes in the Pacific Ocean. Journal of Geophysical Research-Oceans. 1998, 103: 15853-15868. 10.1029/98JC00146.

    CAS  Google Scholar 

  92. 92.

    Van Dover CL: Community structure of mussel beds at deep-sea hydrothermal vents. MEPS. 2002, 230: 137-158.

    Google Scholar 

  93. 93.

    Marsh AG, Mullineaux LS, Young CM, Manahan DT: Larval dispersal potential of the tubeworm Riftia pachyptila at deep-sea hydrothermal vents. Nature. 2001, 411 (3): 77-80. 10.1038/35075063.

    CAS  PubMed  Google Scholar 

  94. 94.

    Adams DK, Flierl GR: Modeled interactions of mesoscale eddies with the East Pacific Rise: Implications for larval dispersal. Deep-Sea Research Part I. 2010, 57 (10): 1163-1176. 10.1016/j.dsr.2010.06.009.

    Google Scholar 

  95. 95.

    Adams DK, McGillicuddy DJ, Zamudio L, Thurnherr AM, Liang X, Rouxel O, German CR, Mullineaux LS: Surface-generated mesoscale eddies transport deep-sea products from hydrothermal vents. Science. 2011, 332: 580-583. 10.1126/science.1201066.

    CAS  PubMed  Google Scholar 

  96. 96.

    Mullineaux LS, Speer KG, Thurnherr AM, Maltrud ME, Vangriesheim A: Implications of cross-axis flow for larval dispersal along mid-ocean ridges. Cahiers De Biologie Marine. 2002, 43 (3–4): 281-284.

    Google Scholar 

  97. 97.

    Lutz RA, Jablonski D, Rhoads DC, Turner RD: Larval dispersal of a deep-sea hydrothermal vent bivalve from the Galapagos Rift. Mar Biol. 1980, 57 (1100): 127-133.

    Google Scholar 

  98. 98.

    Won Y-J, Jones WJ, Vrijenhoek RC: Absence of co-speciation between deep-sea mytilids and their thiotrophic endosymbionts. J Shellfish Res. 2008, 27 (1): 129-138. 10.2983/0730-8000(2008)27[129:AOCBDM]2.0.CO;2.

    Google Scholar 

  99. 99.

    Won Y-J, Hallam SJ, O’Mullan GD, Pan IL, Buck KR, Vrijenhoek RC: Environmental acquisition of thiotrophic endosymbionts by deep-sea mussels of the genus Bathymodiolus. Appl Environ Microbiol. 2003, 69 (11): 6785-6792. 10.1128/AEM.69.11.6785-6792.2003.

    CAS  PubMed Central  PubMed  Google Scholar 

  100. 100.

    Bergstrom CT, Lachmann M: The Red King effect: When the slowest runner wins the coevolutionary race. Proc Natl Acad Sci USA. 2003, 100 (2): 593-598. 10.1073/pnas.0134966100.

    CAS  PubMed Central  PubMed  Google Scholar 

  101. 101.

    Sachs JL, Essenberg CJ, Turcotte MM: New paradigms for the evolution of beneficial infections. Trends Ecol Evol. 2011, 26 (4): 202-209. 10.1016/j.tree.2011.01.010.

    PubMed  Google Scholar 

Download references


The project was conducted with the expert help of captains, crews and pilots assigned to the R/V Atlantis, R/V Atlantis II, and research submersible Alvin (Woods Hole Oceanographic Institute). We thank Dr. Richard A. Lutz (Rutgers University) for contributing the N9 samples and participating in numerous oceanographic expeditions. W. Joe Jones, Phil Hoos, C. Robbie Young, and E. Anderson contributed to sampling, sequencing, and statistical analyses. Funding was provided by NSF grants OCE-9910799 and OCE-0241613 (to RCV), the David and Lucile Packard Foundation (to MBARI) and the National Research Foundation of Korea Grant funded by the Korean Government NRF-2011-013-C00063 (to YJW).

Author information



Corresponding author

Correspondence to Shannon B Johnson.

Additional information

Competing interest

The authors declare that they have no competing interest.

Electronic supplementary material

Table S1.

Additional file 1: Mean percentage of sequence divergence (K2P corrected) within (d w , boldface on diagonal) and between (d a off-diagonal) eastern Pacific Bathymodiolus mussels. Pairwise comparisons involving the S23 sample are highlighted in blue. The d values between northern and southern regions are highlighted in gray. Table S2. Simulation results for the Structure analysis of eastern Pacific Bathymodiolus mussels under correlated and uncorrelated allelic frequency models. The number of possible genotypic clusters, K, varied from one to eleven. The natural log of the probability of the data for a given value of K was averaged across the number of simulations (runs) per K. Bayes factors were estimated based on these averages. (DOCX 136 KB)

Figure S1.

Additional file 2: IMA2 posterior distributions for demographic parameters in B. thermophilus under a two-population model: a) N13–S23 vs. S31–S38; b) N13–N21 vs. S31–S38. (EPS 1 MB)

Figure S2.

Additional file 3: Recombination analyses for the loci (a) EF1α, 12 recombination events, (b) Col-1, 17 recombination events, (c) SAHH, 11 recombination events, and (d) Cat, 10 recombination events. Open circles represent mutations and stars represent recombination events. Colors on OTUs are blue (NEPR), green (GAR), yellow (SEPR), and red (PAR). (EPS 3 MB)

Figure S3.

Additional file 4: Analysis of isolation-by-distance in B. thermophilus populations from north of the Easter Microplate. Open triangles are φ ST values estimated from mitochondrial sequences in the GAR and EPR (N13–S21) samples. Filled circles are F ST values estimated from the combined nuclear loci in the GAR, and EPR (N13–S17) samples. (EPS 3 MB)

Authors’ original submitted files for images

Rights and permissions

This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Reprints and Permissions

About this article

Cite this article

Johnson, S.B., Won, YJ., Harvey, J.B. et al. A hybrid zone between Bathymodiolusmussel lineages from eastern Pacific hydrothermal vents. BMC Evol Biol 13, 21 (2013).

Download citation


  • Bathymodiolus antarcticus n. sp
  • Hybridization
  • Recombination
  • Linkage disequilibrium
  • Deep-sea
  • Hydrothermal vent
  • Bathymodiolus thermophilus