# Extinction intensity during Ordovician and Cenozoic glaciations explained by cooling and palaeogeography

Erin Saupe, Huijie Qiao, Yannick Donnadieu, Alexander Farnsworth, Alan Kennedy-Asser, Jean-Baptiste Ladant, Daniel Lunt, Alexandre Pohl, Paul Valdes, Seth Finnegan

<figure><img src="image_1.png" /></figure>

## To cite this version:

Erin Saupe, Huijie Qiao, Yannick Donnadieu, Alexander Farnsworth, Alan Kennedy-Asser, et al. Extinction intensity during Ordovician and Cenozoic glaciations explained by cooling and palaeogeography. Nature Geoscience, 2020, 13 (1), pp.65-70. 10.1038/s41561-019-0504-6. hal-02902805

HAL Id: hal-02902805

https://hal.science/hal-02902805v1

Submitted on 7 Oct 2025

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Extinction intensity during Ordovician and Cenozoic glaciations explained by cooling and palaeogeography

Erin E. Saupe<fnref n="1" />,<fnref n="8" />*, Huijie Qiao²,<fnref n="8" />, Yannick Donnadieu<fnref n="3" />, Alexander Farnsworth⁴, Alan T. Kennedy-Asser⁴, Jean-Baptiste Ladant⁵, Daniel J. Lunt⁴, Alexandre Pohl⁵,⁶, Paul Valdes⁴ and Seth Finnegan⁷,<fnref n="8" />

A striking feature of the marine fossil record is the variable intensity of extinction during superficially similar climate transitions. Here we combine climate models and species trait simulations to explore the degree to which differing palaeogeographic boundary conditions and differing magnitudes of cooling and glaciation can explain the relative intensity of marine extinction during greenhouse-icehouse transitions in the Late Ordovician and the Cenozoic. Simulations modelled the response of virtual species to cooling climate using a spatially explicit cellular automaton algorithm. We find that palaeogeography alone may be a contributing factor, as identical changes in meridional sea surface temperature gradients caused greater extinction in Late Ordovician simulations than in Cenozoic simulations. Differences in extinction from palaeogeography are significant, but by themselves are insufficient to explain observed differences in extinction intensity. However, when simulations included inferred changes in continental flooding and interval-specific models of sea surface temperature, predicted differences in relative extinction intensity were more consistent with observations from the fossil record. Our results support the hypothesis that intense extinction in the Late Ordovician is partially attributable to exceptionally rapid and severe cooling compared to Cenozoic events.

The fossil record documents substantial variation in extinction rates over the past 540 million years, punctuated by major extinction events<fnref n="1" />,². The degree to which major extinction events have common trigger mechanisms remains the subject of intense debate, but considerable evidence suggests that superfi- cially similar trigger mechanisms can have very different impacts depending on their magnitude and variation in the local and global boundary conditions under which they occur. For example, although several extinction events are thought to have been caused by emplacement of large igneous provinces (LIPs)<fnref n="3" />, the extreme severity of the Permian-Triassic mass extinction (252 million years ago (Ma)) may be attributable to the particular characteristics of the sedimentary rocks into which Siberian Traps volcanics were intruded⁴. The Late Ordovician mass extinction (LOME, ~445 Ma) has long posed a conundrum because it is the only major mass extinction thought to have been driven by climate cooling (although causes of the Late Devonian mass extinction remain uncertain and may include glaciation⁵). Late Ordovician extinctions cluster in two distinct pulses, with the first and largest appearing to coincide with cooling, expansion of south polar ice sheets, and sea level fall near the Katian-Hirnantian boundary, and the second appearing to coin- cide with warming and deglaciation in the mid-late Hirnantian⁶⁻<fnref n="8" />. Although more recent episodes of cooling and glaciation have been associated with elevated extinction of marine invertebrates⁹⁻¹³, none of these cooling events caused the magnitude of extinction char- acteristic of the LOME. This disparity has led to suggestions that

cooling alone cannot account for the LOME, and that other factors such as eutrophication, anoxia, heavy metal toxicity or eruption of a LIP may be the primary cause of the event¹⁴⁻²³.

Here we explore the degree to which differing boundary conditions may help to explain the very different extinction intensities associated with greenhouse-icehouse transitions in the Late Ordovician and in the Cenozoic. We focus in particular on a relatively unexplored factor—the role of palaeogeography and continental configuration. Previous work suggests that many marine ectotherms have thermal niches that are relatively stable on evolutionary timescales and rarely evolve at coarse scales in response to climate change, especially if change is relatively rapid²⁴⁻²⁶. Instead, the more common response of marine ectotherms to temperature change is habitat tracking²⁷⁻³¹—shifts in the geographic range of species to follow preferred abiotic/biotic conditions. The relative ability of entire faunas to track habitat can therefore provide a measure of faunal resilience to extinction³². We hypothesize that the geographic arrangement of continental crust may be important because it has the potential to facilitate or inhibit habitat tracking during periods of climate change (Fig. 1)³³.

To illustrate this, consider a geologically rapid temperature change of comparable magnitude occurring in a world with (1) a single continent that spans from the north to the south pole and (2) a single circumequatorial continent (Supplementary Fig. 1). In the former world, shallow marine species restricted to shelf and slope depths (which constitute the majority of the marine fossil record) would be

<footnote n="1">¹Department of Earth Sciences, University of Oxford, Oxford, UK. ²Institute of Zoology, Chinese Academy of Sciences, Chaoyang District, Beijing, China.</footnote>

