The rise of grasslands is linked to atmospheric CO2 decline in the late Palaeogene

The grassland biome (steppes, savannas, and prairies) covers vast areas of Earth's surface and accounts for as much as one-third of terrestrial net primary production. While grasses (Poaceae) dominate biomass, other families—especially daisies (Asteraceae)—are equally or more species-rich. The rise of grasslands fostered new landscapes and co-evolutionary radiations (e.g., grazing mammals). Biome ages are often inferred from fossil first appearances (e.g., grass phytoliths, pollen of daisies, grasses, amaranths), but fossil-based inferences can be complemented by fossil-calibrated molecular phylogenies, as done for tropical rainforests. Phylogenetic approaches to grassy biomes have largely focused on the origins of C4 grasslands; a broader view of grassland expansion timing using the dominant families Poaceae and Asteraceae has been lacking. Here, the authors estimate when grasslands first expanded using large time-calibrated phylogenies for Asteraceae (a newly assembled supertree) and Poaceae (the largest available tree), to detect temporal shifts in lineage diversification and to test correlations between diversification shifts and environmental fluctuations (atmospheric CO2 and global paleo-temperature). A major analytical challenge is sparse and taxonomically biased species sampling in these hyper-diverse clades. Existing diversification methods accommodate incomplete sampling but often assume unrealistic uniform missing-species distributions, which can bias rate estimates. The authors develop a Bayesian framework in RevBayes that: (i) models episodic (piecewise-constant) diversification through time with rate-shift episodes; (ii) incorporates an empirically informed, non-uniform distribution of missing species using taxonomic placements; and (iii) tests correlations between diversification rates and environmental variables using environmentally dependent diversification models.

Prior work has used fossil phytoliths and pollen records to infer the timing of open-habitat and grassland evolution, with evidence for Late Cretaceous microfossils and macroscopic fossils from the Eocene, and a substantial Oligocene–Miocene increase in diversity. Phylogenetic dating has illuminated the evolution of other biomes (e.g., Malpighiales, palms, Inga for tropical rainforests), but grassy biomes have been less explored with this approach, with most attention on C4 grassland origins. Earlier studies identified diversification shifts in Asteraceae (e.g., Vernonioid clade, Heliantheae alliance) and Poaceae, often during the Neogene, using methods such as MEDUSA, BAMM, and turboMEDUSA. Theoretical work has highlighted identifiability issues in diversification models for extant-only phylogenies when allowing arbitrarily complex rate functions, underscoring the need for constrained, piecewise-constant models and rigorous model comparison. Additionally, biased species sampling can strongly affect diversification inference, motivating methods that leverage taxonomic information on missing species.

