Phylogenomic analyses of echinoid diversification prompt a re-evaluation of their fossil record

  1. Nicolás Mongiardino Koch  Is a corresponding author
  2. Jeffrey R Thompson
  3. Avery S Hiley
  4. Marina F McCowin
  5. A Frances Armstrong
  6. Simon E Coppard
  7. Felipe Aguilera
  8. Omri Bronstein
  9. Andreas Kroh
  10. Rich Mooi
  11. Greg W Rouse
  1. Department of Earth & Planetary Sciences, Yale University, United States
  2. Scripps Institution of Oceanography, University of California San Diego, United States
  3. Department of Earth Sciences, Natural History Museum, United Kingdom
  4. University College London Center for Life’s Origins and Evolution, United Kingdom
  5. Department of Invertebrate Zoology and Geology, California Academy of Sciences, United States
  6. Bader International Study Centre, Queen's University, Herstmonceux Castle, United Kingdom
  7. Departamento de Bioquímica y Biología Molecular, Facultad de Ciencias Biológicas, Universidad de Concepción, Chile
  8. School of Zoology, Faculty of Life Sciences, Tel Aviv University, Israel
  9. Steinhardt Museum of Natural History, Israel
  10. Department of Geology and Palaeontology, Natural History Museum Vienna, Austria
9 figures, 2 tables and 1 additional file

Figures

Neognathostomate diversity and phylogenetic relationships.

(A) Fellaster zelandiae, North Island, New Zealand (Clypeasteroida). (B) Large specimen: Peronella japonica, Ryukyu Islands, Japan; Small specimen: Echinocyamus crispus, Maricaban Island, Philippines (Laganina: Scutelloida). (C) Large specimen: Leodia sexiesperforata, Long Key, Florida; Small specimen: Sinaechinocyamus mai, Taiwan (Scutellina: Scutelloida). (D) Rhyncholampas pacificus, Isla Isabela, Galápagos Islands (Cassidulidae). (E) Conolampas sigsbei, Bimini, Bahamas (Echinolampadidae). (F)Apatopygus recens, Australia (Apatopygidae). (G) Hypotheses of relationships among neognathostomates. Top: Morphology supports a clade of Clypeasteroida + Scutelloida originating after the Cretaceous-Paleogene (K-Pg) boundary, subtended by a paraphyletic assemblage of extant (red) and extinct (green) ‘cassiduloids’ (Kroh and Smith, 2010). Bottom: A recent total-evidence study split cassiduloid diversity into a clade of extant lineages closely related to scutelloids, and an unrelated clade of extinct forms (Nucleolitoida; Mongiardino Koch and Thompson, 2021d). Divergence times are much older and conflict with fossil evidence. Cassidulids and apatopygids lacked molecular data in this analysis. Scale bars = 10 mm.

Phylogenetic relationships among major clades of Echinoidea.

(A) Favored topology, as obtained using the full supermatrix and a best-fit partitioning scheme in IQ-TREE (Nguyen et al., 2015). With the exception of a single contentious node within Echinacea (marked with a yellow star), all methods supported the same pattern of relationships, and assigned maximum support values to all nodes. Numbers below major clades correspond to the current numbers of described living species (obtained from Kroh and Mooi, 2020). (B) Likelihood-mapping analysis showing the proportion of quartets supporting different resolutions within Echinacea. While the majority of quartets support the topology depicted in A (shown in red), a relatively large number support an alternative resolution that has been recovered in morphological analyses (shown in blue; Kroh and Smith, 2010). (C) Difference in likelihood score (delta likelihood) for the two resolutions of Echinacea most strongly supported in the likelihood-mapping analysis. Genes were sorted based on their inferred phylogenetic usefulness (Mongiardino Koch, 2021b), and gene-wise delta scores were averaged for datasets composed of multiples of 20 loci. Support for a clade of Salenioida + (Camarodonta + Stomopneustoida), as depicted in A, is seen as positive delta scores and is predominantly concentrated among the most phylogenetically useful loci. This signal is attenuated in larger datasets that contain less reliable genes, eventually favoring an alternative resolution (as seen by negative scores for the largest datasets). Only the 584 loci containing data for the three main lineages of Echinacea were considered. The line corresponds to a second-degree polynomial regression. (D) Resolution and bootstrap scores (see color scale) of the topology within Echinacea found using datasets of different sizes and alternative methods of inference.