<footnote n="3">³Aix-Marseille University, CNRS, IRD, INRA, Collège de France, CEREGE, Aix-en-Provence, France. ⁴School of Geographical Sciences, University of Bristol, Clifton, Bristol, UK. ⁵Department of Earth and Environmental Sciences, University of Michigan, Ann Arbor, MI, USA. ⁶Department of Earth Sciences, University of California, Riverside, Riverside, CA, USA. ⁷Department of Integrative Biology, University of California, Berkeley, Berkeley, CA, USA.</footnote>

<footnote n="8">⁸These authors contributed equally: Erin E. Saupe, Huijie Qiao, Seth Finnegan. *e-mail: erin.saupe@earth.ox.ac.uk</footnote>
<figure><img src="image_3.png" /></figure>

**Figure 1** | Per-cell proportional extinction for each greenhouse-icehouse transition. **a–c**, Simulated patterns are shown isolating the effect of palaeogeography (**a**), palaeogeography and sea level change (**b**) and palaeogeography, sea level and temperature change (**c**). The results show the average extinction across all dispersal and niche scenarios for the broad continental shelf. Temperature change was represented by the Fast Ocean Atmosphere Model (FOAM) for the Late Ordovician and Eocene-Oligocene, and by the Robertson HadCM3BL climate model for the Plio-Pleistocene (see Supplementary Table 3 for details). Note that the FOAM ocean-atmosphere general circulation models (AOGCMs) do not incorporate changes in sea level; thus, simulations incorporating estimates of temperature change would probably produce even greater extinction if sea level and temperature effects were modelled jointly in these AOGCMs.

able to disperse along the north-south oriented coastlines to track habitats within their thermal tolerance range (Fig. 2). Shallow marine species in the latter world would have no such option, and, if the new climate state included temperatures outside of their thermal tolerance range, would be driven to extinction. A similar phenomenon would be expected in a world with many island systems, especially if the island chains are dispersed longitudinally more than latitudinally. Thus, the arrangement of continental crust could play an important role in modulating the expected severity of extinction for a given climate perturbation. The arrangement of continents during the Late Ordovician differed strikingly from that during the Cenozoic (Fig. 1). Whereas the Cenozoic world is characterized by long, dominantly north-south oriented coastlines, the Late Ordovician world had only a single, north-south oriented supercontinent (Gondwana) and numerous island palaeocontinents and terranes³⁴.

There is a second potentially important difference in boundary conditions: whereas epicontinental seaways were of relatively limited extent during the Cenozoic, the Late Ordovician world was characterized by extensive continental flooding. Epicontinental seaways were, in many cases, characterized by temperature ranges and environmental conditions quite distinct from those of open oceans, limiting the ability of epicontinental endemics to colonize adjacent ocean-facing shelves³⁵. Thus, the potential effects of glacioeustatic sea level changes on the area and the availability of suitable shallow marine habitat may have been greater in the Late Ordovician than during the Cenozoic. Finally, available proxy data suggest that the Late Ordovician glaciation may have been more severe than those of the Cenozoic, both in terms of the volume of continental ice sheets (and associated eustatic changes) and the degree of cooling of surface waters³⁶⁻⁴².

We evaluate the influence of all three of these factors (palaeogeography, glacioeustasy and temperature) on expected extinctions of shallow marine taxa during three major cooling events in Earth history: Late Ordovician (~445 Ma), Eocene-Oligocene (~34 Ma) and Plio-Pleistocene (~4 Ma). Our goal is not to precisely model

the pattern of extinctions, but rather to evaluate whether the exceptionally high extinction during the first pulse of the LOME can be explained by any one or a combination of these factors, or instead requires additional explanation. We first isolate the expected effects of palaeogeography, holding all other factors constant. We then utilize a series of climate model simulations to examine how glacioeustatically forced sea level change and estimated temperature change can explain the differences in extinction magnitude observed among the three events.

## Results

**Effect of varying palaeogeography.** When applying the same temperature gradient to all models to isolate the effect of paleogeography, proportional extinction was higher by 0.12 using a Late Ordovician palaeogeographic configuration than when using late Eocene or late Pliocene palaeogeographies (Supplementary Table 5). The difference occurs for virtually all simulated dispersal, niche and continental shelf combinations (Fig. 3 and Supplementary Fig. 5). The only group of species that did not exhibit higher extinction in the Late Ordovician comprised those species with broad thermal niches and poor dispersal abilities, but in this case the differences in proportional extinction across the three palaeogeographies were minimal (<0.02; Supplementary Table 5). As expected, overall extinction was highest when assuming narrow shelves and when considering simulated species with narrow niches and poor dispersal abilities. Differences across palaeogeographic scenarios were largely consistent, regardless of whether a narrow or broad marine shelf was assumed (Supplementary Fig. 5 and Supplementary Table 5). Proportional extinction was highest on islands and in equatorial regions (Fig. 1 and Supplementary Fig. 6) and generally lowest along north-south coastlines (Fig. 4).

**Effect of varying palaeogeography and sea level.** When estimated changes in sea level and coastline position were included in simulations, proportional extinction was substantially higher
<figure><img src="image_4.png" /><figcaption>**Figure 2 | Schematic of dispersal and extinction in the simulation.**</figcaption></figure>

