Repeated species radiations in the recent evolution of the key marine phytoplankton lineage Gephyrocapsa

Abstract

Phytoplankton account for nearly half of global primary productivity and strongly affect the global carbon cycle, yet little is known about the forces that drive the evolution of these keystone microscopic organisms. Here we combine morphometric data from the fossil record of the ubiquitous coccolithophore genus Gephyrocapsa with genomic analyses of extant species to assess the genetic processes underlying Pleistocene palaeontological patterns. We demonstrate that all modern diversity in Gephyrocapsa (including Emiliania huxleyi) originated in a rapid species radiation during the last 0.6 Ma, coincident with the latest of the three pulses of Gephyrocapsa diversification and extinction documented in the fossil record. Our evolutionary genetic analyses indicate that new species in this genus have formed in sympatry or parapatry, with occasional hybridisation between species. This sheds light on the mode of speciation during evolutionary radiation of marine phytoplankton and provides a model of how new plankton species form.

Fig. 1

Global distribution of extant Gephyrocapsa. a Distribution of the Gephyrocapsa isolates and their corresponding flora (adapted from ref. and updated with data from ref. ) shown over modelled-satellite SST monthly climatology from Modis Aqua (2002–2017) plotted with SeaDAS. Relative range of Gephyrocapsa flora are denoted with roman numbers. In the legend, species names are ordered by decreasing relative abundance in each assemblage. Sediment sites used in this study are also shown: 709, 994, 1082, 1313 and 1314. b Scanning electron micrographs of the Gephyrocapsa strains used in this study (Scale = 2 μm)

Fig. 2

Evolution of extant Gephyrocapsa species. a Phylogenetic chronogram of extant Gephyrocapsa based on analysis of genome sequence data of 10 strains. Every node in the phylogeny has a posterior probability at 1. The dating of ancestral nodes is based on relaxed molecular clock calibrated with the first appearance of G. huxleyi (node 2; uniform prior 290–350 ka) in the fossil record. 95% Highest posterior density (HPD) intervals for ages are shown as grey bars. b Densitree plot based on 6278 phylogenies constructed for genomic fragments 10 kb long, and consensus species tree of extant Gephyrocapsa inferred from a multi coalescent method. All nodes are supported with 100% bootstrap. c Demographic modelling of Gephyrocapsa speciation history. The sizes of blue circles are proportional to estimated population sizes of extant and ancestral Gephyrocapsa species; the numbers inside blue circles (=θ/θ_Ghux) show population sizes relative to that in extant G. huxleyi. Direction of inferred interspecific gene flow is shown with purple arrows, with their widths and the number next to arrowhead showing gene flow intensity per generation. For estimates with confidence intervals see Supplementary Tables 2 and 3. d Stacked relative abundance of the Noëlaerhabdaceae assemblage including Gephyrocapsa and Reticulofenestra species over the last ~700 ka at site 1082 (Large Gephyrocapsa include G. oceanica; Medium, G. margereli and G. muellerae; Small, G. aperta and G. ericsonii). e Absolute number of the same assemblage in site 1082. f Global Δ sea surface temperature (ΔSST) over 700 ka, blue line is for mean variations, red and blue bands for 95% confidence intervals (respectively variability and jackknife). Red line corresponds to loess smoothing and associated grey band to 95% confidence interval

Fig. 3

Matsuoka–Okada cycles and diversification in Gephyrocapsa. a Size variation of coccoliths in Gephyrocapsa spp. in the last 1.8 Ma at site 709^,, with three coccolith size enlargement events (MO1–MO3) associated with fossil species. Shades of green colour correspond to the relative density of coccolith. A phylogenetic chronogram of extant Gephyrocapsa based on relaxed molecular clock and calibrated with the first appearance of G. huxleyi (node 2) in the fossil record. Branches of the tree are scaled to coccolith length, estimated by ancestral trait reconstruction (Supplementary Fig. 9), and shown for each node. Mean coccolith length (in μm) for modern and ancestral species is shown in white rectangle at each tip and node. b Coccolith traits measured. c Global Δ sea surface temperature (ΔSST) over 1800 ka. Blue line is for mean variations, red and blue bands for 95% confidence intervals (respectively, variability and jackknife). Red line corresponds to loess smoothing and its associated grey band to 95% confidence interval. White arrows correspond to stepwise changes in relative global SST