Estimated branch lengths across different models of molecular evolution.

Different site-homogeneous models (left) infer similar levels of divergence, and the choice between them induces little distortion in the general tree structure. Site-heterogeneous models on the other hand not only infer a larger degree of divergence between terminals relative to site-homogeneous ones (center and right), but they also distort the tree (i.e., impose a non-isometric stretching), with branch lengths connecting outgroup taxa expanding much more than those within the ingroup clade.

Figure 4 with 1 supplement
The 10 most sensitive node dates are found within Cidaroidea, Aulodonta, Neognathostomata, and among outgroup nodes.

For each, the range shown spans the interval between the minimum and maximum ages found among the consensus topologies of the 80 time-calibrated runs performed.

Figure 4—figure supplement 1
Median ages for selected clades across the consensus trees of the 80 time-calibrated experiments performed.
Figure 5 with 5 supplements
Sensitivity of divergence time estimation to methodological decisions.

Between-group principal component analysis (bgPCA) was used to retrieve axes that separate chronograms based on the clock model (A), model of molecular evolution (B), and gene sampling strategy (C) employed. In the latter case, only the first two out of four bgPCA dimensions are shown. The inset shows the centroid for each loci sampling strategy, and the width of the lines connecting them are scaled to the inverse of the Euclidean distances that separates them (as a visual summary of overall similarity). The proportions of total variance explained are shown on the axis labels. The impact of the clock model is such that a bimodal distribution of chronograms can be seen even when bgPCA are built to discriminate based on other factors (as in C).

Figure 5—figure supplement 1
Sensitivity of divergence time estimation to the use of alternate prior distributions on calibrated nodes.

Between-group principal component analysis (bgPCA) was used to retrieve the axes of maximum discrimination between chronograms estimated by enforcing either uniform or Cauchy prior distributions. The proportion of total variance explained by this axis is shown on the label.

Figure 5—figure supplement 2
Distribution of posterior probabilities for node ages that show an average difference larger than 20 Myr depending on the choice of clock prior.
Figure 5—figure supplement 3
Distribution of posterior probabilities for node ages that show a maximum difference larger than 20 Myr depending on the gene sampling strategy.

The largest differences can be seen in the relatively younger ages of Ambulacraria and Echinodermata when using clock-like genes, and in the relatively older ages for some nodes within cidaroids and asteroids when using loci with high occupancy. Other sampling criteria largely agree on inferred node ages, as can also be seen in Figure 5C as short distances between their centroids in the chronospace.

Figure 5—figure supplement 4
Distribution of posterior probabilities for node ages that are the most affected by the choice of model of molecular evolution.

No node showed average differences larger than 20 Myr, so those suffering the biggest changes are shown instead.

Figure 5—figure supplement 5
Distribution of posterior probabilities for node ages that are the most affected by the choice of prior distributions on calibrated nodes.

No node showed average differences larger than 20 Myr, so those suffering the biggest changes are shown instead.

Figure 6 with 2 supplements
Divergence times among major clades of Echinoidea and other echinoderms.

(A) Consensus chronogram of the two PhyloBayes (Lartillot et al., 2013) runs using clock-like genes under a CAT + GTR + G model of evolution, an autocorrelated log-normal (LN) clock, and Cauchy prior distributions. Node ages correspond to median values, and bars show the 95% highest posterior density intervals. (B) Lineage-through-time plot, showing the rapid divergence of higher-level clades following the P-T mass extinction (shown with dashed lines, along with the Cretaceous-Paleogene [K-Pg] boundary). Each line corresponds to an individual consensus topology from among the 80 time-calibrated runs performed. (C) Posterior distributions of the ages of selected nodes (identified in A with numbers). The effects introduced by the use of different models of molecular evolution and node age prior distributions are not shown, as they represent the least important factors (see Figure 5); the posterior distributions obtained under different settings of these were merged for every combination of targeted loci and clock prior. Tick marks = 10 Myr.

Figure 6—figure supplement 1
Number of lineages inferred to have crossed the Permian-Triassic (P-T) boundary.