a–f, Cells represent a hypothetical portion of a gridded shallow marine shelf in which simulations are run for scenario 1 (a–c) and scenario 2 (d–f). Blue and green cells indicate unsuitable and suitable conditions for species, respectively. In a and d, the brachiopod image denotes the location of a seed cell where a simulation is initiated for scenario 1 (a) and scenario 2 (d). Once the simulation begins (b,e), the incipient species searches for a suitable habitat and disperses to other suitable cells under a warm, stable climate. The climate then cools (c,f), altering the cells suitable for the species. In scenario 1, the virtual species is capable of jumping over multiple, unsuitable cells to colonize newly suitable regions: this species survives the climatic change. In scenario 2, the virtual species cannot traverse as many unsuitable cells (that is, it has poor dispersal ability) and therefore becomes extinct.

in the Late Ordovician scenario (by up to 0.24) than the Eocene and Pliocene scenarios for all niche breadth and dispersal combinations, regardless of continental shelf extent (Supplementary Fig. 7 and Supplementary Table 5). Simulations run with broad marine margins returned slightly lower extinction compared to those run under narrow marine margins (Supplementary Fig. 7). Proportional extinction was highest on islands and in equatorial regions (Fig. 1 and Supplementary Fig. 6). Patterns were broadly similar for simulations that allowed species more than one opportunity to disperse to new cells following sea level fall, although extinction magnitudes were generally lower across all scenarios (Supplementary Fig. 8).

<figure><img src="image_5.png" /><figcaption>**Figure 3 | Proportional extinction results from the simulations.**</figcaption></figure>

Proportional extinction (expressed as percentages for ease of reading) resulting from simulations testing the effect of (1) palaeogeography only (first column), (2) palaeogeography and estimated changes in sea level (second column), and (3) palaeogeography, sea level and temperature change from AOGCMs (third column). Simulations were run in a shallow marine area directly adjacent to land, and only the narrow marine margin is represented here. AOGCM results are depicted for AOGCM model combinations (1–2) for the Late Ordovician, (5–8) for the Eocene–Oligocene and (10–12) for the Plio-Pleistocene (see Supplementary Table 3 for a list of the climate model codes). The Late Ordovician and the FOAM Eocene AOGCMs do not incorporate changes in sea level; thus, simulations incorporating estimates of temperature change would probably produce even greater extinction if sea level and temperature effects were modelled jointly in these AOGCMs. For other scenarios, see Supplementary Figs. 5, 7–12 and 15–17. Eoc, Eocene; Ord, Ordovician; Plio, Pliocene.

**Effect of varying temperature change.** When using interval-specific and palaeogeographically explicit surface temperature change estimates from AOGCMs, extinction magnitudes were substantially higher for the Late Ordovician scenarios than for the Eocene and Pliocene scenarios, in most cases by about a factor of 10 (Figs. 1 and 3, Supplementary Figs. 9–12 and Supplementary Table 5). These relative patterns held regardless of dispersal ability, marine margins or the AOGCMs used to characterize the greenhouse–icehouse transitions. Proportional extinction was dampened slightly for simulations run with broad marine margins compared to those run with narrow marine margins, but relative differences across models and types of species were similar.

The high extinction magnitude predicted for the Late Ordovician is largely a function of the large average temperature drop estimated by AOGCMs for this interval (Supplementary Fig. 13). This cooling substantially exceeds the cooling estimated for Cenozoic greenhouse–icehouse transitions. Although constraints on the magnitude of cooling during the Late Ordovician are few, modelled difference in tropical cooling between the Late Ordovician and Cenozoic transitions are generally consistent with available proxy data³⁶,⁴⁰,⁴²⁻⁴⁴ (Supplementary Fig. 14). When comparing the Cenozoic cooling
<figure><img src="image_6.png" /></figure>

**Figure 4** | Proportional extinction on north-south coastlines compared to east-west coastlines and 'islands'. Results are shown for simulations considering palaeogeography alone, for different niche and dispersal combinations. 'Islands' were defined as small isolated shallow marine shelf areas with contiguous cells. Error bars are 95% binomial confidence intervals. For predicted extinction patterns, colour and symbols indicate niche breadths (narrow and broad) and dispersal abilities (poor and good) of simulated species. As expected, extinction is higher on islands and east-west coastlines compared to north-south coastlines. However, many species occupied both north-south and east-west coastlines, and therefore proportional extinction differences are minimal. See Supplementary Methods and Supplementary Fig. 18 for additional details.

events, Robertson Plc models for the Eocene–Oligocene (AOGCM transition codes 4–7; Supplementary Table 3) and Plio–Pleistocene (transition codes 10–12; Supplementary Table 3) resulted in higher rates of extinction, consistent with the higher cooling predicted by these AOGCMs (Supplementary Fig. 13). Although the temperature declined less, on average, for the Eocene–Oligocene Getech model compared to the Eocene–Oligocene FOAM model (Supplementary Fig. 13), proportional extinction was generally lower for the latter. The Eocene–Oligocene FOAM model returned lower extinction, even given greater climate change, because glacioeustatic sea level did not change for this AOGCM (although coastlines were artificially changed for simulations focused on glacioeustacy, see Supplementary Methods).

**The effect of shallow marine area definition.** Patterns were broadly consistent when simulations were run in an area defined by every cell within 0–200 m water depth (Supplementary Fig. 4), accounting for the extent of continental flooding. Extinction magnitude was higher for the Late Ordovician in virtually all simulated combinations when considering (1) the effect of palaeogeography alone (Supplementary Fig. 15), (2) palaeogeography and sea level change (Supplementary Fig. 16), and (3) interval-specific and palaeogeographically explicit surface temperature change estimates from AOGCMs (Supplementary Fig. 17).

Discussion

