Origin and Diversification of Major Clades in Parmelioid Lichens (Parmeliaceae, Ascomycota) during the Paleogene Inferred by Bayesian Analysis

There is a long-standing debate on the extent of vicariance and long-distance dispersal events to explain the current distribution of organisms, especially in those with small diaspores potentially prone to long-distance dispersal. Age estimates of clades play a crucial role in evaluating the impact of these processes. The aim of this study is to understand the evolutionary history of the largest clade of macrolichens, the parmelioid lichens (Parmeliaceae, Lecanoromycetes, Ascomycota) by dating the origin of the group and its major lineages. They have a worldwide distribution with centers of distribution in the Neo- and Paleotropics, and semi-arid subtropical regions of the Southern Hemisphere. Phylogenetic analyses were performed using DNA sequences of nuLSU and mtSSU rDNA, and the protein-coding RPB1 gene. The three DNA regions had different evolutionary rates: RPB1 gave a rate two to four times higher than nuLSU and mtSSU. Divergence times of the major clades were estimated with partitioned BEAST analyses allowing different rates for each DNA region and using a relaxed clock model. Three calibrations points were used to date the tree: an inferred age at the stem of Lecanoromycetes, and two dated fossils: Parmelia in the parmelioid group, and Alectoria. Palaeoclimatic conditions and the palaeogeological area cladogram were compared to the dated phylogeny of parmelioid. The parmelioid group diversified around the K/T boundary, and the major clades diverged during the Eocene and Oligocene. The radiation of the genera occurred through globally changing climatic condition of the early Oligocene, Miocene and early Pliocene. The estimated divergence times are consistent with long-distance dispersal events being the major factor to explain the biogeographical distribution patterns of Southern Hemisphere parmelioids, especially for Africa-Australia disjunctions, because the sequential break-up of Gondwana started much earlier than the origin of these clades. However, our data cannot reject vicariance to explain South America-Australia disjunctions.