The probabilities of each scenario are estimated from the inferred divergence times of the 16,000 chronograms sampled across all of the analyses performed (200 for the two runs of each combination of sampled loci, model of molecular evolution, and clock model and prior distribution on node ages). The probability of three or more crown group lineages surviving the P-T extinction is 59.58%. Names show the combination of clades crossing the P-T boundary for each scenario of number of survivors.

Figure 6—figure supplement 2
Prior distributions of all constrained nodes.

Five-hundred replicates were sampled from the joint prior, showing appropriately broad distributions of node ages. Blue lines show minima and yellow ones maxima (when enforced); dotted lines show the age of the Permian-Triassic (251.9 Ma) and Cretaceous-Paleogene (66 Ma) mass extinction events. Nodes whose ages are of special interest (Echinoidea, Scutelloida, and Clypeasteroida) are shown in pink, revealing large prior probabilities of the divergence occurring at either side of mass extinction events.

Relationship between the root-to-tip variance (a proxy for the clock-likeness of loci) and the rate of evolution.

The most clock-like loci (shown in red), which are often favored for the inference of divergence times (e.g., Smith et al., 2018; Carruthers et al., 2020), are among the most highly conserved and can provide little information for constraining node ages (see also Mongiardino Koch, 2021b). Clock-like genes with a higher information content were used instead by choosing the loci with the lowest root-to-tip variance from among those that were within one standard deviation from the mean evolutionary rate (shown in blue).

Appendix 3—figure 1
Ordering of loci enforced using genesortR (Mongiardino Koch, 2021b) and its relationship to the seven gene properties employed.

High ranking loci (i.e., the most phylogenetically useful) show low root-to-tip variances (or high clock-likeness), low saturation, and low compositional heterogeneity, as well as high average bootstrap and Robinson-Foulds similarity to a target topology (in this case, with the contentious relationship among major lineages of Echinacea collapsed).

Appendix 3—figure 2
Trace plots of the log-likelihood of different time calibration runs.

All runs show evidence of reaching convergence and stationarity before our imposed burn-in fraction of 10,000 generations (dashed lines). For simplicity, only runs under the CAT + GTR + G model and Cauchy priors are plotted. Those run under uniform node age priors behaved identically, while those run under GTR + G converged much faster.

Tables

Appendix 1—table 1
Transcriptomic/genomic datasets added in this study relative to the taxon sampling of Mongiardino Koch et al., 2018 and Mongiardino Koch and Thompson, 2021d.

Taxa with citations were taken from the literature, and details can be found in the corresponding studies and associated NCBI BioProjects. Taxa are shown in alphabetical order; those identified with ‘OG’ are outgroup taxa (i.e., non-echinoids). Voucher specimens are deposited at the Benthic Invertebrate Collection, Scripps Institution of Oceanography (SIO-BIC), and the Echinoderm Collection, Muséum National d’Histoire Naturelle (MNHN-IE). If multiple accession numbers are shown for a given taxon, these datasets were merged before assembly. Similar details for all other specimens can be found in Mongiardino Koch et al., 2018 and Mongiardino Koch and Thompson, 2021d.