Our simulations suggest that palaeocontinental configuration alone could explain a small part of the observed variation in extinction intensity between the Late Ordovician and Cenozoic greenhouse-icehouse transitions: identical changes in sea surface temperature caused higher extinction with a Late Ordovician continental configuration than with Cenozoic continental configurations. This difference reflects both the relative paucity of long north-south oriented coastlines spanning large latitudinal ranges and the relative abundance of small islands and terrains in the Late Ordovician compared to the Cenozoic (Supplementary Fig. 18). As the climate cooled, the proportion of species able to disperse along continental margins without crossing large areas of open ocean may have been higher in the Cenozoic than in the Late Ordovician. Our simulations do not account for the many dispersed island arcs and terranes with poorly constrained locations in the Late Ordovician that are not included in simplified palaeogeographic reconstructions. Many of these terranes experienced very high apparent extinction rates in the Late Ordovician, although it is not always easy to disentangle the effects of the mass extinction and continental amalgamation¹⁷,⁴⁵.

Although interesting, the predicted differences in extinction intensity as a consequence of palaeogeography alone are too small to account for the exceptional severity of the LOME compared to Cenozoic cooling events. Modelling changes in sea level substantially amplified the expected difference in extinction intensity between the Late Ordovician and Cenozoic. Simulations that also varied the magnitude and pattern of sea surface cooling predicted Late Ordovician extinction magnitudes far exceeding those associated with Cenozoic events (Figs. 1 and 3).

It is difficult to compare the geographic pattern of extinctions from our simulations with empirical patterns for a variety of reasons: our simulations make the unrealistic assumption that species with different ecological attributes are uniformly distributed on continental margins and do not account for possible latitudinal diversity gradients. The biogeographic distributions of species with different ecological attributes in the past are largely unknown, and the geographic distributions of extinct taxa may be obscured and distorted by spatiotemporal variation in sedimentary rock deposition, preservation, exposure and sampling. However, proportional extinction of marine invertebrate genera across palaeolatitudes during the first (end-Katian) pulse of the LOME shows some similarity to predicted patterns, notably relative peaks in the tropics and at high southern latitudes (Supplementary Fig. 6).

Although thermal niche breadth cannot be measured easily for extinct species, it is notable in this context that thermal tolerance range and latitudinal range are highly correlated in modern marine ectotherms⁴⁶, and the first pulse of the LOME was characterized by exceptionally high extinction of brachiopod genera with relatively narrow (<35°) absolute palaeolatitudinal ranges⁴⁷. Genera with ranges largely or entirely limited to epicontinental seaways also experienced elevated extinction risk⁴⁷,⁴⁸. When considered in light of these patterns, the results of our simulations bolster the argument that the combined effects of palaeogeographic configuration, exceptional cooling, and draining of epicontinental seaways might help to explain the severity of the first pulse of the LOME relative to Cenozoic greenhouse-icehouse transitions.

Our analysis focused only on palaeogeography, sea level drop and temperature change, but other factors almost certainly worked in concert to elevate extinction intensity during the Late Ordovician, as evidenced by the fact that extinction rates substantially exceeded those of the Cenozoic even in regions that our simulations suggest should have experienced comparatively low extinction (Supplementary Fig. 6). For example, preferential extinction of deeper-water brachiopods⁴⁷,⁴⁹ and other benthic taxa suggests widespread anoxic and euxinic conditions may have played an
important role in the extinctions<sup>14,16,17,20</sup>. Due to the dependence of both organismal metabolic rates and oxygen availability on ambient temperature, changes in water temperature and oxygenation may have interactive effects on the distribution of suitable habitat for many species<sup>30,51</sup>.

Climate changes are linked intimately to changes in ocean and atmospheric circulation, which themselves are influenced by con- tinental configuration. Thus, disentangling the relative contribu- tions of these variables in empirical systems is difficult, which is why we implemented a simulation framework. We reiterate that our goal was not to create the most realistic simulations, but rather to generate models that isolate the effects of various factors. We used multiple AOGCMs to characterize uncertainty generated by model AOGCM choice, obtaining congruent patterns across model combi- nations. Whether surface temperatures predicted by AOGCMs are, in fact, representative of the simulated intervals is an open question pending more and better proxy constraints, but predicted changes in tropical sea surface temperatures are generally consistent with available data (Supplementary Fig. 14).

Another potential source of error is palaeogeography. Simulations were performed at relatively coarse spatial resolution (1°×1°), and therefore removed small islands and smoothed com- plex geographic features, some of which could have functioned as important biogeographic barriers. We do not, however, anticipate that this smoothing affected proportional extinction in a systematic fashion across events, and higher-resolution reconstructions would not reflect true (low) confidence in palaeogeography during these time periods. Although Cenozoic palaeogeography is comparatively well constrained, there are many uncertainties remaining regard- ing the positions of palaeocontinents and terranes during the Late Ordovician<sup>34,45,52</sup>. However, the relative positions of major continen- tal landmasses, especially those with well-sampled fossil records, are similar across most recent reconstructions, and we do not anticipate that using alternative palaeogeographic models would have a major impact on our findings.

Our simulation framework assumes that species thermal niches are highly conservative and do not evolve in response to climate change. Although overly simplistic, we regard this as a reasonable assumption based on previous research that has shown species' tolerances are stable over millions of years<sup>24</sup>. Niche evolution, if included, would probably dampen extinction rates in simulations, but we do not expect it would systematically affect extinction mag- nitudes for any one climate transition over another, unless niches were more or less labile through time. The static niches employed in our simulation framework resulted in a zone of very high tropical extinction for warm-adapted species with narrow thermal niches. These simulated species are doomed to extinction when climate cools as no suitable habitats remain available to them. Our simu- lations suggest that the proportion of species unable to disperse to suitable habitat may have been higher during Late Ordovician cooling than during Cenozoic cooling events, both because of the relative lack of continuous north-south coastlines and the presence of numerous small islands and terranes, and because a greater pro- portion of shallow marine habitat area was located in the tropics. Whether tropical regions have, in fact, consistently experienced exceptionally high extinction rates during greenhouse-icehouse transitions is unclear<sup>10,12,13,48,53</sup>.