Introduction
In traditional, morphology-based concepts, lichenized fungi often have wide distributions spanning over several continents with a number of species being cosmopolitan. This has led to a widespread notion that the distribution of these fungi is primarily shaped by ecological conditions rather than explained by historical events. In contrast, a few authors have invoked plate tectonics and emphasized vicariance to explain distribution patterns of lichens, especially for species occurring in the Southern Hemisphere [1][2][3][4][5][6]. However, within the last decade, molecular data have helped to revolutionize the species delimitation of lichenized fungi and demonstrated that morphology-based concepts largely underestimate the diversity of lichens [7][8][9][10][11][12][13][14][15][16]. As a result of these studies, it became clear that -although the number of widely distributed species in lichenized fungi is generally higher than in plants or animals -lichens have more restricted distribution areas than previously assumed. For example, supposedly cosmopolitan species in Parmeliaceae and Physciaceae were found to represent distinct lineages in different continents [11,13]. Further, recent progress in our knowledge of the phylogeny of some clades of lichenized fungi revealed the presence of clades at generic rank that originated and diversified in the Southern Hemisphere [17][18][19][20][21]. Hence, we have turned our attention to address the issue of the extent of vicariance and long-distance dispersal to understand the current distribution of lichenized fungi anew using molecular phylogenies.
For this purpose a dated phylogeny with the estimated ages of origin and diversification of the parmelioid group is required. A main problem for building dated phylogenies in fungi is the poor fossil record. While our understanding of the divergence time of angiosperms is well established [22][23][24][25],. time estimates for fungi were long disputed based on uncertainties in the interpretation of the few known fossils [26]. Consequently, published dating estimates ranged from 660 million to 2.15 billion years ago for the origin of Fungi and from 390 million to 1.5 billion years for the split of the two crown groups of fungi, Ascomycota and Basidiomycota [27][28][29][30][31]. Re-examination of the morphology of the fossils and re-evaluation of the published dating studies, however, suggested more consistent results with the origin of the Fungi dating back to between 760 million and 1.06 billion years, and the split of Ascomycota and Basidiomycota for about 500-650 Ma [26,32]. This suggests that terrestrial fungi evolved and diversified more or less simultaneously with the evolution of land plants.
Given that the phylogeny of parmelioid lichens is well resolved [18] and the recent advantages in our understanding of the timing of major events in the evolution of fungi, we feel confident that the times of main nodes differentiation can be estimated. Moreover, this dating approach can be applied to address the issue of the extent of vicariance on distribution patterns in lichens and use parmelioids as an example.
The family Parmeliaceae is widely distributed throughout the world from polar to tropical regions and is one of the largest families of lichenized Ascomycota [18,[33][34][35][36]. Most species form foliose or fruticose thalli, but some also have subcrustose, umbilicate, peltate thalli, and even lichenicolous fungi were found to belong here [18,35,37]. The family includes about 2500 species classified in 84 genera and is characterized by cup-shaped apothecia with cupulate exciple, Lecanora-type asci, often with hyaline ascospores [35]. Within the family, six strongly supported major monophyletic groups can be distinguished [35], which are: alectorioid, cetrarioid, hypogymnioid, letharioid, parmelioid and psiloparmelioid. By far, the largest of these groups is the parmelioid group with about 1500 species [38,39] classified in 27 accepted genera [18]. Phenotypically, the parmelioid lichens are characterized by having mostly foliose thalli with thread-like rhizines on the lower surface, cup-shaped apothecia on the thallus upper surface and Lecanora-type asci with hyaline ascospores (Fig. 1). Within the parmelioid lichens, eight major monophyletic clades can be distinguished, which are Cetrelia-, Hypotrachyna-, Melanohalea-, Parmelia-, Parmelina-, Parmeliopsis-, Parmotremaand Xanthoparmelia-groups [18]. During the last decades, phylogenetic studies based on DNA sequence data have greatly advanced our understanding of the evolution of the family including the phylogenetic relationships among major clades [18,[34][35][40][41][42][43][44][45][46], phenotypical evolution [35,45,47], disparity in substitution rates among clades [48], and the geographic origin of certain clades [17,19]. However, so far none of these have aimed at dating major cladogenesis events in Parmeliaceae. Besides the issues discussed above, this is probably also due to the poor fossil record of Parmeliaceae. Only very few fossils have been recorded of Parmeliaceae and the known ones are all preserved in amber. An Alectoria species was described from Baltic amber (35)(36)(37)(38)(39)(40) [49]. C. Melanelixia glabra occurs from southwestern Europe to eastern Russia. D. Parmelia sulcata, a common species from cold to temperate regions of both Hemispheres. E. Flavoparmelia soredians occurs in warm and humid areas of temperate regions of both Hemispheres. F. Xanthoparmelia conspersa, one of the most widespread species of macrolichens growing on acid rocks in temperate areas of both Hemipheres, excluding Australia and South Africa. The distribution areas after Nimis [105]. All photographs were taken in the field. Scale = 1 cm. doi:10.1371/journal.pone.0028161.g001 Two fossil species of Parmelia were described from Dominican amber   [50] and two specimens of Anzia have been found from European amber (35)(36)(37)(38)(39)(40) [51].
The aims of this study are 1) to put a time-scale on the phylogenetic tree of parmelioid lichens and thus identify when the main nodes differentiated, employing calibration points from a recent dating estimate of fungi [26] and available fossil referencepoints, and 2) to address the impact of vicariance and longdistance dispersal processes in the distribution patterns in parmelioid lichens. Specifically we address the question whether plate tectonics can be invoked to explain distribution areas of species groups or species occurring in the Southern Hemisphere. The substitution rates of the three loci (nuLSU, mtSSU and RPB1) are shown for each lineage of the parmelioid lichens in Fig. 2 and Table 1. The different lineages have similar range of variation in each gene but there are clear differences in mean rates between the three genes, the rate of RPB1 being two to four times higher than mtSSU and nuLSU, respectively.