Datasets and phylogenies: The Asteraceae supertree (2723 sampled species) was constructed from a fossil-calibrated plastid backbone chronogram (14 regions, 54 species across 13 subfamilies; Calyceraceae as outgroups) and eleven subfamily trees grafted onto the backbone. Fossil calibrations used: (i) an Eocene macrofossil capitulum plus associated pollen (45.6 Mya) for non-Barnadesioideae Asteraceae; (ii) Late Cretaceous fossil pollen (Tubulifloridites lilliei, 72.1 Mya) as a stem-group calibration for crown Asteraceae. Alignments via MAFFT; substitution model GTR+I with UCLN relaxed clock; tree prior birth–death; inference in RevBayes. The Poaceae supertree (3595 taxa) was taken from Spriggs et al., with two alternative calibration scenarios: #1 (younger; Eocene megafossil) and #2 (older; Cretaceous phytolith). Diversification framework: An episodic birth–death process was used, with speciation (λ) and extinction (µ) rates constant within predefined time intervals (epochs) but allowed to shift instantaneously at episode boundaries. Rates are homogeneous across lineages within a time interval. Time was discretized into equal-length epochs (tested: 4, 10, 20, 50, 100, 200). Three prior models for time-varying rates were compared: (i) UCLN (uncorrelated lognormal) where ln(λi), ln(µi) ~ Normal(m, σ); (ii) GMRF (Gaussian Markov random field) autocorrelated prior ln(λi) ~ Normal(ln(λi−1), σ); (iii) HSMRF (horseshoe MRRF) with local adaptivity ln(λi) ~ Normal(ln(λi−1), σγi), γi ~ half-Cauchy(0,1). These provide null models without environmental dependence. Empirical taxon sampling: A new sampling model incorporated taxonomic information by assigning known numbers of unsampled species to named clades and integrating over plausible times of their missing speciation events (between the clade stem age and present), extending prior approaches to the episodic model. This reduces bias from non-uniform sampling across clades. Environmentally dependent diversification: The authors modeled correlations between diversification rates and environmental variables (atmospheric CO2; global paleo-temperature). For CO2, they adapted the exponential dependence to the episodic framework: ln(λi) = ln(λ0) + β × ΔCO2i (and analogously for µi). Four models were evaluated: (1) fixed environmental linkage (no additional variation); (2) environmental linkage plus UCLN epoch-wise independent variation; (3) environmental linkage plus GMRF autocorrelated variation (Brownian motion with trend); (4) environmental linkage plus HSMRF autocorrelated, locally adaptive variation. If β = 0, each collapses to its corresponding episodic model without environmental dependence. Initial rate λ0 ~ Uniform(0,100). Temperature used analogous formulations. Environmental data: CO2 series from Beerling & Royer (updated by Dana Royer); temperature from Zachos et al. Environmental variables were averaged over 1-, 2-, and 5-Myr bins; sensitivity to epoch size assessed. Model selection and inference: For time-varying rate models and environmentally dependent models, marginal likelihoods were estimated via stepping-stone sampling (128 stones; 2000 iterations per stone; ~1374 moves/iteration) in RevBayes. Bayes factors (BFs) compared models. Support for the sign of β used posterior odds given a symmetric prior P(β<0) = P(β>0) = 0.5. Two taxon-sampling schemes (uniform vs empirical) were not directly compared via marginal likelihoods because they use different data; instead, a simulation study assessed bias and power. Simulation study: (a) Environment-dependent diversification: trees simulated under UCLN and GMRF with σ ∈ {0, 0.02, 0.04} and β ∈ {0, −0.005, −0.01} using TESS, spanning constant-rate, time-varying without environmental dependence, fixed environmental dependence, and combined time-varying with environmental dependence. Each setting: 10 trajectories and trees; analyzed under the four environmental models. (b) Empirical taxon sampling: missing species randomly added to the daisy tree, divergence times drawn under (i) constant-rate and (ii) time-varying episodic processes (rates from empirical estimates), then pruned to mimic empirical sampling; 100 trees per setting; analyzed with GMRF prior under both empirical and uniform sampling to assess false positives and power.

- Diversification timing: Both Asteraceae and Poaceae experienced their most dramatic increase in diversification rates from the late Oligocene (~28 Mya) to early Miocene (~20 Mya), with rates peaking between ~20–15 Mya and a brief decline ~13–10 Mya, followed by a later increase from the late Miocene (~10 Mya). - Calibration sensitivity: Using Poaceae calibration scenario #2 (Cretaceous phytolith), an earlier Poaceae diversification peak (~35–30 Mya) was detected; otherwise overall patterns were consistent. - Model robustness: Autocorrelated prior models (GMRF for daisies; HSMRF for grasses) were strongly favored over UCLN in time-varying analyses (per Bayes factors), and inferred the same major shifts regardless of epoch number. - Importance of sampling model: Assuming uniform taxon sampling biases diversification estimates. The empirical taxon-sampling approach, which assigns missing species to clades, improved power to detect true time variation and reduced false positives (supported by simulations). - Environmental correlations: There is decisive support for a negative correlation between atmospheric CO2 and diversification rates in both families across all four environmentally dependent models. Bayes factors supporting β < 0 were extremely strong: ~37,501 for the fixed, UC, and GMRF models; 49 for HSMRF. Temperature–diversification correlations were ambiguous and model-dependent, ranging from significantly positive to negative across models, with no consistent support in autocorrelated models. - Paleontological concordance: Low inferred diversification rates before ~35 Mya align with the scarcity of Asteraceae and Poaceae fossils; elevated rates during/after the Oligocene match increased fossil diversity. - LTT limitations: The principal increase in diversification for both families is not apparent in lineage-through-time (LTT) plots, highlighting limitations of LTT when sampling is highly non-uniform across clades. - Mid-Miocene dip (~13–10 Mya): A short-term decline in diversification may relate to a brief CO2 increase (not captured in smoothed curve) and/or the explosive radiation of hypsodont grazers (e.g., horses), potentially impacting grassland diversity and composition.