All of our simulations assumed the same rate of transition from greenhouse to icehouse states—simulated species were forced to track suitable habitat as climate changed from warm to cold in only a set number of ‘time steps’. In reality, these climatic transitions occurred over timescales of 10<sup>3</sup>–10<sup>6</sup> years, which may have allowed species to more effectively track suitable habitat. Simulating more gradual climate transitions would, in general, reduce proportional extinction by permitting species to more effectively track preferred temperatures via small dispersal steps, but whether this reduction

would apply equally to all simulated transitions is unclear. At pres- ent it is difficult to constrain the rate of cooling during the Late Ordovician greenhouse-icehouse transition, and hence a mean- ingful comparison of rates across all three events is not possible. For this reason, it is also difficult to account for the potential effects of changes in palaeogeographic configuration (for example, accretion of peri-Laurentian terranes in the Late Ordovician<sup>45</sup>, uplift of the Isthmus of Panama in the Plio-Pleistocene<sup>54</sup>) that may have occurred over the same timeframe as greenhouse-icehouse climate transitions.

Despite these caveats, our simulation framework provides clear evidence that the expected effects of greenhouse-icehouse tran- sitions differ depending on palaeogeographic configuration, sea level change and magnitude of temperature change. We suggest that this framework represents a promising way forward to pro- vide further constraints on the relative contribution of boundary conditions to extinction and speciation rates, particularly as more sophisticated models and better proxy/palaeogeographic data become available.

Online content

Any methods, additional references, Nature Research reporting summaries, source data, extended data, supplementary informa- tion, acknowledgements, peer review information; details of author contributions and competing interests; and statements of data and code availability are available at https://doi.org/10.1038/s41561- 019-0504-6.

Received: 12 July 2019; Accepted: 13 November 2019;

Published online: 16 December 2019

References

1. Bambach, R. K. Phanerozoic biodiversity mass extinctions. *Annu. Rev. Earth Planet. Sci.* **34**, 127–155 (2006).

2. Raup, D. M. & Sepkoski, J. J. Mass extinctions in the marine fossil record. *Science* **215**, 1501–1503 (1982).

3. Bond, D. P. G. & Wignall, P. B. Large igneous provinces and mass extinctions: an update. *Geol. Soc. Am. Spec. Pap.* **505**, 29–55 (2014).

4. Svensen, H. et al. Siberian gas venting and the end-Permian environmental crisis. *Earth Planet. Sci. Lett.* **277**, 490–500 (2009).

5. Joachimski, M. M. & Buggisch, W. Conodont apatite δ¹⁸O signatures indicate climatic cooling as a trigger of the Late Devonian mass extinction. *Geology* **30**, 711–714 (2002).

6. Brenchley, P. J. et al. High-resolution stable isotope stratigraphy of Upper Ordovician sequences: constraints on the timing of bioevents and environmental changes associated with mass extinction and glaciation. *Geol. Soc. Am. Bull.* **115**, 89–104 (2003).

7. Melchin, M. J., Mitchell, C. E., Holmden, C. & Šorch, P. Environmental changes in the Late Ordovician–Early Silurian: review and new insights from black shales and nitrogen isotopes. *Geol. Soc. Am. Bull.* **125**, 1635–1670 (2013).

8. Sheehan, P. M. The Late Ordovician mass extinction. *Annu. Rev. Earth Planet. Sci.* **29**, 331–364 (2001).

9. Edie, S. M., Huang, S., Collins, K. S., Roy, K. & Jablonski, D. Loss of biodiversity dimensions through shifting climates and ancient mass extinctions. *Integr. Comp. Biol.* **58**, 1179–1190 (2018).

10. Hansen, T. A. Extinction of Late Eocene to Oligocene molluscs: relationship to shelf area, temperature changes and impact events. *Palaios* **2**, 69–75 (1987).

11. Ivany, L. C., Patterson, W. P. & Lohmann, K. C. Cooler winters as a possible cause of mass extinctions at the Eocene/Oligocene boundary. *Nature* **407**, 887–890 (2000).

12. Powell, M. G. Timing and selectivity of the Late Mississippian mass extinction of brachiopod genera from the Central Appalachian Basin. *Palaios* **23**, 525–534 (2008).

13. Stanley, S. M. Anatomy of a regional mass extinction: Plio-Pleistocene decimation of the Western Atlantic bivalve fauna. *Palaios* **1**, 17–36 (1986).

14. Bartlett, R. et al. Abrupt global-ocean anoxia during the Late Ordovician– Early Silurian detected using uranium isotopes of marine carbonates. *Proc. Natl Acad. Sci. USA* **115**, 5896–5901 (2018).

15. Ghienne, J.-F. et al. A Cenozoic-style scenario for the End-Ordovician glaciation. *Nat. Commun.* **5**, 4485 (2014).

16. Hammarlund, E. U. et al. A sulfidic driver for the End-Ordovician mass extinction. *Earth Planet. Sci. Lett.* **331–332**, 128–139 (2012).
17. Harper, D. A. T., Hammarlund, E. U. & Rasmussen, C. M. Ø. End Ordovician extinctions: a coincidence of causes. *Gondwana Res.* **25**, 1294–1307 (2014).