Estimating divergence times
A chronogram based on the analysis of the combined matrix of three loci is shown in Fig. 3, showing the relationships of Lecanoromycetes and Parmeliaceae (highlighted in dark grey color). The detailed chronogram for Parmeliaceae is depicted in Fig. 4. This analysis used as calibration points the age of the split of Lecanoromycetes-Chaethothyriomycetidae (C1); the age of the fossil (P. ambra) assigned to the Parmelia s.s. crown node (C2), and the age of the fossil A. succini assigned to the Alectoria crown (C3). The node ages, 95% highest posterior density intervals (HPD) of ages, substitution rates and 95% HPD of substitution rates for the main clades, divergence points and diversification point of the genera are shown in Table 2.
The alternative analyses gave similar results ( Table 3). The small differences are as follows:1) using two calibration points (C1 and C2) and excluding the Alectoria fossil, node ages and 95% HPD intervals are slightly younger for some of the parmelioid clades; 2) using C1 and C3 and excluding Parmelia, the node ages and 95% HPD intervals had a small decrease in older clades; and 3) using C1, C2 and the Alectoria fossil constraining the age of the alectorioid clade (Alectoria, Bryoria and Pseudephebe; C3*) the node ages are similar with very small fluctuations. The 95% HPD intervals for all the clades largely agree and are similar to those obtained in the first analysis using the calibration points C1, C2 and C3.
The split of the Parmeliaceae core was estimated to be around 109 Ma (85.52-136.55 Ma) when the crustose genus Protoparmelia separated from the rest of the Parmeliaceae (Fig. 4, Table 2). The basal radiation of Parmeliaceae core took place between 60 and Our analyses suggest seven separate major divergence events that led to the evolution of the main clades of parmelioid lichens (Fig. 4, marked with dots and numbers). The earliest divergence,