TaxonCitationTissue typeCollectorLocality (depth)VoucherGenBank accession number
Amphiura filiformis (OG)Dylus et al., 2016SRX2255774
Apatopygus recensThis studyMixedOwen AndersonFoveaux Strait, South Island, New ZealandSIO-BIC E7142SRR16134561
Aspidodiadema hawaiienseThis studyMixedGreg Rouse, Avery HileySeamount 4, Costa Rica, Pacific Ocean (1003 m)SIO-BIC E7363SRR16134560
Asterias rubens (OG)Reich et al., 2015SRX445860
Astrophyton muricatum (OG)Janies et al., 2016; Linchangco et al., 2017SRX1391908
Bathysalenia phoinissaThis studyTube feet, spine musclesBIOMAGLO Cruise TeamMozambique Channel, Mayotte (295–336 m)MNHN-IE-2016–23SRR15130003
Clypeaster japonicusThis studyEggsFrances ArmstrongMisaki Marine Biological Station, Kanagawa, JapanSRR16134552
Echinothrix calamarisThis studyMale gonadPhilippinesSRR16134551
Encope emarginataThis studyEggsGulf Specimen Marine LabApalachee Bay, FL, USASRR16134550
Eucidaris metulariaThis studyMixedGreg RouseAl-Fahal Reef, Red Sea, Makkah, Saudi ArabiaSIO-BIC 2017–008SRR16134549
Fellaster zelandiaeThis studySpines, tube feet, gut, joining tissue of Aristotle’s lanternWilma BlomWestern end of Cornwallis Beach, Auckland, New ZealandSIO-BIC E7920SRR16134548
Histocidaris variabilisThis studyGonadGreg Rouse, Avery HileyLas Gemelas Seamount, near Isla del Coco, Costa Rica, Pacific Ocean (571 m)SIO-BIC E7350SRR16134547
Lamberticrinus messingi (OG)Clouse et al., 2015; Janies et al., 2016SRX1397823
Leodia sexiesperforataThis studyEggsFrances ArmstrongBocas del Toro, PanamaSRR16134546
Metacrinus rotundus (OG)Koga et al., 2016DRX021520
Micropyga tuberculataThis studyTube feet, spine musclesBIOMAGLO Cruise TeamMozambique Channel, Îles Glorieuses (385–410 m)MNHN-IE-2016–39SRR15130004
Ophiocoma echinata (OG)Reich et al., 2015SRX445856
Peronella japonicaThis studyEggsFrances ArmstrongKanazawa, JapanSRR16134545
Rhyncholampas pacificusThis studyMixedCarlos A. Conejeros-VargasPanteón Beach, Puerto Ángel Bay, Oaxaca, MexicoSIO-BIC E7918SRR16134558
Scaphechinus mirabilisSimakov et al., 2015DRX187887
DRX187888
Sinaechinocyamus maiThis studyGregory’s diverticulumJih-Pai LinMiaoli County, TaiwanSIO-BIC E7919SRR16134557
Sterechinus neumayeriCollins et al., 2021ERX3498697
ERX3498698
ERX3498699
ERX3498700
ERX3498701
Stereocidaris nascaensisThis studyMuscle surrounding Aristotle’s lanternCharlotte SeidOff San Ambrosio, Desventuradas Islands, Chile (215 m)SIO-BIC E7154SRR16134559
Tetrapygus nigerThis studyFemale gonad, tube feetFelipe AguileraDichato Bay, ChileSRR16134553
SRR16134554
Tripneustes gratillaThis studyPedicellariaePhilippinesSRR16134556
Tromikosoma hispidumThis studyTube feetLisa Levin, Todd LitkeQuepos Plateau, Costa Rica, Pacific Ocean (2067 m)SIO-BIC E7251SRR16134555
Appendix 1—table 2
Details of molecular datasets and supermatrix.

Statistics for raw reads and assemblies are shown for datasets incorporated in this study relative to the sampling of Mongiardino Koch et al., 2018 and Mongiardino Koch and Thompson, 2021d, where similar statistics can be found for the other datasets. Taxa are shown in alphabetical order; those identified with ‘OG’ are outgroup taxa (i.e., non-echinoids). Novel datasets correspond to Illumina pair-end transcriptomes, except for two draft genomes (Bathysalenia phoinissa and Micropyga tuberculata). Mean insert size is expressed in number of base pairs. For transcriptomes, read pairs (initial) shows numbers input into Agalma (Dunn et al., 2013), that is, after processing with Trimmomatic (Bolger et al., 2014) or UrQt (Modolo and Lerat, 2015), while read pairs (retained) show those that passed the Agalma curation checks (including ribosomal, adapter, quality, and compositional filters), and represent the final number of read pairs used for assembly. For genomes, see information in the bioinformatic details above. Assemblies were further sanitized with either alien_index (Ryan, 2014) alone or in combination with CroCo (Simion et al., 2018), and the number of assembled transcripts/contigs shows the size of datasets after these curation steps. The number of loci shows the occupancy of terminals in the supermatrix output by Agalma (1346 loci at 70% occupancy), after which loci were further removed by TreeShrink (Mai and Mirarab, 2018), resulting in the final occupancy scores.