18. Vandenbroucke, T. R. A. et al. Metal-induced malformations in Early Palaeozoic plankton are harbingers of mass extinction. *Nat. Commun.* **6**, 7966 (2015).

19. Zhou, L. et al. Changes in marine productivity and redox conditions during the Late Ordovician Hirnantian glaciation. *Palaeogeogr. Palaeoclimatol. Palaeoecol.* **420**, 223–234 (2015).

20. Zou, C. et al. Ocean euxinia and climate change 'double whammy' drove the Late Ordovician mass extinction. *Geology* **46**, 535–538 (2018).

21. Melott, A. L. et al. Did a gamma-ray burst initiate the Late Ordovician mass extinction? *Int. J. Astrobiol.* **3**, 55–61 (2004).

22. Gong, Q. et al. Mercury spikes suggest volcanic driver of the Ordovician–Silurian mass extinction. *Sci. Rep.* **7**, 5304 (2017).

23. Jones, D. S., Martini, A. M., Fike, D. A. & Kaiho, K. A volcanic trigger for the Late Ordovician mass extinction? Mercury data from South China and Laurentia. *Geology* **45**, 631–634 (2017).

24. Saupe, E. et al. Macroevolutionary consequences of profound climate change on niche evolution in marine molluscs over the past three million years. *Proc. R. Soc. Lond.* **B 281**, 20141995 (2014).

25. Peterson, A. T. Ecological niche conservatism: a time-structured review of evidence. *J. Biogeogr.* **38**, 817–827 (2011).

26. Stigall, A. L. When and how do species achieve niche stability over long time scales? *Ecography* **37**, 1123–1132 (2014).

27. Poloczanska, E. S. et al. Global imprint of climate change on marine life. *Nat. Clim. Change* **3**, 919–925 (2013).

28. Sunday, J. M., Bates, A. E. & Dulvy, N. K. Thermal tolerance and the global redistribution of animals. *Nat. Clim. Change* **2**, 686–690 (2012).

29. Clarke, A. & Crame, J. A. The Southern Ocean benthic fauna and climate change: a historical perspective. *Philos. Trans. R. Soc. B 338*, 299–309 (1992).

30. Barnes, D. K., Griffiths, H. J. & Kaiser, S. Geographic range shift responses to climate change by Antarctic benthos: where we should look. *Mar. Ecol. Prog. Ser.* **393**, 13–26 (2009).

31. Valentine, J. W. *Paleoecologic Molluscan Geography of the Californian Pleistocene* 134 (Univ. California Press, 1961).

32. Brett, C. E., Hendy, A. J., Bartholomew, A. J., Bonelli, J. R. Jr & McLaughlin, P. I. Response of shallow marine biotas to sea-level fluctuations: a review of faunal replacement and the process of habitat tracking. *Palaios* **22**, 228–244 (2007).

33. Valentine, J. W. *Evolutionary Paleocology of the Marine Biosphere* (Prentice Hall, 1973).

34. Cocks, L. R. M. & Torsvik, T. H. Earth geography from 500 to 400 million years ago: a faunal and palaeomagnetic review. *J. Geol. Soc. Lond.* **159**, 631–644 (2002).

35. Stanley, S. M. Thermal barriers and the fate of perched faunas. *Geology* **38**, 31–34 (2010).

36. Lear, C. H., Bailey, T. R., Pearson, P. N., Coxall, H. K. & Rosenthal, Y. Cooling and ice growth across the Eocene–Oligocene transition. *Geology* **36**, 251–254 (2008).

37. Pusz, A. E., Thunell, R. C. & Miller, K. G. Deep water temperature, carbonate ion and ice volume changes across the Eocene–Oligocene climate transition. *Paleoceanogr. Paleoclimatol.* **26**, PA2205 (2011).

38. Miller, K. G. et al. Climate threshold at the Eocene–Oligocene transition: Antarctic ice sheet influence on ocean circulation. *Geol. Soc. Am. Spec. Pap.* **452**, 169–178 (2009).

39. Coxall, H. K., Wilson, P. A., Pälike, H., Lear, C. H. & Backman, J. Rapid stepwise onset of Antarctic glaciation and deeper calcite compensation in the Pacific Ocean. *Nature* **433**, 53–57 (2005).

40. Finnegan, S. et al. The magnitude and duration of Late Ordovician–Early Silurian glaciation. *Science* **331**, 903–906 (2011).

41. Trotter, J. A., Williams, I. S., Barnes, C. R., Lécuyer, C. & Nicoll, R. S. Did cooling oceans trigger Ordovician biodiversification? Evidence from conodont thermometry. *Science* **321**, 550–554 (2008).

42. Liu, Z. H. et al. Global cooling during the Eocene–Oligocene climate transition. *Science* **323**, 1187–1190 (2009).

43. Herbert, T. D., Peterson, L. C., Lawrence, K. T. & Liu, Z. Tropical ocean temperatures over the past 3.5 million years. *Science* **328**, 1530–1534 (2010).

44. Katz, M. E. et al. Stepwise transition from the Eocene greenhouse to the Oligocene icehouse. *Nat. Geosci.* **1**, 329–334 (2008).

45. Rasmussen, C. M. Ø. & Harper, D. A. T. Did the amalgamation of continents drive the End Ordovician mass extinctions? *Palaeogeogr. Palaeoclimatol. Palaeoecol.* **311**, 48–62 (2011).