Discussion
In this study we put a time-scale on the phylogenetic tree of parmelioid lichens using three DNA regions with different evolutionary rates. The estimated ages allow addressing how vicariance and long-distance dispersal shaped the current distribution patterns of parmelioid lichens, specially the disjunct distributions of species groups in the Southern Hemisphere. Moreover the dated phylogeny provides a general picture of the palaeoclimatic conditions prevalent on Earth when the main lineages differentiated.
The three DNA regions used in this study to build the phylogenetic trees have different evolutionary rates. We found higher substitution rates in the protein coding gene RPB1 than in the nuLSU and mtSSU ribosomal DNA (Fig. 2, Table 2). Ribosomal DNA has been frequently used in molecular studies of Parmeliaceae and other lichenized fungi [34,[41][42][44][45][46][47][52][53][54] but so far few molecular studies have used protein-coding genes to infer phylogenetic relationships in Parmeliaceae [18,35,43,[55][56]. In our analysis the RPB1 gene (with high substitution rates) provided better resolution of the terminal lineages of the tree while the more conserved genes with lower substitution rate (nuLSU, mtSSU) better supported the backbone of the tree topology.
Substitution rates have been used to estimate divergence times due to the lack of fossil records. Takamatsu & Matsuda [57] calculated a substitution rates for Erysiphales (2.52610 29 s/s/y for nuLSU). Our estimation of the substitution rate for the nuLSU of the Parmeliaceae (1.12610 29 s/s/y; Table 2) is in the same order of magnitude but slower than in the phytopathogenous fungi of Erysiphales. This indicates that caution should be used when applying substitution rates of a group of fungi to estimate divergence times of an unrelated group [58][59].
The scarcity and uncertainty of the fossil record was a major obstacle for estimation of dates of radiation in most groups of fungi. We have estimated divergence times of the major lineages of Parmeliaceae and times of the main radiations of the genera using a comprehensive phylogenetic hypothesis of the family, calibrated with the available fossil evidence, a root time inferred from [26], and allowing a relaxed-clock model for the rates of evolution of the main clades.
Berbee & Taylor [32] estimated the age of the parmeliaceae crown node at about 60 Ma based on substitution rates, and a minimum age of the family was given at about 40 Ma according to fossil records [50][51]. However, the age estimates obtained herein suggested that Parmeliaceae evolved much earlier. Our analyses indicate that the Parmeliaceae core originated rather recently with a stem node age estimate of 108 Ma (Fig. 3, 4, Table 2) and a crown node age estimate of 74 Ma (Fig. 4, Table 2). For other major families of lichenized ascomycetes much older crown node ages have been estimated (Rivas Plata & Lumbsch, pers. com.), including Graphidaceae (156 Ma), Physciaceae (153 Ma), and Ramalinaceae (126 Ma). The parmelioid clade is the largest and most strongly supported monophyletic group of the family. Within the parmelioid crown a total of seven major divergence events at different times have been found (Fig. 4), ranging from early Eocene (separation of Xanthoparmelia from the Parmotrema-clade) to Oligocene (separation of Flavoparmelia from Parmotrema and Canoparmelia crozalsiana-clade). The radiations of parmelioid genera were estimated to start at the end of Eocene (Melanelixia) and occurred during the Oligocene-Miocene for most of the genera.
Major clades of parmelioid lichens either show distinct distribution patterns of the clade or include numerous species with disjunct distributions, as in the genera Xanthoparmelia, Austroparmelina, Melanohalea, Parmelina, and Remototrachyna [11,17,19,41,[60][61][62][63][64]. Examples of disjunct distributions include Xanthoparmelia and Austroparmelina. Xanthoparmelia with ca. 800 species, the most speciose clade of Parmeliaceae, occurs worldwide, although in some cases the species delimitation has recently been challenged [65], and in other cases they have been shown much more restricted than previously known [61][62][63][64]. Despite being cosmopolitan, Xanthoparmelia has two main areas of distribution in arid regions of the Southern Hemisphere (Australia and Africa). A distribution pattern spanning over Australia and Africa cannot be explained by vicariance since these landmasses separated much earlier ( Fig. 5; 110-135 Ma) than the split of the Xanthoparmelia lineage from the Parmotrema-clade (37-58 Ma). The data, however, cannot reject vicariance as the reason for distribution ranges of species occurring in South America and Australia. The separation of these continents (35-52 Ma) happened at a similar time as the origin of this lineage.
Austroparmelina includes 13 species that occur in southern and eastern Australia, Tasmania and New Zealand. Two species have wider distribution areas, occurring also in South Africa (A. labrosa, A. pseudorelicina; [18]); and one species also is known from South America (A. labrosa, Chile; [66]). According to our estimates the genus separated from its sister lineages in the late Eocene. Thus the presence of Austroparmelina in South Africa is most plausibly explained as a result of long distance dispersal because separation of Africa took place much earlier (Fig. 5). As in the case of Xanthoparmelia, our data cannot reject vicariance events as an explanation for the presence of A. labrosa in South America and Australia.
In general the estimated ages of diversification events found in our analyses indicate that disjunct distribution patterns of Southern Hemisphere lineages cannot be explained by vicariance. Transoceanic long distance dispersal is the most plausible cause to explain these distribution patterns. This is especially true for taxa occurring in Africa and Australia. This is consistent with recent studies employing molecular clock approaches interpreting distribution patterns in other groups of fungi and plants [58][59][67][68][69][70][71][72][73][74][75][76], suggesting that long-distance dispersal is an important factor to explain the current distribution patterns of plants and Two calibration points based on fossil assignation are marked as C2 (crown node of Parmelia) and C3 (crown node of Alectoria). C3* is the alternative calibration with the Alectoria fossil assigned to the alectorioid crown. The Parmelioid group is highlighted by the grey box. Numbers inside circles refer to divergence nodes between the main clades. Grey bars show the 95% highest posterior density intervals (HPD). Detailed ages are given in Table 2. doi:10.1371/journal.pone.0028161.g004 Table 2. Estimated ages and substitution rates of the most recent common ancestors (MRCA) for the main clades obtained with partitioned BEAST analyses using three calibration points. fungi. In the case of lineages distributed at present in South America, Antarctica and Australia the estimated ages do not discard that vicariance could have resulted from the break-up of continents that had occurred 35-52 Ma. For genera with more restricted distribution ranges and species groups occurring mainly in the Holarctic, additional phylogeographical data are necessary to test biogeographical hypotheses. This is the case e.g., for Parmelina, a genus that occurs in areas with a Mediterranean climate in the Northern Hemisphere [11]. Its separation from Myelochroa, a genus with center of distribution in eastern Asia [77] but extends further into temperate regions, was estimated as having happened in the late Eocene. Similarly the recently described genus Remototrachyna, a Southeast Asian element with only one pantropical species [19], diverged from its sisterlineage Bulbothrix also in the late Eocene.
Our analyses suggest a complex relationship of diversification events and palaeoclimatic conditions. The origin of Parmeliaceae was estimated in the late Cretaceous when the climate was warmer than today, and when subtropical to temperate fauna and flora extended well into polar latitudes [78][79]. The first radiation of the parmelioid group (60 Ma) occurred just before the Early Eocene Climatic Optimum (52 Ma) (Fig. 6) when temperature and atmospheric CO 2 reached maximum levels [80]. During this time period, the general climate was warm and humid, associated with tectonic changes and volcanism [81][82]. Most of the parmelioid lineages, however originated during the Eocene cooling and Oligocene glaciation (Fig. 6). During the Eocene-Oligocene transition (33.5 Ma) a profound global climate shift took place changing the Cretaceous/early Palaeogene ''Green House'' conditions to ''Ice House'' conditions, with the growth of Antarctic ice sheets to approximately their modern size, and the appearance of Northern Hemisphere glacial ice [80]. The radiation of the genera took place at different times from the early Oligocene to the Pliocene. Melanelixia started to radiate during the early Oligocene. During this time, the general climate was characterized for cold conditions of the Oligocene glaciation (B, Fig. 6). However Melanohalea, Parmelia, Remototrachyna, Flavoparmelia, Xanthoparmelia and Myelochroa started to radiate between the Late Oligocene Warming. and the Mid-Miocene Climatic Optimum (C-D, Fig. 6). The radiation of other genera (Punctelia, Parmotrema, Cetrelia, Parmelina, Canoparmelia crozalsiana-clade) is estimated as having occurred at the transition from the middle Miocene to the Early Pliocene, before the climate became cooler, drier and seasonal at the end of Pliocene [80]. Table 3. Estimated ages of the most recent common ancestor (MRCA) of the main clades obtained in the alternative BEAST analyses.

Conclusions
Using three calibration points (one at the split of Lecanoromycetes and Eurotiomycetes, and the ages of two fossil lichens) we obtained the first dated phylogeny of parmelioid lichens and estimated the ages of divergence of the well-resolved lineages and main genera. The radiation of the Parmelioid occurred near the Cretaceous-Tertiary (K/T) boundary, before the climatic optima. These age estimates indicate that long-distance dispersal has played a major role in shaping the current distribution of the Southern Hemisphere parmelioid lichens and that continental drift of Gondwana landmasses and vicariance cannot explain the Africa-Australia disjunct patterns. The major genera originated during Eocene and Oligocene, and radiated during cooling periods at different times from the late Oligocene to early Pliocene.

Materials and Methods
Taxon sampling, sequence alignment and selection of substitution model A dataset of 225 specimens of Lecanoromycetes with complete sequences of nuLSU rDNA, mtSSU rDNA and the protein coding RPB1 gene, generated in previous studies [18,[83][84], was compiled for this study. GenBank accession numbers are listed in Table S1. The parmelioid lineages (the main focus of this work) were represented by 96 OTUs including 24 genera. The major lineages of Lecanoromycetes [85][86] and outgroups (Chaethothyriomycetidae and Dothideomycetes), represented by 129 OTUs were included in the analyses to prevent that uneven sampling across the tree could distort the apparent trend in speciation through time [87]. The sequences of each locus were aligned separately using Muscle 3.6 [88] and the ambiguous positions removed using Gblocks with default settings [89][90].
The general time reversible model including estimation of invariant sites (GTR+I+G) was selected by jModelTest v 0.1.1 [91][92] as the most appropriate nucleotide substitution model for the three separated loci.

Calibration of nodes and dating analysis
The divergence time analyses were performed using BEAST v.1.6.1 [93]. For the dating analysis it is recommended to use a user starting tree instead of the random starting tree built by BEAST. The latter is very likely to violate the temporal and/or topological constraints specified to calibrate divergence times, and cause an error when attempting to initiate the MCMC. For building this tree we checked the phylogenetic signal of our matrix running preliminary ML and Bayesian analyses using Garli 0.96 [94] and MrBayes 3.1.1 [95]. ML analyses were carried out with the default settings and Bayesian analyses were performed assuming a GTR+I+G model, run for 5 million of generations with 4 chains and every 100 th tree sampled. The first 5000 generations were burned in and a majority rule consensus tree was calculated with the sumt option. The topologies generated separately for each locus by ML and Bayesian analyses were congruent with the topology of the three loci concatenated, and with the general phylogeny obtained by Crespo et al. [18].
Three points of calibration were used for this study. The principal calibration point (C1) was the divergence time of 280-330 Ma for the stem of Lecanoromycetes following [26]. In addition we used the ages of two fossil lichens: the diversification node (C2) of Parmelia was calibrated with fossils from the Dominican amber (Parmelia ambra, 15-45 Ma, [50]), and the crown node (C3) of Alectoria with a fossil from the Baltic amber (35)(36)(37)(38)(39)(40) Ma, [49]). The assignation of fossils to extant groups is a crucial matter in dating analyses, and in the case of the lichens the fossil record is so sparse that this becomes particularly important. Parmelia ambra is a fossil from the Dominican amber resembling Parmelia saxatilis and similar species [50]. It presents unclear terminal pseudocyphellae, elongate isidia, plane to concave upper surface, and simple to dichotomously branched rhizines. All these features are characteristic of the Parmelia s. s. clade, and thus the fossil could also be related to the phylogenetically close genus Relicina (sister group of the Parmelia s.s.), the 'Parmelia signifera' group, or to the morphologically close Nipponoparmelia [18]. Relicina was discarded because P. ambra does not have cilia, a feature present on the Relicina species. The species of the 'Parmelia signifera' group have subsquarrose rhizines different from the simple to dichotomously branched rhizines of the fossil. On the other hand, Nipponoparmelia presents lobes rolled upwards [96] different from the flat lobes of the fossil. Moreover, only N. isidioclada has isidia but the rhizines are much branched, not simple or bifurcate, than those of the fossil. Thus, the Parmelia ambra fossil was used to calibrate the Parmelia s. s. crown node in the starting tree.
The Alectoria fossil from the Baltic amber [49] is morphologically related to the alectorioid clade (Bryoria, Pseudophebe and Alectoria; [35]). The fossil has abundant apothecia, a character that it shares with Alectoria, while the related Bryoria and Pseudephebe genera rarely have apothecia. Thus this fossil was used to constrain the crown node of Alectoria. Nevertheless, due to the inevitable uncertainty in Figure 6. Chronogram of parmelioid genera during Cenozoic time and its relationship to global temperature changes. The global temperature changes obtained from the deep-sea oxygen and carbon isotype proxies after Zachos et al [107].The major diversification events of parmelioid lichens are mapped onto the temperature curve: the square represents the parmelioid most recent common ancestor (MCRA), circles indicate splits of major lineages (numbers in the tree refer to those in Table 2 placement of fossil taxa we assessed the impact of individual fossil calibration on divergence time estimates using alternative analysis: 1) using a single fossil for calibration, either Parmelia (C2 node) or Alectoria (C3 node); and 2) using the Alectoria fossil to calibrate the whole alectorioid clade, including Alectoria, Bryoria and Pseudephebe (C3*).
The divergence time analyses were performed using BEAST v.1.6.1 [93]. We used as starting tree the ML tree obtained with Garli 0.96 [94], made ultrametric using nonparametric rate smoothing (NPRS) implemented in TreeEdit v.10a10 [97] with the divergence between Lecanoromycetes and Chaethothyriomycetidae set at 305 Ma. We constrained the position of Chaethothyriomycetidae as sister clade of Lecanoromycetes based on Schoch et al. [98]. Previously to the analysis, we test that this constraint is not significantly worse than the unconstrained topology, using the Shimodaira-Hasegawa test (SH) [99] and Expected Likelihood Weights test (ELW) [100]. Both tests were run on Tree-Puzzle 5.2 [101].
The final dating analysis was performed with a partitioned BEAST analysis with unlinked substitutions models (GTR+I+G) across the loci, a Birth-Death process tree prior, and a relaxed clock model (uncorrelated lognormal) for each partition. Calibration points were defined as prior distributions: 1) the split of Lecanoromycetes and Chaethothyriomycetidae (C1) was calibrated with a uniform distribution (280-330 Ma). 2) The calibrations points with fossils were considered as minimal ages and calibrated with a lognormal distribution [102]. The Parmelia crown node (C2) at log-normal mean = 2.77, offset = 14, lognormal standard deviation = 0.5. The Alectoria crown node (C3) at lognormal mean = 3.61, offset = 34, lognormal standard deviation = 0.75.
The final analysis was run for 10 million generations, with parameter values sampled every 1000 generation. We checked the stationary plateau with Tracer v. 1.4.1 [103]. We discarded 10% of the initial trees as burn in and the consensus tree was calculated using Tree Annotator v 1.6.1 [96]. The results were visualized with FigTree v. 1.3.1 [104]. The substitution rates for each locus were obtained running independent BEAST analyses for each dataset using the same parameters as in the partitioned analysis. Ages and rates were estimated for all the nodes with more than 0.95 of posterior probability both in the BEAST runs and in the previous Bayesian analysis.

Supporting Information
Table S1 Specimens used in this study with GenBank accession numbers. (DOC)