Figures
Abstract
Background
Bivalves are very ancient and successful conchiferan mollusks (both in terms of species number and geographical distribution). Despite their importance in marine biota, their deep phylogenetic relationships were scarcely investigated from a molecular perspective, whereas much valuable work has been done on taxonomy, as well as phylogeny, of lower taxa.
Methodology/Principal Findings
Here we present a class-level bivalve phylogeny with a broad sample of 122 ingroup taxa, using four mitochondrial markers (MT-RNR1, MT-RNR2, MT-CO1, MT-CYB). Rigorous techniques have been exploited to set up the dataset, analyze phylogenetic signal, and infer a single final tree. In this study, we show the basal position of Opponobranchia to all Autobranchia, as well as of Palaeoheterodonta to the remaining Autobranchia, which we here propose to call Amarsipobranchia. Anomalodesmata were retrieved as monophyletic and basal to (Heterodonta + Pteriomorphia).
Conclusions/Significance
Bivalve morphological characters were traced onto the phylogenetic trees obtained from the molecular analysis; our analysis suggests that eulamellibranch gills and heterodont hinge are ancestral characters for all Autobranchia. This conclusion would entail a re-evaluation of bivalve symplesiomorphies.
Citation: Plazzi F, Ceregato A, Taviani M, Passamonti M (2011) A Molecular Phylogeny of Bivalve Mollusks: Ancient Radiations and Divergences as Revealed by Mitochondrial Genes. PLoS ONE 6(11): e27147. https://doi.org/10.1371/journal.pone.0027147
Editor: Carles Lalueza-Fox, Institut de Biologia Evolutiva-Universitat Pompeu Fabra, Spain
Received: June 23, 2011; Accepted: October 11, 2011; Published: November 1, 2011
Copyright: © 2011 Plazzi et al. 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 author and source are credited.
Funding: This work was supported by the University and Research Italian Ministry (MIUR PRIN07, grant number 2007NSHJL8_002 to MP) and the the “Canziani Bequest” fund (University of Bologna, grant number A.31.CANZELSEW to MP). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.
Competing interests: The authors have declared that no competing interests exist.
Introduction
The impressive biological success of bivalves is a perfect example of evolutionary potentials embedded in a clear-cut modification of the already successful molluscan body plan. The major distinguishing features of the bivalve mollusk are the peculiar architecture of their shell, the lateral compression (and general reduction) of the foot and the complete loss of the radula. First bivalves appeared in the Cambrian period [1]–[7], but the oldest genera are poorly known and tough to interpret. Pojeta [5] retained only five Cambrian genera as actually bivalves: two of them are rather studied (Fordilla and Pojetaia), whereas Arhouriella, Camya and Tuarangia are much more controversial [1], [4]. Different factors triggered the Ordovician bivalve radiation: the evolution of a feeding gill [4], the presence of a byssus gland in the adult [4], [8], the development of an infaunal way of life [6], [7] linked to the so-called “Cambrian substrate revolution” (see [7]; and reference therein). Actually, from rare, pedal-feeding surface-dwellers of the early Cambrian, all the principal clades of extant bivalves evolved in the middle Ordovician [3], [6], [7] in a “two-pulse process” [4], [7]. Since then, bivalve phylogeny was a flourishing of branches on a wide tree.
Most probably today's protobranchs resemble those first bivalve species, showing a well-developed foot and true molluscan ctenidia, principally devoted to gas exchange [9], [10]. Nevertheless protobranchs developed long palp proboscides to bring food to the mouth, that were probably lacking in the earliest forms Pojetaia and Fordilla [3], [7]. The modification of gills for filter feeding, with the consequent reduction and loss of palp proboscides, the gain of byssus, allowing epifaunal life, the mantle margin fusion, with the emergence of siphons, triggered bivalves' adaptive radiations along geological eras [3], [8], [11], [12].
Nowadays, bivalves are arranged into four big clades, which are given the status of subclasses. Protobranchs were divided in two clusters. Species belonging to order Nuculoida are considered among the most primitive bivalves and were included in the subclass Palaeotaxodonta by Newell [13]. The order Solemyoida was described as unrelated to nuculoids, and was included in the subclass Lipodonta (sensu [2]), but more recently other authors preferred to merge both orders in the subclass Protobranchia [3], [14]–[16]; indeed, molecular analyses supported a sister group relationship between the two orders [17], [18]. Furthermore, the superfamily Nuculanoidea was removed from Protobranchia [19]–[22], and Giribet [12] proposed the name Opponobranchia referring to the subclass-rank clade Nuculoida + Solemyoida.
Sister group of the Opponobranchia is the Autobranchia ( = Autolamellibranchiata sensu [12]), i.e. bivalves with modified ctenidia, without palp proboscides, generally filibranch or eulamellibranch. Some authors, like Waller [16], treated Autobranchia as a subclass itself. Following the most widely accepted taxonomy, however, three subclasses, substantially identical to the definition in Newell [13], belong to Autobranchia: Heterodonta, Palaeoheterodonta, and Pteriomorphia.
Relationships within Autobranchia are still contentious: many studies retrieved a monophyletic clade called Heteroconchia, joining Palaeoheterodonta and Heterodonta [6], [12], [16], [19], [21], [23]. Conversely, several phylogenetic analyses resulted in a close relationship between Pteriomorphia and Heterodonta [2], [20], [22], [24]–[26].
Eventually, species belonging to subclass Anomalodesmata (order Pholadomyoida) are generally eulamellibranch, siphonate burrowers, which developed some remarkable adaptations: some of them are septibranch and deep-water carnivors. Formerly ascribed to their own subclass [13], [27], they are currently considered as a basal, monophyletic clade among Heterodonta ([12], [19], [20], [28]–[31]; but see [22]).
As mentioned, only a handful of comprehensive molecular phylogenetic studies have been released to date. After some pioneering analyses [25], [32], [33], and the extensive effort of Campbell [34], most recent deep phylogenies concentrate on single subclasses: Pteriomorphia [17], [35], Anomalodesmata [29], [30], and particularly Heterodonta, the most diverse group [31], [36]–[38]. Direct optimization [39] was used for wide scale phylogenetic reconstructions, as Giribet and Wheeler [19] and Giribet and Distel [20] assembled a thorough total evidence matrix, the broadest ever assembled on bivalve evolution.
Finally, our previous study [22] was the first attempt to infer a complete evolutionary tree of the class with a robust, two-steps phylogenetic analysis. The aim of that work was to develop a sound pipeline to approach bivalve molecular phylogenetics: the present paper follows this pipeline by adding more bivalve taxa, to obtain an in-depth survey of the evolutionary tree of Bivalvia. This study represents the biggest dataset of bivalve mollusks to date, which has been characterized by four mitochondrial genes. Thanks to this improved dataset, we will address all those issues that were not possible to discuss in detail in our previous paper [22], mainly focusing on deep relationships linking bivalve subclasses.
Results
Sequence data
A total of 60 sequences from 29 species were obtained for this study and deposited in GenBank under Accession Numbers JF496737-JF496786. Sequences of MT-RNR1 (12s), MT-RNR2 (16s), MT-CO1 (cox1) and MT-CYB (cytb) were 19, 9, 17, and 15, respectively. Details of the concatenated alignment (Dataset S1) are listed in Table 1; final alignment (as well as the phylogenetic tree) was deposited in TreeBASE under the Study Accession URL http://purl.org/phylo/treebase/phylows/study/TB2:S11320. After removal of ambiguously aligned positions and related indel characters, 2260 nucleotides and 735 indels were left for phylogenetic analyses, for a total of 2995 characters. The complete dataset includes 436 sequences from 122 bivalves and five outgroup species. Interestingly, we found four protein-coding gene (PCG) sequences (Neopycnodonte cochlear, Spondylus gaederopus, and Talochlamys multistriata for MT-CO1; Laevicardium crassum for MT-CYB) where single-site gaps have to be included to obtain a correct alignment. In our previous work, we noted the same for Hyotissa hyotis and Barbatia cfr. setigera cytochrome b sequences [22]. The alignment, both at the nucleotide and amminoacid level, is otherwise good, therefore it is unlikely we are facing a NUMT (i.e., a mitochondrial pseudogene; [40]), inasmuch that no NUMTs have been reported for bivalves to date [41], [42]. The simplest explanation involves a sequencing error and cannot be dismissed. Anyway, even if we do not have empirical data on this account, single nucleotide indels in apparently functional mitochondrial genes–MT-CYB being one of them–have been reported and discussed elsewhere [43–45; and reference therein], which in turn would surely deserve further investigation. For phylogenetic purposes, we inserted missing data instead of single-site gaps whenever they mapped in a region of the gene included in the alignment.
Evaluating phylogenetic signal
Phylogenetic Representativeness test aims to measure the degree of representativeness of a sample with respect to the group it should represent in a phylogenetic analysis (Fig. 1; see [46]). The measured Average Taxonomic Distinctness (AvTD) of our sample of 86 bivalve genera fell within the 95% confidence interval of AvTD computed from 100 random subsample of the same dimension. However, the Variance in Taxonomic Distinctness (VarTD) was clearly higher than its 95% confidence interval (Fig. 1A). Moreover, the AvTD of our sample was within the range of 95% lower confidence limit yielded by the shuffling test (Fig. 1B). Most probably, the little sampling among Anomalodesmata taxa (which are indeed hard to obtain) is the main reason of the border-line AvTD and the high VarTD we found.
A, Average Taxonomic Distinctness (AvTD) and Variation in Taxonomic Distinctness (VarTD) computed for the sample used for this study. AvTD is plotted on left axis: the circle represents the value obtained from the present sample, whereas continuous lines indicate the lower 95% confidence limit, the maximum value for that sample dimension (thick lines), and the mean AvTD (thin line). VarTD is plotted on the right axis: the diamond represents the value obtained from the present sample, whereas dotted lines indicate the minimum value for that sample dimension, the upper 95% confidence limit (thick lines), and the mean VarTD (thin line). B, shuffling test with 100 randomly shuffled master lists (see text for details). Mean VarTD (thin dotted lines), upper 95% VarTD confidence limit (upper thick dotted lines), lower 95% AvTD confidence limit (lower thick continue lines), and mean AvTD (thin continue lines) are shown as the 95% confidence intervals across the replicates. Axes, circle, and diamond as above.
Uncorrected pairwise distances plotted on Maximum Likelihood (ML) pairwise distances showed some saturation in substitutions along our dataset (Fig. 2), which is expected given the depth of this phylogeny. Accounting for multiple hits in aligned sequences, Maximum Likelihood distances are greater than the corresponding uncorrected distances. Such a degree of saturation deserves some caution in analyzing this dataset, by implementing non-trivial evolutionary models and carefully assessing results' statistical support: as a matter of facts, saturation plots are compatible with a high-level phylogenetic reconstruction (see, f. i., [47]), as uncorrected distances only partially level off on Maximum Likelihood distances and a statistically significant positive trend is present in all gene partitions (Fig. 2 and Table S1), with the expected exception of third codon positions. Therefore, best models for bivalve mitochondrial phylogenetic inference will have to discard these sites, or analyze them in a codon-based context, thus confirming our previously proposed pipeline for bivalves phylogenetic analysis [22].
Uncorrected (“p-”) distances plotted on Maximum Likelihood distances to estimate saturation in our dataset for each of six partitions separately. See text for details on partition names. Units of x-axis in substitutions/position. Linear interpolation is also shown, whose details are listed in Table S1.
Likelihood Mapping (LM; Fig. 3) allowed the estimation of the amount of signal present in our data; first of all, 1000 random quartets were drawn without constraints. They are evenly (P>0.05) distributed in the simplex, but only 8.6% of them do fall into the star-like tree area, while 85.2% map near one of the three vertices, indicating that in most cases a topology is strongly favored over alternative hypotheses. The concatenated alignment as well as single genes and partitions (data not shown) were examined, and in all cases a preferred topology was individuated (Fig. 3). 8 out of 13 analyses indicated the unrooted topology ((Palaeoheterodonta + Heterodonta) + (Anomalodesmata + Pteriomorphia)) as the most supported; the second most supported topology was ((Palaeoheterodonta + Anomalodesmata) + (Heterodonta + Pteriomorphia)), which was retrieved for 3 partitions. The same Maximum Likelihood model was used as in the saturation analysis; results from all 13 analyses were significantly different from the null hypothesis (P<0.005) and more than 60% of them pointed towards the same backbone tree–therefore, we are confident that this approach is able to overcome multiple hits-linked flaws.
Each analysis was performed on 1,000 random quartets; the left simplex shows point distribution; the central one the subdivision among the three corners; the right one the subdivision among Voronoi cells [96], [115]. A, Likelihood Mapping for the concatenated alignment without grouping. B, Likelihood Mapping for the concatenated alignment with Opponobranchia excluded and remaining taxa subdivided into Palaeoheterodonta (a), Anomalodesmata (b), Heterodonta (c), and Pteriomorphia (d). The three possible topologies are shown at vertices.
Preliminary phylogenetic reconstructions
Neighbornet networks of the complete alignment were produced for single genes and for the concatenated alignment, based on uncorrected and LogDet distances. All networks are essentially similar, varying only in the positions of some taxa, like Lucinella, Loripes, Cuspidaria, Nuculana, Astarte, and Cardita. Figure 4 shows the LogDet neighbornet network for the complete alignment: all genera and families are retrieved as well-defined clades, with the exception of mytilids and Chlamys. Although the network is less clearly tree-like in deep relationships, some sharp signal is present for major groups, like Palaeoheterodonta (the Unionidae are very well distinct in all networks). Opponobranchia often cluster together with Haliotis and other outgroups. The position of anomalodesmatans is unstable among different genes and distance methods: under LogDet model, they cluster together next to part of Mytilidae (Lithophaga lithophaga, Modiolula phaseolina, Modiolus sp.), whereas under the uncorrected method Cuspidaria is found close to Loripes and Lucinella between Opponobranchia and Heterodonta, and Pandora and Thracia are found in a star-like region of the tree with Cardita, Astarte and Nuculana. These last three genera are found among pteriomorph species under the LogDet model. Single-gene networks are generally consistent with this topology, with local resolution decreasing in some part of the graph. Long branches were found in some single-gene networks only (mostly those of ribosomal markers), whereas, in the concatenated alignment, this was only the case for the scaphopod outgroup Siphonodentalium lobatum.
Results of molecular evolution models for each partition are extensively listed in Table S2. For ML analysis, the model selected for the partition all was implemented with PAUP* [48]. The heuristic search with 150 bootstrap replicates yielded a consensus tree with generally high support values (Fig. S1). Bivalves did not cluster in a supported monophyletic clade: the scaphopod Siphonodentalium lobatum was found to be the sister group of a polytomy with Katharina, Haliotis, Thais, genus Nucula, Solemya, and all remaining bivalves (the Autobranchia), whose monophyly has a bootstrap proportion (BP) value of 65. The first split separates Palaeoheterodonta (BP = 100) and a broad assemblage of species belonging to Heterodonta, and Pteriomorphia: neither was retrieved as monophyletic, nor were anomalodesmatans. Families and genera are generally monophyletic, with some exceptions, like Arcidae.
Model-decision and Bayesian Inference
As expected, results from Akaike Information Criterion (AIC) and Bayes Factor (BF) tests (Table 2 and S3) were straightforward in distinguishing between 4by4 and codon models: all partitioning schemes implementing the M3 codon model (i.e., p14–p17) outperformed those implementing the classical 4by4 analysis (i.e., p01–p13). The AIC test selected p14 as the best model for our dataset, with an Estimated Marginal Likelihood (EML) of −118,729.10, whereas BF selected p17 (EML = −118,205.79). It has to be noted that codon-based analyses are extremely demanding in terms of computational power: therefore, as detailed in Methods section, we used single MC3 analyses with half generations with respect to 4by4 models. Four of such analyses were run to estimate convergence within and among runs, and parameters and trees were finally summarized given the convergence evidence. In all cases, we could compute final statistics and consensus tree from 2 runs, with the exception of p17, where we could use only 3,416 generations from a single run, which is an order of magnitude lower than we did for models p14–p16. Therefore, the preference of BF for model p17 could be an effect of the low and different sample size of this specific run; moreover, AIC should be more conservative whenever these concerns are present, in that it accounts for overparametrization in the model by penalizing a high number of free parameters K (see [22]; and reference therein for further details). In conclusion, we regard to p14 as the most supported tree of our study, and it is shown in Figure 5.
Shown is p14 tree, computed partitioning our dataset into ribosomal and protein coding genes; these were analyzed using the M3 codon model (see text for details). Values at the nodes are Posterior Probabilities (PP); nodes were collapsed if PP<0.95. Color code as follows: violet, Opponobranchia; blue, Palaeoheterodonta; brown, Anomalodesmata; green, Heterodonta; red, Pteriomorphia. Astartoidea and Carditoidea, previously classified as heterodonts and included hereafter among pteriomorphians, are shown in orange and marked with asterisks.
Five monophyletic clusters with Posterior Probabililty (PP) equal to 1 were obtained, corresponding to the five traditional subclasses. Opponobranchia (here Nucula and Solemya) was retrieved as monophyletic and basal to Autobranchia, whose topology was found to be (Palaeoheterodonta + (Anomalodesmata + (Heterodonta + Pteriomorphia))). Nodes are robustly supported along the whole tree, and most have PP = 1.00.
Subclass Palaeoheterodonta is divided into two extant orders, Trigonioida and Unionoida. Cristaria plicata is basal to remaining palaeoheterodonts in our tree. A polytomy separates Lanceolaria grayana, the genus Unio, the genus Anodonta, the cluster Pyganodon + Psaudanodonta, and a cluster with remaining unionids with Alathyria jacksoni (family Hyriidae). Therefore, family Unionidae is paraphyletic because of Alathyria, subfamily Anodontinae is paraphyletic as well, because of Cristaria, and subfamily Unioninae is polyphyletic. On the other hand, subfamily Ambleminae is monophyletic, and 3 out of 4 tribes are represented in our tree: only the tribe Lampsilini is represented with more than one genus (Epioblasma, Lampsilis, Venustaconcha), and it is monophyletic. No specimen from order Trigonioida was included in this study.
Only one order, Pholadomyoida, belongs to subclass Anomalodesmata. Although the subclass is monophyletic, internal relationships are unresolved. However, Thracia and Pandora cluster together as sister group of Cuspidaria with PP = 0.85 in p14 and this relationship is present in all trees, being also supported with PPs>0.95 in some of them. Therefore, a signal, albeit weak, is present for the monophyly of Pandoroidea (suborder Pholadomyina).
Superfamily Lucinoidea (Loripes + Lucinella) is basal to all remaining heterodonts. These (excluding Astarte + Cardita, see below) are arranged as a polytomy separating two big clusters and two small clades, (Abra + Donax) and (Ensis + Sinonovacula). The first big cluster can be described as ((Dreissenoidea + Myoidea) + (Mactroidea + (Corbiculoidea + Glossoidea + Veneroidea))). Genera Dreissena and Mactra are monophyletic, as are families Mactridae and Veneridae. Relationships within venerids are well resolved, and subfamily Tapetinae and Meretricinae are monophyletic; only subfamily Chioninae is not monophyletic, because of the sister group relationship between Clausinella and Venus. The second big cluster can be described as (Hiatelloidea + Cardioidea). Subfamily Tridacninae is basal to a polytomy with Fraginae (Lunulicardia + Corculum), Laevicardiinae (Laevicardium), and a cluster with Cardiinae (Acanthocardia) and Cerastodermatiinae (Cerastoderma).
Two clades are basal to the core of Pteriomorphia. The first is the cluster (Astarte + Cardita), which is generally ascribed to Heterodonta as composed by superfamilies Astartoidea and Carditoidea. The second is the monophyletic family of Mytilidae, which is divided in two sister groups: on one side, (Lithophaginae + Modiolinae); on the other side, (Mytilaster sp. + (Crenellinae + Mytilinae)). Relationships within the core of Pteriomorphia are well resolved: it is subdivided into three clusters, one of them represented by Nuculana commutata alone, which was formerly ascribed to Protobranchia. The second cluster has Anomia as basal to Limoidea and Pectinoidea (both monophyletic). Genus Acesta is monophyletic and sister group of the cluster (Lima pacifica galapagensis + (Lima sp. + Limaria sp.)), therefore genus Lima is paraphyletic. Spondylus (family Spondylidae) and Parvamussium (family Propeamussiidae) are basal to a heterogeneous clade of intermingled Pectinidae and Propeamussiidae (Adamussium, Amusium), where many lower taxa are found as polyphyletic: Chlamydinae, Pectininae, genus Chlamys. Conversely, subfamily Patinopectininae is monophyletic due to the sister group relationship between Patinopecten and Mizuhopecten. The third cluster is composed by order Arcida as sister group of (Pteriida + Ostreina). With minor exceptions, like the polyphyly of Barbatia, and the paraphyly of Pteriida, Pteriidae, Arcidae, and Arcinae, most taxa were recovered as monophyletic: namely, we could retrieve as highly supported clusters subfamilies Pycnodonteinae, Ostreinae, families Gryphaeidae, Ostreidae, superfamilies Ostreoidea, Arcoidea, suborders Ostreina, Pteriina, Arcina, and order Arcida.
Morphological features
Six morphological characters were traced and optimized on p14 tree: gill type, shell microstructure [13], hinge [19], gill cilia [49]–[51], stomach type [52], and labial palps [10]. Parsimony reconstructions of ancestral states are shown in Figure 6; ML was also implemented for all those characters where multiple states were not used, and results were in complete agreement with parsimony.
Optimization of six major morphological characters on bivalve phylogeny as retrieved in this work. Each tree shows the parsimony reconstruction of ancestral state given the p14 topology and a matrix of morphological characters compiled following [13] and [19]; see text for more detail. A, gill grade; B, hinge; C, gill cilia; D, stomach type; E, labial palps; F, shell microstructure.
Discussion
Phylogenetic signal
All the evidence from our dataset points towards the conclusion that phylogenetic signal is available through the combined use of these four mitochondrial markers, but it is not trivial to detect it correctly.
This is expected because of the phylogenetic depth of this study: bivalves arose 530 million years ago, in the earliest Cambrian ([3], [22], [53]; and reference therein). The saturation profile (see Fig. 2) is compatible with the old age of the class; expectedly, repeated substitution events at the same site (multiple hits) were found. Nevertheless, given the significantly positive trend that was found in most case (and especially for MT-RNR2, MT-CO1, and MT-CYB; see Fig. 2 and Table S1), irrespective of gene properties, we may conclude that the use of complex evolutionary models should account for the saturation occurred in the four analyzed genes. This is further corroborated by Likelihood Mapping and neighbornet networks (see Fig. 3 and 4). Evidence of monophyly were found for all the major groups of bivalve systematics, with special reference to pteriomorph radiation and minor exceptions. Some groups appear to be particularly well-defined in our dataset, like Ostreidae, Pectinidae, Unionidae, and Veneridae. Moreover, the method of Likelihood Mapping implements precise and statistically tested evolutionary models, which are able to account for multiple hits along genes and for rate mutation heterogeneity. Indeed, the use of Likelihood Mapping simplex could finally ascertain the presence of strong phylogenetic signal in our dataset (see Fig. 4A) and also the first emergence of one or two preferred topologies (see Fig. 4B).
For this reason, we think that the use of mitochondrial genes should not be discarded a priori to address deep phylogenies, but rather they have to be rigorously tested before the analysis. It is well known that the deeper is an evolutionary relationship, the more refined a technique must be to unveil and exploit it. This is especially the case for the general backbone of bivalve tree, which had to be targeted with advanced Bayesian inference. In this study, as in our previous preliminary analysis [22], selected models tend to merge homogeneous markers in single partitions (i.e. ribosomal genes on one side and PCGs on the other), indicating that this is most likely the best trade-off between a detailed, realistic model and overparametrization.
Bivalve phylogeny
The p14 Bayesian tree (Fig. 5) is very well resolved; the high number of taxa it included allow to address many evolutionary issues about bivalves.
Opponobranchia was confirmed as separate to all Autobranchia; the reduced length of branches leading to Nuculoidea and Solemyiodea constitutes an evidence that these species tend to retain most ancestral characters, as widely hypothesized (see, f.i., [3], [9], [54]; and reference therein).
Phylogenetic relationships within Palaeoheterodonta are unclear, especially for subfamily Unioninae and family Hyriidae. Possibly, this is due to the widespread presence of Doubly Uniparental Inheritance (DUI) among Unionidae, which hampered traditional phylogenetic reconstructions. Therefore, we refer to most recent works on palaeoheterodont evolution ([55]–[60]; and reference therein) and, above all, to the recent work of Breton et al. [61] on the DUI-related comparative mitochondrial genomics of freshwater mussels. In any case, the monophyly of the subclass is not challenged in our study, given the high PP value (1.00) and the length of the branch separating Palaeoheterodonta from their sister group.
Palaeoheterodonta is confirmed to be the sister group of all remaining Autobranchia, as found in our previous study [22]. This is not in agreement with other molecular and morphological studies [12], [16], [19], [21], [23], which considered Palaeoheterodonta more related to Heterodonta than to Pteriomorphia, thus constituting a monophyletic group called Heteroconchia. However, other molecular studies retrieved Palaeoheterodonta as basal to (Heterodonta + Pteriomorphia): Canapa et al. [25] obtained this result on the basis of the RN18S1 nuclear gene, whereas Giribet and Distel [20] used a big dataset and four molecular markers (RN18S1, RN28S1, MT-CO1, and histone H3). Actually, it is unclear why Giribet preferred in his thorough review [12] the Heteroconchia hypothesis when his most recent experimental work was not supporting it [20]. Moreover, a very recent study exploiting complete mitochondrial genomes obtained Palaeoheterodonta to be basal to remaining Autobranchia [26]. Interestingly, the same relationship has been proposed also on morphological grounds: Cope [2], for instance, showed that parsimonious analysis of shell microstructural types leads to similar conclusions.
The Amarsipobranchia
We here contend the monophyly of Heteroconchia sensu [12] and therefore we propose the node-based name “Amarsipobranchia” for the clade comprising Anomalodesmata, Heterodonta, and Pteriomorphia, as it never got a formal name. This term derives from the Greek “marsipos” (µάρσιπος) for “pouch” and means “gills not inserted into a pouch”, in reference to the relationships between anterior filaments of the inner demibranch and the oral groove. In Nuculoidea, Solemyidae, Unionoidea, and possibly Trigonioidea at least the first few anterior filaments are inserted unfused into a distal oral groove, whereas in other bivalves they are fused or not inserted at all ([9], [10], [13]; and reference therein). Although this is not a universal feature of all extant Anomalodesmata, Heterodonta, and Pteriomorphia (for example, inserted unfused anterior filaments are found also in Mytiloidea and Astartidae), this character has most probably to be considered as a symplesiomorphy of this group (see below and Fig. 6E). Although further work is needed to confirm the validity of this taxon, we feel the usefulness of giving a taxonomic name to a clade that is receiving growing support in recent analyses (see, f. i., [2], [20], [22], [24]–[26]).
Anomalodesmata appears to be basal to Heterodonta and Pteriomorphia. In our previous study [22], we obtained anamalodesmatans to be basal to Pteriomorphia, but not monophyletic. In some other studies, anomalodesmatans were found to be a monophyletic clade among Heterodonta [19], [20], [28]–[31] and their subclass status was questioned ([12]; and reference therein). Given our mitochondrial dataset, we can here suggest anomalodesmatans as a monophyletic subclass of Bivalvia, but it is clear that more taxa have to be sampled to completely unravel this point. This is similar to some results of Giribet and Wheeler [19]. Within the subclass, we could not completely confirm the sister group relationship between Pholadomyina and Cuspidariina. Actually, they are also very distinguishable from a morphological point of view, given the eulamellibranch gills of Pandoroidea and the septibranch condition of Cuspidariina [13].
As Astarte and Cardita have been included within Pteriomorphia (see below), the subclass Heterodonta corresponds here to the Euheterodonta sensu [20]. The basal position is occupied by Lucinoidea, confirming the work of John Taylor and colleagues [31], [36]–[38]. Few conclusions can be drawn from this study on Tellinoidea and Donacoidea sensu [62], as the clusters (Abra + Donax) and (Ensis + Sinonovacula) were not completely resolved in the p14 tree. Generally speaking, we tentatively recommend a superfamily Tellinoidea comprising Psammobiidae, Semelidae, and Donacidae, as proposed by Vokes [63]. Our tree shows three more big clusters of Heterodonta, which could correspond to three orders. An order Cardiida sensu novo would contain Hiatelloidea as sister group of Cardioidea, whose only family here represented is the family Cardiidae. Subfamily Tridacninae is basal to remaining subfamilies (Fragine, Laevicardiinae, Cardiinae, Cerastodermatiinae), confirming recent studies on cardiids evolution ([64]–[66]; and reference therein). We retrieved the monophyletic group that Taylor et al. [31] called Neoheterodontei; we recommend the definition of two sister orders Myida and Veneroida sensu novo, which are represented here as (Myoidea + Dreissenoidea) and (Mactroidea + (Glossoidea + Corbiculoidea + Veneroidea)), respectively. The subfamiliar taxonomy of Veneridae needs to be assessed further, as already suggested by Kappner and Bieler [67] and Taylor et al. [31].
Pteriomorphians and their relationships with Astarte and Cardita
Pteriomorphia is robustly monophyletic in our analysis, as repeatedly demonstrated [17], [35]; in this study, however, we present the unexpected result of the inclusion of Astarte cfr. castanea and Cardita variegata within this subclass as sister species. This cluster is consistent with previous molecular and morphological work [19], [20], [31], [68]. Superfamilies Astartoidea, Carditoidea, as well as Crassatelloidea, have generally been regarded as the most primitive heterodonts [19], [20], [34], [69], but also different positions have been proposed [70], [71]. Specifically, Giribet and Distel [20] also proposed Carditoidea (including Astarte castanea) and Crassatelloidea to be the sister group of Nuculanoidea. This is not confirmed since in our study Nuculana commutata is among basal Pteriomorphia (see also [19], [20]), which is commonly accepted nowadays [12], [21]. We prefer the ordinal name Carditoida sensu [21] to indicate this clade, even if they essentially correspond to the Archiheterodonta sensu [31], because this name could lead to confusion if this topology is confirmed.
Deeper inside the pteriomorphian clade, the basal position of Mytilidae is not new, as shown by Waller [16], Carter et al. [72], Steiner and Hammer [17], Giribet and Wheeler [19], and Matsumoto [35] with morphology and molecules (but see [2], [3]). We also agree with Distel [73], who found some concerns about the monophyly of some subfamilies of Mytilidae, namely Mytilinae. We also note that the well known, even if not universally accepted, classification of Ostreina and Pectinina, as suborders of the order Ostreoida, is no longer sustainable, as already noted by Canapa et al. [25], nor is the order Pterioida sensu [63]. We propose to erect an order Nuculanoida for the superfamily Nuculanoidea (see above) and then to regard to pteriomorph systematics in terms of two big clusters. In the first, Anomioidea are basal to Limida sensu [62] as sister group to Pectinoidea, comprising Spondylidae, Propeamussiidae, and Pectinidae in our tree, although further investigations are deserved here, with special reference to Anomiidae (traditionally classified as Pectinina) and pectinid relationships (see, f.i., [74]). Given our tree, we suggest to consider an order Pectinida sensu novo which would include Anomioidea, Limoidea and Pectinoidea.
In the second cluster, we individuate the group (Arcida + (Pinnina + Pteriina + Ostreoida sensu novo)); this leaves unresolved the relationships within the order Pteriida, and it would exclude the possibility to elevate the suborder Pinnina sensu [62] to the ordinal rank. In such a scenario about pteriomorph evolution, Arcida would occupy a somewhat different position with respect to results of Distel [73] and Steiner and Hammer [17], albeit maintaining its basal condition.
Finally, Striarca lactea has been generally classified as member of the subfamily Striarcinae within family Noetiidae; however, several authors have also appraised both subfamilies Striarcinae and Noetiinae as members of the family Arcidae [75]–[77], which would render Arcidae monophyletic in our tree. Moreover, the genus Asperarca Sacco 1898 has been occasionally considered as a synonym of Barbatia Gray 1840 (see, f.i., [62]; but see also [63], [78]), which would render the genus Barbatia paraphyletic in our tree.
Tracing and optimizing major morphological characters on the evolutionary tree
Given the phylogenetic reconstruction we discussed above, major morphological features of bivalve shell and soft parts should be re-evaluated.
Quite surprisingly, the two most used characters for bivalve taxonomy, i.e. gills and shell hinge, do not follow the evolutionary scenarios commonly accepted so far. Protobranch gills (true ctenidia) should be considered the ancestral state among Bivalvia; this is not surprising since most mollusks do have true ctenidia. The question is more puzzling when the “feeding gill” arose among Autobranchia: commonly the filibranch gill has been considered as ancestral, while the eulamellibranch one as derived. The situation, according to our tree, should be exactly the opposite: eulamellibranch gills appear to be the plesiomorphic (ancestral) state in Autobranchia (see Fig. 6A). This is mainly due to the fact that all palaeoheterodonts and most anomalodesmatans, the two groups that arose first among Autobranchia according to our tree, do have an eulamellibranchiate condition (except some anomalodesmatans, which are derived septibranchs). If we accept this, then the filibranch condition of pteriomorphians seems to have evolved from an eulamellibranchiate one. Moreover, according to our tree, the filibranch condition might be occurred at least five times among Pteriomoprhia (Anomioidea, Pectinoidea, Pterioidea, Arcoidea, and Mytiloidea), but there are three unresolved polytomies in this portion of the tree and a better resolution could result in a more parsimonious reconstruction of filibranch condition. Furthermore, filibranchs are already not considered as a natural group [51].
Even more surprisingly, the eulamellibranch condition seems to have reverted to the ancestral protobranchiate state in the superfamily Nuculanoidea. However, the respiratory apparatus of Nuculidae (and Solemyidae, recall the symbiosis with chemoautotrophic bacteria) seems more adherent to the ancestral protobranchiate architecture. Conversely, Nuculanidae possess pumping gills that are very different in function: these are specialized filaments that work as a single septum, ctenidial filaments are alternate and not opposite, the blood space within them is greatly increased, and there is no connection between the palps and the ctenidia [9], [49]. Actually, Yonge [9] stated that in this group “the ctenidia have evolved along lines of their own”. Moreover, siphons are lacking in nuculids and solemyids, whereas they are present in nuculanids, who also possess a posterior unpaired tentacle, a marginal sense organ and three wide digestive diverticula instead of two [9].
Therefore, it is logic to conclude from these anatomical data and from our molecular phylogenetic reconstruction that this peculiar type of gills should not be considered as a variation on the protobranchiate grade, but as an autapomorphy of this family within Pteriomorphia. In this scenario, the words “protobranch” and “ctenidium” would be misleading, as the gills of Nuculanidae are quite unrelated to those of Opponobranchia. We suggest to abandon this terminology and prefer the word “antliobranch” to define this characteristic type of gills, with reference to the pumping action of the filaments (from the Greek áντλεω meaning “ to pump”, “to bail out”). Of course, more studies are needed to better reconsider gills morphology in the light of molecular phylogenies; nevertheless, it has to be noted that what we commonly call protobranch, filibranch or eulamellibranch gills might be artifactual assemblies of different gills types, and maybe these unexpected results will trigger further morphological studies on gills anatomy.
Similarly to gills, the heterodont hinge (once considered more derived) seems to be again the basal condition of Autobranchia (Fig. 6B), so that Nuculanoidea and Arcoidea independently evolved their own taxodont hinge: therefore, taxodont hinges of Nucula, Nuculana, and arks should not be considered as homologous characters. Interestingly, Fang [7] described as basal to all bivalves a “pretaxodont dentition” consisting of 1 or few umbonal teeth, orthomorphodont, that seems very similar to the ancient condition suggested in Figure 6B. Teeth were lost in four cases: Solemyoidea, Dreissenoidea, Hiatelloidea, and most Pteriomorphia, with the exception of Astartoidea and Carditoidea, which retained the ancestral condition of Autobranchia (heterodont hinge), and Nuculanoidea and Arcoidea, which evolved a taxodont dentition on their own. This, as above, needs further studies, once again because different kind of hinges of different origin might possibly hide under the terms heterodont, taxodont and edentate.
On the other hand, the other characters we investigated (gill cilia, stomach type, labial palps and shell microstructure) fit well in the proposed phylogeny. F.i., Type 1 gill cilia are the plesiomorphic condition among bivalves, while Type 2 arose only once in a pteriomorphian clade, excluding Carditoidea + Astartoidea and Mytiloidea, which are therefore supported as basal among pteriomorphians (Fig. 6C). Stomach type (Fig. 6D) again follow quite well the obtained tree and only Type 3 stomach seems to appear twice independently. Labial palps of Type 1 are shared between Opponobranchia and Palaeoheterodonta, thus supporting the basal condition of the latter. Labial palps of Type 3 sensu [10] are synapomorphic for Amarsipobranchia (Fig. 6E), and they mutated into Type 2 in three lineages: Cardioidea, Carditoidea, and Veneroida. Finally, nacreous shell microstructure (Fig. 6F) seems to be the ancestral state of all Bivalvia, while cross lamellar shells appeared once at the arose of Amarsipobranchia.
Finally, as already mentioned, categories we used to map key features on the p14 phylogenetic tree must be taken as just broad umbrella-terms, since most of these character states may hide different discrete conditions. Given the unquestionable interest of a possible re-interpretation of bivalve evolutionary morphology, we hope to trigger further evolutionary work on these issues. More outgroups have also to be included in order to infer the correct characters polarization, meaning more work on mollusks has to be done.
Conclusions
The phylogenetic hypothesis on bivalve evolution we extensively described in the previous paragraph is shown in Figure 7. Its major outcomes and new proposals are: i) mitochondrial genomes are informative for bivalve phylogeny, given a proper phylogenetic approach; ii) the congruence to the established taxonomy of clades obtained from our phylogeny is a further evidence that our tree inference is rather driven by historical signal than homoplasy; iii) the basal subdivision in Opponobranchia and Autobranchia is confirmed; iv) Palaeoheterodonta was retrieved as sister group of a cluster comprising all remaining Autobranchia, which we propose to term Amarsipobranchia; v) Anomalodesmata is apparently monophyletic and maintains a basal status among Amarsipobranchia; vi) three ordinal categories are proposed, namely Cardiida (Hiatelloidea and Cardioidea), Carditoida (Astartoidea and Carditoidea), and Pectinida (Anomioidea, Limoidea, and Pectinoidea); finally, vii) the heterodont hinge and eulamellibranch gills may possibly be re-interpreted as ancestral character states in Autobranchia, and a revision of gill and hinge structures and evolution should be undertaken.
The evolutionary tree was sketched as outlined in this paper (see text for details). Superfamilial relationships are shown, with proposed ordinal classification; for anomalodesmatans, we used the nomenclature from Newell [13] and Vokes [63]. Color code as in Figure 5. Asterisks mark newly-proposed ordinal categories; Neoheterodontei sensu [31] and Amarsipobranchia (comprising subclasses Anomalodesmata, Heterodonta, and Pteriomorphia) are also shown.
In our study, morphological characters and molecular phylogenies are generally in agreement, but sometimes do not. This is not surprising, being different kind of data under different kind of evolutionary histories. Nevertheless, an effort should be taken to better fit both kind of data in Bivalvia, and more integrated work is needed. Finally, especially for deep nodes, the outcomes of molecular phylogenetics should always be compared with, and eventually validated by, all the expertise in the field, merging to a widely accepted phylogenetic hypothesis, encompassing the whole evidence from DNA to morphology.
Further improvements of the present work will increase the available dataset either by exploiting more mitochondrial and nuclear markers or by further enlarging the sample, with special reference to some underrepresented groups (like Anomalodesmata, Anomiidae, Nuculanoidea, Solenoidea, Tellinoidea, Trigonioida): the investigation of deep bivalve phylogeny is as just as started.
Materials and Methods
Taxon Sampling, PCR Amplification, and Sequencing
Species whose sequences were added to the bivalve mitochondrial dataset are listed in Table S4, along with the specimen voucher number of the MoZoo Lab collection (www.mozoolab.net). Specific permits were not needed for the described field studies, as specimens were always collected where such specific permissions were not required: sampling localities were not privately-owned or protected in any way. Moreover, field studies did not involve endangered or protected species. DNA was extracted by Qiagen (Valencia, USA) DNeasy® Blood and Tissue kit. PCR amplification and cloning were carried out as described in [22]: briefly, the Invitrogen (Carlsbad, USA) or ProMega (Madison, USA) Taq DNA polymerase kits were used following manufacturers' instructions to amplify target sequences (the mitochondrial genes MT-RNR1, MT-RNR2, MT-CO1, MT-CYB); a wide range of reaction conditions were used, as different species and markers needed different PCR settings. Typically, a denaturation step of 2′ at 94°C was followed by 35 cycles composed by denaturation of 1′ at 94°C, annealing of 30″−1′ at 46–56°C, and extension of 1′ at 72°C. A final extension step of 5′ at 72°C was added before cooling amplicons at 4°C. We used the same primers as in [22]; specific PCR conditions are available from FP upon request. Sequencing reactions were performed through Macrogen (World Meridian Center, Seoul, South Korea) facility. We put special care into avoiding paralogous sequences due to the presence of the DUI mechanism in some bivalve mollusks: as extensively described in [22], we extracted DNA, whenever possible, from pedal muscle, and avoided to clone amplicons if not strictly requested.
Assembling the dataset
Electropherograms were read through MEGA 4 [79]: sequencer files were manually checked and edited when necessary. The CLC Sequence Viewer 6.4 software (CLC bio, Aarhus, Denmark) was used to organize and to download sequences from GenBank (at December 2010). We then retrieved those taxa for which at least three on four markers were present. Four alignments were prepared with CLC Sequence Viewer and aligned with ClustalW [80] at the EBI server (http://www.ebi.ac.uk/Tools/msa/clustalw2/) [81]. For ribosomal genes, the IUB matrix was used with a 25 penalty for gap opening and a 5 penalty for gap extension; for PCGs, penalties were set to 50 and 10, respectively. When a sequence was not available for a given species, it was replaced with a stretch of missing data in that alignment; Hartmann and Vision ([82]; and reference therein) showed that a large amount of missing data does not lead to incorrect phylogeny in itself, as long as sufficient data are available. 37,579 out of 286,893 (∼13.10%) missing data were inserted for this reason; in the whole matrix, they sum up to 52,152 (∼18.18%), and to 67,232 out of 380,238 (∼17.68%) taking into account also indel presence/absence data (see below). In many cases, we lumped together sequences of different congeneric species to represent the genus they belong to: this is a rather common practice in deep phylogenetic studies and does not lead to inconsistent results at the class level, which is targeted in this study (see, f.i., [22], [83]). Five outgroups were selected for this study: the polyplacophoran Katharina tunicata, two scaphopods (Graptacme eborea and Siphonodentalium lobatum) and two gastropods (Haliotis tuberculata and Thais clavigera). Table S5 lists all sequences used for this study, both downloaded from GenBank and produced in our laboratory.
Region of ambiguous alignment for ribosomal genes were detected by GBlocks [84], [85] with the following parameters: minimum number of sequences for a conserved position, half + 1; minimum number of sequences for a flanking position, half + 1; maximum number of contiguous nonconserved positions, 50; minimum length of a block, 10; allowed gap positions, all. Finally, gaps were coded following the simple indel method of Simmons and Ochoterena [86] as described in [22]; this was carried out with the software GapCoder [87].
Evaluating phylogenetic signal
Taxon sampling was investigated through the method described in [46], which has the property of involving only preexistent taxonomic knowledge about the target group, and does not need any preliminary genetic analysis: for this reason, this is a truly a priori test on taxonomic coverage. All analyses were carried out through the software PhyRe [46] and the bivalve checklist compiled by Millard [62], with 100 random resamplings in all cases. Shuffling test was performed at the family level: 100 master lists were generated and the number of splits, merges, and moves was set to 12, 8, and 4, respectively. We empirically showed in our previous paper [22] that a sample size of about 30 species is sufficient to correctly estimate all molecular evolutionary parameters from a bivalve dataset (given the four mitochondrial markers we employ here); therefore, we did not use any a posteriori test for taxon sampling, as the sample size is more than four times in this study.
A saturation analysis was conducted as in [47], [88] through the program PAUP* 4.0b10 [48] using PaupUp graphical interface [89]. The uncorrected (“p-”) distances were plotted on Maximum Likelihood distances given the proper molecular evolution model (see below). Linear interpolation and its significance were computed with the software PaSt 1.89 [90]. The saturation test was conducted independently for the four markers and, about PCGs, for third codon positions only.
We used SplitsTree 4.6 [91], [92] to obtain phylogenetic networks in which more splits leading to specific clades are shown than in a strictly bifurcating tree. This method aimed to evaluate phylogenetic signal in raw data, therefore the neighbornet network was chosen [93], [94], based on either uncorrected (“p-”) or Log-Det distances.
Presence and properties of phylogenetic signal were also tested with the LM approach [95], [96] as implemented in the software TreePuzzle 5.2 [97], [98]. The complete alignment was used as a dataset, while outgroups were excluded. Molecular evolutionary parameters were given as computed by ModelTest [99] and 1000 random quartets were drawn to produce the final result. Four-cluster Likelihood-Mapping [99] analysis was also conducted: in this case, we excluded outgroups and Opponobranchia (given the stable basal position in all analyses) and subdivided all remaining taxa between four subclasses (Anomalodesmata, Heterodonta, Palaeoheterodonta, and Pteriomoprhia). Significance of results was tested with a Chi-Square test assuming as a null distribution an even presence of observations in each of the three regions of the triangle.
Model decision tests and tree inference
Our dataset was arranged, according to [22], in 26 different partitions: the complete alignment (all), the concatenated ribosomal genes (rib), the concatenated PCGs (prot), individual genes (12s, 16s, cox1, cytb), individual codon positions among the prot partition and single PCGs (prot_1, prot_2, prot_3, cox1_1, cox1_2, cox1_3, cytb_1, cytb_2, cytb_3), the concatenated first and second codon positions (prot_12, cox1_12, cytb_12), and the corresponding indel characters coded as 0/1, irrespective of codon positions (all_indel, rib_indel, prot_indel, 12s_indel, 16s_indel, cox1_indel, cytb_indel). These partitions were assembled in 13 different partitioning schemes, as shown in Table S6. The best-fitting evolutionary model was selected with ModelTest 3.7 using the graphical interface provided by MrMTgui [100]; the Bayesian Information Criterion (BIC) was preferred as a model decision criterion ([101]; and reference therein).
ML analysis was carried out with PAUP* 4.0b10. The alignment was not partitioned and molecular evolutionary parameters computed by ModelTest 3.7 were used for likelihood calculations. Gaps were treated as missing data and binary characters were excluded from the analysis. The outgroups were forced to be paraphyletic with respect to the ingroup. Bootstrap consensus tree using full heuristic ML searches with stepwise additions and TBR branch swapping was constructed to assess nodal support. As described in [22], 150 input files were sent to the University of Oslo Bioportal facility (http://www.bioportal.uio.no) in a parallel run, each computing the maximum likelihood tree for a single bootstrap replicate. Random seed were generated according to PAUP* recommendations with Microsoft Excel® and the consensus tree was computed with Phyutility [102].
All the 13 partitioning schemes were investigated in a Bayesian Analysis with the software MrBayes 3.1.2 [103], [104] hosted at the University of Oslo Bioportal. Initially, the so-called “4by4” nucleotide model (i.e., a traditional 4×4 substitution matrix) was used for all partitioning schemes. For 4 partitioning schemes (see Table 2), namely those containing PCG (prot, cox1, or cytb) partitions, we implemented for PCGs a codon model [105], [106], the M3 model.
10,000,000 generations of two parallel MC3 analyses of 4 chains each were run for each 4by4 partitioning scheme. Since in this analysis we are focusing on the relationships among subclasses, Bivalves were constrained to be monophyletic with respect to the five molluscan outgroups. Nucleotide partitions were treated according to ModelTest results; binary partitions were treated with the default model for restriction data enforcing the coding = variable option and a gamma heterogeneity in substitution distribution. Convergence was estimated by PSRF [107] and by plotting standard deviation of average split frequencies sampled every 1,000 generations. For each M3 analysis 4 independent run of 5,000,000 generations of one single MC3 algorithm were run and convergence among and within runs was estimated via the AWTY tools (http://king2.scs.fsu.edu/CEBProjects/awty/awty_start.php) [108]. A tree was sampled every 100 (4by4 models) or every 125 (M3 models) generations and the consensus was computed at convergence after burnin removal.
The EML computed by MrBayes 3.1.2 was used for AIC [109] and BF [110], as described in ([22]; and reference therein). Briefly, the AIC provides an estimate of the Kullback-Leibler distance [111], i.e. the distance of the model from the reality, considering a penalty computed on the number of free parameters; therefore, smaller values are preferable. On the other hand, the BF involves pairwise comparisons among models through the EML ratio: the larger is the BF value, the more the first model overcomes the second one.
All trees were graphically edited by PhyloWidget [112] and Dendroscope [113] softwares. Optimization of morphological characters on the best evolutionary topology was carried out with Mesquite 2.74 [114]: matrix was taken from [13], with the exception of hinge type, which was coded following [19]. The parsimony method was chosen, as in two cases multiple state characters were coded; in other cases, we tested parsimony results with ML approach, using the MK1 model as implemented by Mesquite.
Supporting Information
Figure S1.
Maximum Likelihood tree. Shown is the consensus tree of 150 bootstrap replicates, using the concatenated alignment as a single partition. Values at the nodes are Bootstrap Proportions (BP); nodes were collapsed if BP<60.
https://doi.org/10.1371/journal.pone.0027147.s001
(EPS)
Table S2.
Molecular evolution models selected by ModelTest 3.7.
https://doi.org/10.1371/journal.pone.0027147.s003
(RTF)
Table S4.
Species used in our laboratory for this study.
https://doi.org/10.1371/journal.pone.0027147.s005
(RTF)
Table S5.
GenBank accession number of sequences used for this study.
https://doi.org/10.1371/journal.pone.0027147.s006
(RTF)
Table S6.
Partitioning schemes adopted for this study.
https://doi.org/10.1371/journal.pone.0027147.s007
(RTF)
Dataset S1.
Complete alignment in NEXUS format. Raw data from 127 taxa (122 bivalves+5 outgroups) are provided for a total of 2,789 molecular characters, without indel coding and removal of ambiguously aligned position.
https://doi.org/10.1371/journal.pone.0027147.s008
(TXT)
Acknowledgments
We would like to thank the captain and the crew of the r/v Urania, Lorenzo Angeletti and volunteers that were essential for efficient sampling during the CNR ARCO and ARCADIA cruises (E.U. Hermes and Hermione projects involving AC and MT) that provided a substantial number of bivalve samples used in the present paper. We are also very thankful to Mirco Bergonzoni, Gianluigi Bini, Valter Gasperi, Diego Viola, and Gianluca Zuffi for providing specimens used for this study. Special thanks are due to Valeria Abram and Lorenzo Donati for their invaluable help during lab work. This is ISMAR, Bologna scientific contribution n. 1715. The present work was greatly improved by comments and suggestions of three anonymous referees.
Author Contributions
Conceived and designed the experiments: FP MP. Performed the experiments: FP. Analyzed the data: FP. Contributed reagents/materials/analysis tools: FP MP. Wrote the paper: FP AC MT MP. Organized the two Adriatic oceanographic cruises ARCO and ARCADIA, which made specimens' collection possible: AC MT.
References
- 1. Runnegar B, Pojeta J Jr (1992) The earliest bivalves and their Ordovician descendants. Am Malacol Bull 9: 117–122.
- 2.
Cope JCW (1996) The early evolution of the Bivalvia. In: Taylor JD, editor. Origin and Evolutionary Radiation of the Mollusca. Oxford: Oxford University Press. pp. 361–370.
- 3.
Morton B (1996) The evolutionary history of the Bivalvia. In: Taylor JD, editor. Origin and Evolutionary Radiation of the Mollusca. Oxford: Oxford University Press. pp. 337–359.
- 4. Cope JCW, Babin C (1999) Diversification of bivalves in the Ordovician. Geobios 32: 175–185.
- 5. Pojeta J Jr (2000) Cambrian Pelecypoda (Mollusca). Am Malacol Bull 15: 157–166.
- 6.
Cope JCW (2002) Diversification and biogeography of bivalves during the Ordovician Period. In: Crame JA, Owen AW, editors. Palaeobiogeography and Biodiversity Change: the Ordovician and Mesozoic-Cenozoic Radiations. London: Geological Society. pp. 25–52.
- 7. Fang Z-J (2006) An introduction to Ordovician bivalves of southern China, with a discussion of the early evolution of the Bivalvia. Geol J 41: 303–328.
- 8. Tsubaki R, Kameda Y, Kato M (2011) Pattern and process of diversification in an ecologically diverse epifaunal bivalve group Pterioidea (Pteriomorphia, Bivalvia). Mol Phylogenet Evol 58: 97–104.
- 9. Yonge CM (1939) The protobranchiate Mollusca: a functional interpretation of their structure and evolution. Phil Trans R Soc B 230: 79–147.
- 10. Stasek CR (1963) Synopsis and discussion of the association of ctenidia and labial palps in the bivalve Mollusca. Veliger 6: 91–97.
- 11. Stanley SM (1968) Post-Paleozoic adaptive radiation of infaunal bivalve mollusks: a consequence of mantle fusion and siphon formation. J Paleontol 42: 214–229.
- 12.
Giribet G (2008) Bivalvia. In: Ponder WF, Lindberg DR, editors. Phylogeny and Evolution of the Mollusca. Berkeley: University of California Press. pp. 105–142.
- 13. Newell ND (1965) Classification of the Bivalvia. Am Mus Novit 2206: 1–25.
- 14. Starobogatov YI (1992) Morphological basis for phylogeny and classification of Bivalvia. Ruthenica 2: 1–25.
- 15.
von Salvini-Plawen L, Steiner G (1996) Synapomorphies and plesiomorphies in higher classification of Mollusca. In: Taylor JD, editor. Origin and Evolutionary Radiation of the Mollusca. Oxford: Oxford University Press. pp. 29–51.
- 16.
Waller TR (1998) Origin of the molluscan class Bivalvia and a phylogeny of major groups. In: Johnston PA, Haggart JW, editors. Bivalves: An Eon of Evolution. Calgary: University of Calgary Press. pp. 1–45.
- 17.
Steiner G, Hammer S (2000) Molecular phylogeny of the Bivalvia inferred from 18S rDNA sequences with particular reference to the Pteriomorphia. In: Harper EM, Taylor JD, Crame JA, editors. The Evolutionary Biology of the Bivalvia. London: The Geological Society of London. pp. 11–29.
- 18. Passamaneck YJ, Schander C, Halanych KM (2004) Investigation of molluscan phylogeny using large-subunit and small-subunit nuclear rRNA sequences. Mol Phylogenet Evol 32: 25–38.
- 19. Giribet G, Wheeler W (2002) On bivalve phylogeny: A high-level analysis of the Bivalvia (Mollusca) based on combined morphology and DNA sequence data. Invert Biol 121: 271–324.
- 20.
Giribet G, Distel DL (2003) Bivalve phylogeny and molecular data. In: Lydeard C, Lindberg DR, editors. Molecular systematics and phylogeography of mollusks. Washington: Smithsonian Books. pp. 45–90.
- 21. Bieler R, Mikkelsen PM (2006) Bivalvia – a look at the branches. Zool J Linn Soc 148: 223–235.
- 22. Plazzi F, Passamonti M (2010) Towards a molecular phylogeny of Mollusks: Bivalves' early evolution as revealed by mitochondrial genes. Mol Phylogenet Evol 57: 641–657.
- 23.
Waller TR (1990) The evolution of ligament systems in the Bivalvia. In: Morton B, editor. The Bivalvia. Hong Kong: Hong Kong University Press. pp. 49–71.
- 24. Cope JCW (1997) The early phylogeny of the class Bivalvia. Palaeontology 40: 713–746.
- 25. Canapa A, Marota I, Rollo F, Olmo E (1999) The small-subunit rRNA gene sequences of venerids and the phylogeny of Bivalvia. J Mol Evol 48: 463–468.
- 26. Doucet-Beaupré H, Breton S, Chapman EG, Blier PU, Bogan AE, et al. (2010) Mitochondrial phylogenomics of the Bivalvia (Mollusca): searching for the origin and mitogenomic correlates of doubly uniparental inheritance of mtDNA. BMC Evol Biol 10: 50.
- 27.
Myra Keen A (1963) Marine molluscan genera of Western North America: an illustrated key. Stanford: Stanford University Press. 126 p.
- 28.
Harper EM, Hide EA, Morton B (2000) Relationships between the extant Anomalodesmata: a cladistic test. In: Harper EM, Taylor JD, Crame JA, editors. The Evolutionary Biology of the Bivalvia. London: The Geological Society of London. pp. 129–143.
- 29. Harper EM, Dreyer H, Steiner G (2006) Reconstructing the Anomalodesmata (Mollusca: Bivalvia): morphology and molecules. Zool J Linn Soc 148: 395–420.
- 30. Dreyer H, Steiner G, Harper EM (2003) Molecular phylogeny of Anomalodesmata (Mollusca: Bivalvia) inferred from 18S rRNA sequences. Zool J Linn Soc 139: 229–246.
- 31. Taylor JD, Williams ST, Glover EA, Dyal P (2007) A molecular phylogeny of heterodont bivalves (Mollusca: Bivalvia: Heterodonta): new analyses of 18S and 28S rRNA genes. Zool Scr 36: 587–606.
- 32. Steiner G, Müller M (1996) What can 18S rDNA do for bivalve phylogeny? J Mol Evol 43: 58–70.
- 33. Adamkewicz SL, Harasewych MG, Blake J, Saudek D, Bult CJ (1997) A molecular phylogeny of the bivalve mollusks. Mol Biol Evol 14: 619–629.
- 34.
Campbell DC (2000) Molecular evidence on the evolution of the Bivalvia. In: Harper EM, Taylor JD, Crame JA, editors. The Evolutionary Biology of the Bivalvia. London: The Geological Society of London. pp. 31–46.
- 35. Matsumoto M (2003) Phylogenetic analysis of the subclass Pteriomorphia (Bivalvia) from mtDNA COI sequences. Mol Phylogenet Evol 27: 429–440.
- 36. Williams ST, Taylor JD, Glover EA (2004) Molecular phylogeny of the Lucinoidea (Bivalvia): non-monophyly and separate acquisition of bacterial chemosymbiosis. J Moll Stud 70: 187–202.
- 37. Taylor JD, Williams ST, Glover EA (2007) Evolutionary relationships of the bivalve family Thyasiridae (Mollusca: Bivalvia), monophyly and superfamily status. J Mar Biol Ass U K 87: 565–574.
- 38. Taylor JD, Glover EA, Williams ST (2009) Phylogenetic position of the bivalve family Cyrenoididae – removal from (and further dismantling of) the superfamily Lucinoidea. Nautilus 123: 9–13.
- 39. Wheeler W (1996) Optimization alignment: The end of multiple sequence alignment in phylogenetics? Cladistics 12: 1–9.
- 40. Sorenson MD, Quinn TW (1998) Numts: A Challenge for Avian Systematics and Population Biology. The Auk 115: 214–221.
- 41. Bensasson D, Zhang D-X, Hartl DL, Hewitt GM (2001) Mitochondrial pseudogenes: Evolution's misplaced witnesses. TRENDS Ecol. Evol. 16: 314–321.
- 42. Zbawicka M, Burzyński A, Wenne R (2007) Complete sequences of mitochondrial genomes from the Baltic mussel Mytilus trossulus. Gene 406: 191–198.
- 43. Mindell DP, Sorenson MD, Dimcheff DE (1998) An Extra Nucleotide Is Not Translated in Mitochondrial ND3 of Some Birds and Turtles. Mol Biol Evol 15: 1568–1571.
- 44. Beckenbach AT, Robson SKA, Crozier R (2005) Single Nucleotide +1 Frameshifts in an Apparently Functional Mitochondrial Cytochrome b Gene in Ants of the Genus Polyrhachis. J Mol Evol 60: 141–152.
- 45. Rosengarten RD, Sperling EA, Moreno MA, Leys SP, Dellaporta SL (2008) The mitochondrial genome of the hexactinellid sponge Aphrocallistes vastus: Evidence for programmed translational frameshifting. BMC Genomics 9: 33.
- 46. Plazzi F, Ferrucci RR, Passamonti M (2010) Phylogenetic Representativeness: a new method for evaluating taxon sampling in evolutionary studies. BMC Bioinformatics 11: 209.
- 47. Tanabe Y, Saikawa M, Watanabe MM, Sugiyama J (2004) Molecular phylogeny of Zygomycota based on EF-1α and RPB1 sequences: limitations and utility of alternative markers to rDNA. Mol Phylogenet Evol 30: 438–449.
- 48.
Swofford DL (2002) PAUP*. Sunderland: Sinauer Associates. Phylogenetic analysis using parsimony (*and other methods) version 4.0b10.
- 49. Atkins D (1936) On the ciliary mechanisms and interrelationships of lamellibranchs. Part I. Some new observations on sorting mechanisms in certain Lamellibranchia. Quarterly Journal of Microscopical Science 79: 181–308.
- 50. Atkins D (1937) On the ciliary mechanisms and interrelationships of lamellibranchs. Part III. Types of lamellibranch gills and their food currents. Quarterly Journal of Microscopical Science 79: 375–421.
- 51. Atkins D (1938) On the ciliary mechanisms and interrelationships of lamellibranchs. Part VII. Latero-frontal cilia of the gill filaments and their phylogenetic value. Quarterly Journal of Microscopical Science 80: 345–436.
- 52. Purchon RD (1958) Phylogeny in the Lamellibranchia. Proc Cent And Bicent Congr Biol Singapore 69-82:
- 53. Brasier MD, Hewitt RA (1978) On the Late Precambrian- Early cambrian Hartshill Formation of Warwickshire. Geol Mag 115: 21–36.
- 54.
Morton B, Yonge CM (1964) Classification and structure of the Mollusca. In: Wilbur KM, Yonge CM, editors. Physiology of Mollusca, 1. New York: Academic Press. pp. 1–58.
- 55. Graf DL, Ó Foighil D (2000) The evolution of brooding characters among the freshwater pearly mussels (Bivalvia: Unionoidea) of North America. J Molluscan Stud 66: 157–170.
- 56. Machordom A, Araujo R, Erpenbeck D, Ramos M-Á (2003) Phylogeography and conservation genetics of endangered European Margaritiferidae (Bivalvia: Unionoidea). Biol J Linn Soc 78: 235–252.
- 57.
Roe KJ, Hoeh WR (2003) Systematics of freshwater mussels (Bivalvia: Unionoida). In: Lydeard C, Lindberg DR, editors. Molecular Systematics and Phylogeography of Mollusks. Washington: Smithsonian Books.. pp. 91–122.
- 58. Serb JM, Buhay JE, Lydeard C (2003) Molecular systematic of the North American freshwater bivalve genus Quadrula (Unionidae: Ambleminae) based on mitochondrial ND1 sequences. Mol Phylogenet Evol 28: 1–11.
- 59. Huff SW, Campbell D, Gustafson DL, Lydeard C, Altaba CR, et al. (2004) : Investigations into the phylogenetic relationships of the threatened freshwater pearl-mussels (Bivalvia, Unionoidea, Margaritiferidae) based on molecular data: implications for their taxonomy and biogeography. J Molluscan Stud 70: 379–388.
- 60. Khalloufi N, Toledo C, Machordom A, Boumaïza M, Araujo R (2011) The unionids of Tunisia: taxonomy and phylogenetic relationships, with redescription of Unio ravoisieri Deshayes, 1847 and U. durieui Deshayes, 1847. J Molluscan Stud 0: 1–13.
- 61. Breton S, Doucet-Beaupré H, Stewart DT, Piontkivska H, Karmakur M, et al. (2009) Comparative Mitochondrial Genomics of Freshwater Mussels (Bivalvia: Unionoida) With Doubly Uniparental Inheritance of mtDNA: Gender-Specific Open Reading Frames and Putative Origins of Replication. Genetics 183: 1575–1589.
- 62.
Millard V (2001) Classification of Mollusca: A classification of world wide Mollusca. 2nd edition. South Africa. 1931 p.
- 63.
Vokes HE (1980) Genera of the Bivalvia. Ithaca: Palaeontological Research Institution. 307 p.
- 64. Maruyama T, Ishikura M, Yamazaki S, Kanai S (1998) Molecular phylogeny of zooxanthellate bivalves. Biol Bull 195: 70–77.
- 65. Schneider JA, Ó Foighil D (1999) Phylogeny of giant clams (Cardiidae: Tridacninae) based on partial mitochondrial 16S rDNA gene sequences. Mol Phylogenet Evol 13: 59–66.
- 66. Kirkendale L (2009) Their Day in the Sun: molecular phylogenetics and origin of photosymbiosis in the ‘other’ group of photosymbiotic marine bivalves (Cardiidae: Fraginae). Biol J Linnean Soc 97: 448–465.
- 67. Kappner I, Bieler R (2006) Phylogeny of venus clams (Bivalvia: Venerinae) as inferred from nuclear and mitochondrial gene sequences. Mol Phylogenet Evol 40: 317–331.
- 68. Healy JM (1995) Sperm ultrastructure in the marine bivalve families Carditidae and Crassatellidae and its bearing of unification of the Crassatelloidea with the Carditoidea. Zool Scr 24: 21–28.
- 69. Park JK, Ó Foighil D (2000) Sphaeriid and corbiculid clams represent separate heterodont bivalve radiations into freshwater environments. Mol Phylogenet Evol 14: 75–88.
- 70. Yonge CM (1969) Functional morphology and evolution within the Carditacea (Bivalvia). Proceedings of the Malacological Society of London 38: 493–527.
- 71. Purchon RD (1987) Classification and evolution of the Bivalvia: An analytical study. Phil Trans R Soc B 316: 277–302.
- 72.
Carter JG, Campbell DC, Campbell MR (2000) Cladistic perspectives on early bivalve evolution. In: Harper EM, Taylor JD, Crame JA, editors. The Evolutionary Biology of the Bivalvia. London: The Geological Society of London. pp. 47–79.
- 73. Distel DL (2000) Phylogenetic relationships among Mytilidae (Bivalvia): 18S rRNA data suggest convergence in mytilid body plans. Mol Phylogenet Evol 15: 25–33.
- 74. Puslednik L, Serb JM (2008) Molecular phylogenetics of the Pectinidae (Mollusca: Bivalvia) and effect of increased taxon sampling and outgroup selection on tree topology. Mol Phylogenet Evol 48: 1178–1188.
- 75. Reinhart PW (1935) Classification of the pelecypod family Arcidae. Bulletin du Musée Royal d'Histoire Naturelle de Belgique 11: 5–68.
- 76. Rost H (1955) A report on the family Arcidae. Allan Hancock Foundation publications. Series 1. Allan Hancock Pacific expeditions 20: 177–249.
- 77.
Myra Keen A (1971) Sea shells of Tropical West America: Marine mollusks from Baja California to Peru. Stanford: Stanford University Press. 1080 p. 2 ed.
- 78. La Perna R (1998) On Asperarca Sacco, 1898 (Bivalvia, Arcidae) and two new Mediterranean species. Bollettino Malacologico 33: 11–18.
- 79. Tamura K, Dudley J, Nei M, Kumar S (2007) MEGA4: Molecular Evolutionary Genetics Analysis (MEGA) software version 4.0. Mol Biol Evol 24: 1596–1599.
- 80. Thompson JD, Higgins DG, Gibson TJ (1994) CLUSTAL W: improving the sensitivity of progressive multiple sequence alignment through sequence weighting, position-specific gap penalties and weight matrix choice. Nucleic Acids Res 22: 4673–4680.
- 81. Chenna R, Sugawara H, Koike T, Lopez R, Gibson TJ, et al. (2003) Multiple sequence alignment with the Clustal series of programs. Nucleic Acids Res 31: 3497–3500.
- 82. Hartmann S, Vision TJ (2008) Using ESTs for phylogenomics: Can one accurately infer a phylogenetic tree from a gappy alignment? BMC Evol Biol 8: 95.
- 83. Li B, Dettaï A, Cruaud C, Couloux A, Desoutter-Meniger M, et al. (2009) RNF213, a new nuclear marker for acanthomorph phylogeny. Mol Phylogenet Evol 50: 345–363.
- 84. Castresana J (2000) Selection of Conserved Blocks from Multiple Alignments for Their Use in Phylogenetic Analysis. Mol Biol Evol 17: 540–552.
- 85. Talavera G, Castresana J (2007) Improvement of Phylogenies after Removing Divergent and Ambiguously Aligned Blocks from Protein Sequence Alignments. Syst Biol 56: 564–577.
- 86. Simmons MP, Ochoterena H (2000) Gaps as characters in sequence-based phylogenetic analyses. Syst Biol 49: 369–381.
- 87. Young ND, Healy J (2003) GapCoder automates the use of indel characters in phylogenetic analysis. BMC Bioinformatics 4: 6.
- 88. Pereyra S, García G (2008) Patterns of genetic differentiation in the Gymnogeophagus gymnogenys species complex, a neotropical cichlid from South American basins. Environ Biol Fish 83: 245–257.
- 89.
Calendini F, Martin J-F (2009) PaupUP Home Page. Available: http://www.agro-montpellier.fr/sppe/Recherche/JFM/PaupUp/main.htm. Accessed 2009 Jan 21.
- 90. Hammer Ø, Harper DAT, Ryan PD (2001) PAST: paleontological statistics software package for education and data analysis. Palaeontol Electron 4: 1–9.
- 91. Dress A, Huson D, Moulton V (1996) Analyzing and Visualizing Sequence and Distance Data Using SplitsTree. Discrete Appl Math 71: 95–109.
- 92. Huson D, Bryant D (2006) Application of phylogenetic networks in evolutionary studies. Mol Biol Evol 23: 254–267.
- 93. Bryant D, Moulton V (2004) Neighbor-Net: An Agglomerative Method for the Construction of Phylogenetic Networks. Mol Biol Evol 21: 255–265.
- 94. Wägele JW, Letsch H, Klussmann-Kolb A, Mayer C, Misof B, et al. (2009) Phylogenetic support values are not necessarily informative: the case of the Serialia hypothesis (a mollusk phylogeny). Frontiers in Zoology 6: 12.
- 95. Strimmer K, von Haeseler A (1996) Quartet Puzzling: A Quartet Maximum-Likelihood Method for Reconstructing Tree Topologies. Mol Biol Evol 13: 964–969.
- 96. Strimmer K, von Haeseler A (1997) Likelihood-mapping: A simple method to visualize phylogenetic content of a sequence alignment. Proc Natl Acad Sci U S A 94: 6815–6819.
- 97. Schmidt HA, Strimmer K, Vingron M, von Haeseler A (2002) TREE-PUZZLE: maximum likelihood phylogenetic analysis using quartets and parallel computing. Bioinformatics 18: 502–504.
- 98.
Schmidt HA, von Haeseler A (2003) Maximum-Likelihood Analysis Using TREE-PUZZLE. In: Baxevanis AD, Davison DB, Page RDM, Stormo G, Stein L, editors. Current Protocols in Bioinformatics. New York: Wiley and Sons. pp. 6.6.1–6.6.25.
- 99. Posada D, Crandall KA (1998) Modeltest: Testing the model of DNA substitution. Bioinformatics 14: 817–818.
- 100.
Nuin P (2007) MrMTgui: Cross-platform interface for ModelTest and MrModeltest. Available: http://www.genedrift.org/mtgui.php. Accessed 2007 Jul 11.
- 101. Luo A, Qiao H, Zhang Y, Shi W, Ho SYW, et al. (2010) Performance of criteria for selecting evolutionary models in phylogenetics: a comprehensive study based on simulated datasets. BMC Evol Biol 10: 242.
- 102. Smith SA, Dunn CW (2008) Phyutility: a phyloinformatics tool for trees, alignments, and molecular data. Bioinformatics 24: 715–716.
- 103. Huelsenbeck JP, Ronquist F (2001) MRBAYES: Bayesian inference of phylogeny. Bioinformatics 17: 754–755.
- 104. Ronquist F, Huelsenbeck JP (2003) MRBAYES 3: Bayesian phylogenetic inference using mixed models. Bioinformatics 19: 1572–1574.
- 105. Goldman N, Yang Z (1994) A codon-based model of nucleotide substitution for protein coding DNA sequences. Mol Biol Evol 11: 725–736.
- 106. Muse SV, Gaut BS (1994) A likelihood approach for comparing synonymous and nonsynonymous substitution rates, with application to the chloroplast genome. Mol Biol Evol 11: 715–724.
- 107. Gelman A, Rubin DB (1992) Inference from iterative simulation using multiple sequences. Stat Sci 7: 457–511.
- 108. Nylander JAA, Wilgenbusch JC, Warren DL, Swofford DL (2008) AWTY (are we there yet?): A system for graphical exploration of MCMC convergence in Bayesian phylogenetics. Bioinformatics 24: 581–583.
- 109.
Akaike H (1973) Information theory and an extension of the maximum likelihood principle. In: Petrox BN, Caski F, editors. Second International Symposium on Information Theory. Budapest: Akademiai Kiado.. 267 p.
- 110. Kass RE, Raftery AE (1995) Bayes Factors. J Am Stat Assoc 90: 773–795.
- 111. Kullback S, Leibler RA (1951) On information and sufficiency. Ann Math Stat 22: 79–86.
- 112. Jordan GE, Piel WH (2008) PhyloWidget: Web-based visualizations for the tree of life. Bioinformatics 15: 1641–1642.
- 113. Huson DH, Richter DC, Rausch C, Dezulian T, Franz M, et al. (2007) Dendroscope – An interactive viewer for large phylogenetic trees. BMC Bioinformatics 8: 460.
- 114.
Maddison WP, Maddison DR (2011) Mesquite. Available: http://mesquiteproject.org. Accessed 2011 Mar 2.
- 115. Nieselt-Struwe K, von Haeseler A (2001) Quartet-Mapping, a Generalization of the Likelihood-Mapping Procedure. Mol Biol Evol 18: 1204–1219.