46. Sunday, J. M., Bates, A. E. & Dulvy, N. K. Global analysis of thermal tolerance and latitude in ectotherms. *Proc. R. Soc. Lond. B* **278**, 1823–1830 (2010).

47. Finnegan, S., Rasmussen, C. M. Ø. & Harper, D. A. T. Biogeographic and bathymetric determinants of brachiopod extinction and survival during the Late Ordovician mass extinction. *Proc. R. Soc. B* **283**, 20160007 (2016).

48. Finnegan, S., Heim, N. A., Peters, S. E. & Fischer, W. W. Climate change and the selective signature of the Late Ordovician mass extinction. *Proc. Natl Acad. Sci. USA* **109**, 6829–6834 (2012).

49. Finnegan, S., Rasmussen, C. M. Ø. & Harper, D. A. Identifying the most surprising victims of mass extinction events: an example using Late Ordovician brachiopods. *Biol. Lett.* **13**, 20170400 (2017).

50. Deutsch, C., Ferrel, A., Seibel, B., Pörtner, H.-O. & Huey, R. B. Climate change tightens a metabolic constraint on marine habitats. *Science* **348**, 1132–1135 (2015).

51. Penn, J. L., Deutsch, C., Payne, J. L. & Sperling, E. A. Temperature-dependent hypoxia explains biogeography and severity of end-Permian marine mass extinction. *Science* **362**, eaat1327 (2018).

52. Torsvik, T. H. & Cocks, L. R. M. *Earth History and Palaeogeography* (Cambridge Univ. Press, 2016).

53. Reddin, C. J., Kocsis, Á. T. & Kiessling, W. Climate change and the latitudinal selectivity of ancient marine extinctions. *Paleobiology* **45**, 70–84 (2019).

54. O'Dea, A. et al. Formation of the Isthmus of Panama. *Science* **2**, e1600883 (2016).

**Publisher's note** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

© The Author(s), under exclusive licence to Springer Nature Limited 2019
## Methods

**Simulation overview.** A full description of our methods is provided in the Supplementary Information. We used a gridded simulation framework⁵⁵⁻⁵⁷ to examine the potential effects of palaeogeography, glacioeustatically forced sea level fall and temperature change on extinction magnitude across three greenhouse to icehouse transitions: Late Ordovician, Eocene–Oligocene and Plio–Pleistocene. The gridded world in these simulations consisted of shallow marine margins around continents (both broad and narrow; Fig. 1 and Supplementary Fig. 2). Virtual species occupied cells in this gridded world if the average surface temperature of the cell fell within the species' niche (that is, was suitable) and searched for suitable cells based on assigned dispersal ability (Supplementary Methods). We simulated four types of virtual species characterized by different conditions of thermal niche breadth (narrow or broad) and dispersal ability (poor and good). Sea surface temperatures are derived from either simplified hypothetical meridional gradients (Supplementary Methods) or interval-specific AOGCMs (Supplementary Section 3). Each simulation begins in warm 'greenhouse' conditions, during which virtual species are seeded and allowed time to disperse and colonize suitable habitat. Following this burn-in period, cold 'icehouse' conditions are imposed, changing the extent and distribution of suitable cells for each species (Fig. 2 and Supplementary Methods). Species unable to disperse to suitable habitat under the colder conditions become extinct, whereas species with suitable habitat within their dispersal range survive (Fig. 2).

**Palaeogeography.** To isolate the influence of palaeogeography on extinction potential, we used a single set of simplified 'greenhouse' and 'icehouse' surface temperature gradients (Supplementary Methods and Supplementary Fig. 3) and varied only palaeogeography across the three greenhouse-icehouse transitions (Late Ordovician, Eocene–Oligocene and Plio-Pleistocene). Palaeogeography, specifically the configuration of shallow marine margins (Fig. 1 and Supplementary Fig. 2), was estimated from reconstructions in ref. ⁵⁸ for the Late Ordovician, ref. ⁵⁹ for the late Eocene (Ypresian) and Robertson Plc for the Pliocene⁶⁰. Virtual species were generated (that is, simulations initiated) at every possible shallow marine cell bordering land in each palaeogeographic reconstruction, for a total of 2,233, 2,258 and 2,892 virtual species in the Ordovician, Eocene and Pliocene, respectively. We tabulated the number of virtual species that became extinct when the climate cooled and compared proportional extinction for each palaeogeographic configuration (Supplementary Methods). This process was repeated for each possible combination of ecological and shelf parameters (niche breadth × dispersal ability × marine shelf extent).

**Eustatic sea level.** The second set of simulations focused on the effect of changing sea level on habitat reachability and hence extinction magnitude across the three greenhouse-icehouse transitions. When climate transitioned to 'icehouse' conditions, we changed the position of shallow marine margins to simulate the effect of glacioeustatic drawdown (Supplementary Methods). We employed the same hypothetical meridional temperature gradients used in the first set of simulations, such that palaeogeography and sea level change, but not climate change, varied across the three climate transitions. Shallow marine margins from the first set of simulations were used to model greenhouse conditions, and marine margins during icehouse conditions were estimated by either assuming a globally uniform sea level drop (for the Ordovician and Eocene) or based on specific palaeogeographic reconstructions (Pleistocene)⁶¹ (Supplementary Methods). Aside from shifting the position of shallow marine habitat (that is, sea level change), simulations differed from the previous round only in that the probability of a species dispersing to the new coastline following sea level drop was specified as a function of both its dispersal ability and the distance from the prior coastline (Supplementary Table 1). For simplicity and because we have few constraints on rates of sea level drop, we invoked 'jump' dispersal to new coastlines, rather than attempting to estimate rates of sea level change; that is, species were allowed only one opportunity to jump to the changed marine shelf position. These simulations used the same virtual species and ecological and shelf parameters as the first set of simulations, and extinction magnitudes across simulations were compared using the same framework.