The study addresses when and why grassland-dominant lineages (Asteraceae, Poaceae) diversified by integrating fossil-calibrated phylogenies with Bayesian diversification models that account for non-uniform sampling and test environmental drivers. The simultaneous, tree-wide acceleration in diversification during the late Oligocene–early Miocene coincides with a major Cenozoic atmospheric CO2 decline to near-modern levels and associated global cooling (Coolhouse state). The decisive negative correlation between CO2 and diversification indicates that declining CO2 was associated with increased diversification of these grassland families, consistent with physiological and ecological theory: lower CO2 can limit plant performance and induce ecophysiological drought, favoring open habitats (grasslands and shrublands) over forests. Paleo-vegetation modeling during the Last Glacial Maximum similarly supports grassland expansion under low CO2. Temperature effects were inconclusive across models, suggesting that CO2 may have been a more direct or stronger macroevolutionary driver for these clades than temperature alone, or that intertwined variables (aridity, seasonality, fire) and biotic interactions complicate temperature signals. The mid-Miocene diversification dip aligns with the rapid expansion of hypsodont grazers, whose grazing can restructure communities, potentially reducing diversification in particular intervals, followed by later increases coincident with C4 grass expansion in warm, fire-prone environments and radiations of hyper-diverse Asteraceae lineages (e.g., Senecio). The findings underscore that global environmental changes—particularly CO2 decline—were key macroevolutionary drivers shaping the origin and expansion of grasslands. Methodologically, the use of an episodic, identifiable diversification model, rigorous Bayesian model comparison, and an empirical taxon-sampling scheme proved critical to recover robust temporal patterns and environmental correlations, mitigating biases from sparse, uneven sampling and over-flexible models.

This work shows that the rise and early diversification of grassland-dominant plant families (Asteraceae and Poaceae) were synchronous and tightly linked to the late Paleogene–early Neogene decline in atmospheric CO2, with major diversification increases beginning after ~34 Mya and peaking around 20–15 Mya. Correlations with paleo-temperature were not consistently supported. A brief mid-Miocene slowdown may reflect biotic impacts from hypsodont grazer radiations and/or short-term CO2 fluctuations, with subsequent increases associated with C4 grass expansion and radiations within Asteraceae. Methodologically, the study contributes a Bayesian, episodic, environmentally dependent diversification framework with an empirical taxon-sampling model implemented in RevBayes, enabling robust inference in hyper-diverse, sparsely sampled clades. Given projections of rising atmospheric CO2, the findings imply substantial future shifts in grassland composition and diversity, including potential declines in C4 lineages. Future directions include extending models to jointly assess multiple environmental drivers (e.g., CO2, temperature, aridity, fire regimes), integrating fossil occurrence data more directly, testing lineage-specific rate heterogeneity, and incorporating explicit biotic interactions (e.g., herbivore diversification) to parse causal mechanisms.

- Sampling and coverage: Hyper-diverse families necessitated sparse and uneven sampling; although the empirical taxon-sampling model mitigates bias, uncertainty remains about the exact placement and timing of missing species. - Model assumptions: The episodic model assumes tree-wide (homogeneous) diversification rates within epochs, not lineage-specific rate heterogeneity; true processes may include clade-specific dynamics. - Identifiability: While piecewise-constant models are identifiable under constraints, extant-only phylogenies can be consistent with multiple diversification histories; inference depends on model structure and priors. - Environmental proxies and resolution: CO2 and temperature reconstructions carry uncertainties and were averaged over million-year bins; short-lived fluctuations may be smoothed out. - Calibration uncertainty: Alternative Poaceae fossil calibration scenarios yield different timing of peaks; fossil placement debates (e.g., Cretaceous phytoliths) affect absolute timing. - Temperature effects: Ambiguous model-dependent results for paleo-temperature indicate limited power or confounding with other variables (aridity, fire, seasonality).