SpeciesMean insert sizeRead pairs (initial)Read pairs (retained)Assembled transcripts/contigsNumber of lociRemoved by TreeShrinkFinal occupancy
Abatus agassizii522538.4
Abatus cordatus497236.8
Acanthaster planci (OG)864863.6
Amphiura filiformis (OG)180.161,558,17352,206,727416,946761855.9
Apatopygus recens415.374,315,68756,898,850274,3801152585.2
Apostichopus japonicus (OG)849962.4
Araeosoma leptaleum1184287.8
Arbacia lixula1234191.6
Aspidodiadema hawaiiense228.7109,716,219104,518,714311,03210601177.9
Asterias rubens (OG)230.431,890,61325,495,009103,090805959.1
Asthenosoma varium1254792.6
Astrophyton muricatum (OG)195.925,191,95422,281,829149,146478535.1
Bathysalenia phoinissa154,120605644.5
Brissus obesus773057.4
Caenopedina hawaiiensis1094281.1
Clypeaster japonicus198.610,505,5209,298,31674,743829161.5
Clypeaster rosaceus1242292.1
Clypeaster subdepressus1233391.4
Colobocentrotus atratus1158086.0
Conolampas sigsbei1001574.0
Cystechinus c.f. giganteus661049.1
Dendraster excentricus337125.0
Diadema setosum305122.6
Echinarachnius parma1187388.0
Echinocardium cordatum1012275.0
Echinocardium mediterraneum1185188.0
Echinocyamus crispus709652.2
Echinometra mathaei1055078.4
Echinothrix calamaris203.968,871,08760,443,904251,7881252492.7
Encope emarginata206.212,015,88810,461,87066,2411076179.9
Eucidaris metularia321.638,802,15929,047,40183,318412230.5
Eucidaris tribuloides414030.8
Evechinus chloroticus1196188.8
Fellaster zelandiae313.962,619,79151,742,748168,5981269194.2
Heliocidaris erythrogramma931069.2
Histocidaris variabilis329.284,103,67273,884,180144,0441189588.0
Holothuria forskali (OG)743854.6
Lamberticrinus messingi (OG)195.022,989,82020,234,59759,049351325.9
Leodia sexiesperforata203.816,211,18214,174,47568,9641064179.0
Leptosynapta clarki (OG)403429.6
Lissodiadema lorioli770057.2
Loxechinus albus1265193.9
Lytechinus variegatus1257293.2
Mellita tenuis1002174.4
Meoma ventricosa1053078.2
Mesocentrotus nudus1247092.6
Metacrinus rotundus (OG)198.523,791,83220,900,34583,718642747.2
Micropyga tuberculata170,514415430.5
Ophiocoma echinata (OG)278.528,427,02624,836,025130,153712752.4
Paracentrotus lividus1266094.1
Patiria miniata (OG)740854.4
Peronella japonica208.317,696,28715,707,931106,1101043677.0
Prionocidaris baculosa794358.8
Psammechinus miliaris1173187.1
Rhyncholampas pacificus346.863,116,41352,160,741234,9101070279.3
Saccoglossus kowalevskii (OG)668649.2
Scaphechinus mirabilis401.710,169,1958,796,067127,664974472.1
Sinaechinocyamus mai297.260,208,17252,265,201164,6461233191.5
Sphaerechinus granularis1214190.1
Sterechinus neumayeri207.326,376,27921,308,840122,0011097181.4
Stereocidaris nascaensis327.069,962,49560,316,050121,5081200688.7
Stomopneustes variolaris908067.5
Strongylocentrotus purpuratus1266593.7
Tetrapygus niger202.563,040,54152,761,671163,0841287295.5
Tripneustes gratilla154.662,962,23957,015,692130,3761309197.2
Tromikosoma hispidum300.457,208,35747,909,038234,5021238791.5

Additional files

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Open citations (links to open the citations from this article in various online reference manager services)

Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)

  1. Nicolás Mongiardino Koch
  2. Jeffrey R Thompson
  3. Avery S Hiley
  4. Marina F McCowin
  5. A Frances Armstrong
  6. Simon E Coppard
  7. Felipe Aguilera
  8. Omri Bronstein
  9. Andreas Kroh
  10. Rich Mooi
  11. Greg W Rouse
(2022)
Phylogenomic analyses of echinoid diversification prompt a re-evaluation of their fossil record
eLife 11:e72460.
https://doi.org/10.7554/eLife.72460