**Temperature.** The third set of simulations evaluated the potential role of differences in degree of cooling on extinction magnitude across the three greenhouse-icehouse transitions (Supplementary Methods). These simulations employed surface temperature layers from interval-specific AOGCMs instead of the hypothetical meridional gradients used in previous simulations, and therefore modelled the combined effects of variation in palaeogeography, magnitude of eustatic changes, and magnitude and distribution of sea surface cooling. Note, however, that the FOAM AOGCMs for the Ordovician and Eocene did not incorporate change in eustatic sea level. Temperature changes for the Late Ordovician were represented by one AOGCM combination, for the Eocene-Oligocene by three AOGCM combinations, and for the Plio-Pleistocene by two AOGCM combinations, based on stage-level model availability (models and parameters are described in Supplementary Section 3 and Supplementary Tables 2 and 3). These simulations used the same virtual species and combinations of

ecological and shelf parameters as the first two sets of simulations, and extinction magnitudes across simulations were compared using the same framework.

**Sensitivity to definition of shallow marine areas.** To evaluate the potential effect of eustatically driven changes in areas of shallow marine habitat on extinction magnitude across the three climate transitions, we ran additional simulations in which every cell within 0–200 m water depth was considered as potentially inhabitable if thermal conditions were within the virtual species' niche (Supplementary Methods and Supplementary Fig. 4). Shallow marine cells were defined based on coarsely resolved bathymetric layers and virtual species were seeded in each cell (see Supplementary Table 4 for numbers). Species could disperse to cells within the extent of this shallow marine area based on their dispersal ability. The same simulation framework described above was then used to test for the effect of (1) palaeogeography, (2) palaeogeography and sea level drop and (3) palaeogeography, sea level drop and temperature change on extinction magnitude across the three greenhouse-icehouse transitions. Unlike the above simulations, both the position and the extent of suitable habitat changed with sea level drop in these simulations.

## Data availability

Data from simulations are provided as Supplementary Data.

## Code availability

Simulation code is provided as Supplementary Software.

## References

55. Qiao, H., Saupe, E. E., Soberón, J., Peterson, A. T. & Myers, C. E. Impacts of niche breadth and dispersal ability on macroevolutionary patterns. *Am. Nat.* **188**, 149–162 (2016).

56. Rangel, T. F. L., Diniz-Filho, J. A. F. & Colwell, R. K. Species richness and evolutionary niche dynamics: a spatial pattern-oriented simulation experiment. *Am. Nat.* **170**, 602–616 (2007).

57. Saupe, E. E. et al. Non-random latitudinal gradients in range size and niche breadth predicted by spatial patterns of climate. *Glob. Ecol. Biogeogr.* **28**, 928–942 (2019).

58. Blakey, R. C. *Colorado Plateau Geosystems, 1–1* (Deep Time Maps, 2016); https://deptimemaps.com/

59. Scotese, C. *Digital Paleogeographic Map Archive*, CD-ROM (PALEOMAP Project, Scotese, 2001).

60. Harris, J. et al. in *Petroleum Systems Analysis—Case Studies: AAPG Memoir Vol. 114* (eds AbuAli, M. A., Moretti, I. & Nordgård Bolås, H. M.) 37–60 (American Association of Petroleum Geologists, 2017).

## Acknowledgements

We thank C. Myers, A. Townsend Peterson and J. Soberón for help with developing the initial simulation framework, from which simulations were launched. H.Q. was supported by the Natural Science Foundation of China (31772432). E.E.S. was supported by an EAR NSF Postdoctoral Fellowship and Leverhulme Grant #DGR01020. A.P., J.-B.L. and Y.D. thank the CEA/CCRT for providing access to HPC resources of TGCC under allocation 2014-012212 made by GENCI. This is a contribution to IGCP Project-653, 'The Onset of the Great Ordovician Biodiversification Event'. A.F. and D.J.L. acknowledge funding from NERC through NE/K014757/1, NE/1005722/1, NE/1005714/1 and (P.J.V. also) NE/ P013805/1. D.J.L. and P.J.V. acknowledge funding through ERC grant 'The Greenhouse Earth System' (T-GRES, project reference 340923). A.F., A.T.K.-A., D.J.L. and P.J.V. are thankful for the use of the computational facilities of the Advanced Computing Research Centre, University of Bristol (http://www.bris.ac.uk/acrc) (Bluecrystal). S.F. acknowledges funding from the David and Lucile Packard Foundation. This is PBDB publication no. 353. We thank all contributors to the PBDB, with special thanks to the top 10 for our data download: A. Hendy, W. Kiessling, P. Wagner, S. Holland, M. Uhen, B. Kröger, J. Alroy, A. Miller, M. Clapham and J. Sessa.

## Author contributions

E.E.S. and S.F. designed the study, collected the data and wrote the manuscript. A.F., A.P., A.T.K.-A., D.J.L., J.-B.L., P.V. and Y.D. contributed climate modelling data. E.E.S. and H.Q. performed analyses. All authors read and commented on the manuscript.

## Competing interests

The authors declare no competing interests.

## Additional information

**Supplementary information** is available for this paper at https://doi.org/10.1038/s41561-019-0504-6.

**Correspondence and requests for materials** should be addressed to E.E.S.

**Peer review information** Primary Handling Editor: James Super.

**Reprints and permissions information** is available at www.nature.com/reprints.