Reader TTS

ABSTRACT The hyperthermal events of the Cenozoic, including the Paleocene-Eocene Thermal Maximum, provide an opportunity to investigate the potential effects of climate warming on marine

ecosystems. Here, we examine the shallow benthic marine communities preserved in the late Cretaceous to Eocene strata on the Gulf Coastal Plain (United States). In stark contrast to the

ecological shifts following the end-Cretaceous mass extinction, our data show that the early Cenozoic hyperthermals did not have a long-term impact on the generic diversity nor composition

of the Gulf Coastal Plain molluscan communities. We propose that these communities were resilient to climate change because molluscs are better adapted to high temperatures than other taxa,

as demonstrated by their physiology and evolutionary history. In terms of resilience, these communities differ from other shallow-water carbonate ecosystems, such as reef communities, which

record significant changes during the early Cenozoic hyperthermals. These data highlight the strikingly different responses of community types, i.e., the almost imperceptible response of

molluscs versus the marked turnover of foraminifera and reef faunas. The impact on molluscan communities may have been low because detrimental conditions did not devastate the entire Gulf

Coastal Plain, allowing molluscs to rapidly recolonise vacated areas once harsh environmental conditions ameliorated. SIMILAR CONTENT BEING VIEWED BY OTHERS CARBONATE COMPENSATION DEPTH

DRIVES ABYSSAL BIOGEOGRAPHY IN THE NORTHEAST PACIFIC Article Open access 24 July 2023 SPATIAL VARIATIONS IN THE STABLE ISOTOPE COMPOSITION OF THE BENTHIC ALGAE, _HALIMEDA TUNA_, AND

IMPLICATIONS FOR PALEOTHERMOMETRY Article Open access 01 October 2020 MARKED SPATIAL HETEROGENEITY OF MACRO-BENTHIC COMMUNITIES ALONG A SHALLOW-MESOPHOTIC DEPTH GRADIENT IN REUNION ISLAND

Article Open access 30 December 2024 INTRODUCTION Human activities are drastically changing conditions in coastal marine ecosystems by polluting, destroying habitats, overexploiting

resources, enabling invasive species, and driving climate warming. Increased greenhouse gas emissions associated with these activities harm marine communities by expanding hypoxic dead

zones, increasing ocean acidity, and causing thermal stress1,2. In order to understand how communities will respond to climate-related stressors, we look to potential deep-time analogues and

the community shifts recorded by fossils. The Eocene witnessed two separate long-term warming trends of ~6 °C culminating in the late Ypresian and Bartonian, known as the Early Eocene

Climatic Optimum (EECO) and the Middle Eocene Climatic Optimum (MECO), respectively (Fig. S1). Superimposed on these long-term trends are many short-lived intervals of increased carbon

injection into the atmosphere and increased sea surface temperatures, known as hyperthermals, and of these the Paleocene-Eocene Thermal Maximum (PETM) and the Eocene Thermal Maximum 2 have

the highest magnitude and pace3 (Fig. S1). It has been argued that these hyperthermal events represent the best analogues for projected climatic change, as they were caused by rapid

increases in _p_CO2 and involved various environmental consequences, such as ocean acidification and intensification of the hydrological cycle2,4. In addition, the peak of the MECO saw rapid

warming5 and could be considered another hyperthermal. Similar hyperthermal events of smaller magnitude have also been recorded from the Paleocene, e.g., the latest Danian Event6. The PETM

was the most rapid warming event of the early Cenozoic and had the largest ecological impact on marine ecosystems of any hyperthermal during that time7. In shallow-water carbonates there was

a substantial decline in reef volume8 that manifested as a shift from coral-algae reefs to large foraminiferal carbonate ramps9; in the deep-sea, there was a major extinction of deep-sea

benthic foraminifera10, dwarfing of both benthic foraminifera and ostracods11; a rapid diversification of pteropods12, and poleward shifts in the distribution of planktic foraminifera,

dinoflagellates, and radiolarians13. Supposed causes of PETM ecological changes include ocean deoxygenation, rising temperatures, shoaling of the calcium carbonate compensation depth, and

variations in food supply7,13. Nevertheless, the response of non-reefal shallow marine ecosystems to the PETM remains unclear, and there are few studies that quantitatively investigate

changes in the composition of macrobenthic assemblages. Along the United States (US) Atlantic Coastal Plain (South Carolina to New Jersey), the PETM interval contains few or no adult

molluscs, potentially due to low oxygen conditions and/or ocean acidification14,15. Conversely, the US Gulf Coastal Plain (Texas to Georgia) shows no evidence that the PETM resulted in a

diversity decline or body size decline in molluscan communities16,17. Furthermore, Ivany _et al_.16 do not report a faunal turnover at the family-level, but their data do suggest changes in

the dominant genera and species. To improve our understanding of the impact of the early Cenozoic hyperthermals on shallow marine benthic communities, we quantitatively investigated changes

in their diversity and composition along the Gulf Coastal Plain. This study tests the hypothesis that early Cenozoic hyperthermal events were associated with significant, long-term changes

in community diversity and composition. We compiled a dataset of species abundance and richness from the Late Cretaceous through Eocene (Maastrichtian-Priabonian) interval and quantitatively

assessed the molluscan communities for changes in (i) taxon richness, (ii) taxonomic composition, and (iii) functional composition. Although the faunal record does not allow for the

assessment of the short-term (up to millennial-scale) responses of molluscs to the early Cenozoic hyperthermals, our comprehensive analysis shows that the early Cenozoic hyperthermals did

not significantly impact the evolutionary history of benthic molluscan communities. RESULTS DIVERSITY CHANGES Analyses of the raw and Shareholder Quorum Sampling (SQS) data yield similar

trends in the species, generic, and functional richness of benthic assemblages respectively (Fig. 1). Our analyses provide evidence for two substantial changes in generic richness: a

significant (Kruskall-Wallis test (KW): _p_ < 0.01), temporary drop at the Cretaceous/Paleogene (K/Pg) boundary (richness quickly recovered by the P1b biozone, ~700 kyr later) and high

values of richness in the middle Eocene (Lutetian-Bartonian). The only early Cenozoic hyperthermal with enough samples for statistical inference is the PETM, which does not record a

significant change in taxonomic diversity. The late Ypresian samples, which coincide with the Early Eocene Climatic Optimum (EECO), are less diverse than those of other time bins; however,

some outliers have a greater richness than other time bins (Fig. 1a,b). The turnover rates show that the decline of genus richness at the K/Pg boundary is associated with a high generic

extinction rate, whereas the late Ypresian records a relatively low turnover rate suggesting that the decrease in sample diversity is unlikely to reflect a biotic crisis (Fig. 2). Other than

the K/Pg boundary, high turnover rates associated with both elevated extinction and origination of species and genera occur in the late Danian, coincident with the hypothesised late Danian

hyperthermal. Thirty-one different benthic modes of life were recognised in this study, with no more than 21 functional groups in any one sample and no more than six modes of life following

SQS (Fig. 1c). In general, the different time bins fall within categories of high, medium, and low functional richness. The intervals that show significant declines in functional diversity

belong to the early Danian, late Ypresian, and Priabonian samples (KW: _p_ < 0.01) (Fig. 1c). The low functional richness in the early Danian represents a consequence of the

end-Cretaceous mass extinction and has been recognised in other datasets18. In contrast, low functional richness in the late Ypresian, coincident with the EECO, has not been documented

before. COMPOSITIONAL SHIFTS Network analyses of fossil data illustrate compositional differences among samples of different ages, and specifically, how taxa are distributed among the

samples19 (Fig. 3). We structured the data into unipartite and bipartite networks, each consisting of nodes and links. In the unipartite networks (Fig. 3a,b), each node is a sample; a link

indicates that two samples share one or more taxa; and the links are weighted equal to the Bray-Curtis similarity scores of the connected samples20. The bipartite networks (Fig. 3c,d), in

contrast, contain both sample and taxon nodes where samples are connected by non-weighted links to their taxa, and vice versa, but no two samples (or two taxa) are directly connected to each

other19. We partitioned these networks into modules (clusters of nodes) using ‘community detection algorithms’19,20,21,22, which seek to identify the most densely connected nodes separated

by the regions with fewest links. The unipartite network was partitioned into non-overlapping modules using weighted random walks, so the results reflect relative abundances of taxa.

Conversely, the bipartite network was partitioned into overlapping modules using a non-weighted modularity-based algorithm focused on the presence/absence of taxa20. The graphs of these

networks illustrate that, although the Cretaceous and Paleogene samples differ with regard to both species and genera, the various time bins of the Paleogene are, by and large, only

distinguishable at the species-level. Overall, the algorithms identified numerous clusters at the species level—20 in the unipartite network (Fig. 3a) and 12 in the corresponding bipartite

network (Fig. 3c)—but only two modules in the genus-level networks (Fig. 3b,d). These results show that the species-level networks consist of distinct clusters of various Cretaceous,

Paleocene, and Eocene ages; the boundaries between these modules correspond to changes in communities around the K/Pg boundary, early Danian, Ypresian, and the late Bartonian. Conversely,

each genus-level network has two modules: one composed of Cretaceous and early Danian samples and a second consisting of all remaining Paleocene and Eocene samples (Fig. 3b,d). In terms of

taxonomic composition, the amount of overlap between these modules is negligible (Fig. 3d), as the two overlapping modules in the bipartite network share only 15 genera of 797 total (1.9% of

the data). The similarity of results for the unipartite and bipartite networks, which are based on relative abundance and presence/absence data, respectively, suggest that the differences

between the time intervals primarily reflect changes in the presence/absence of taxa. Several samples of early Danian age are assigned to the first module with the Cretaceous nodes because

they evidently consist of genera that survived the K/Pg transition. Beginning at some point in the Danian, the origination/immigration of genera caused significant change in the generic

composition of marine communities, as demonstrated by the clustering of all remaining Paleogene samples. Although the Paleocene and Eocene differ to a large degree in terms of species, their

samples are indistinguishable in terms of genera. Therefore, network analysis shows that mollusc communities experienced greater compositional change at the species-level than the

genus-level throughout the stratigraphic interval of study and that, unlike the end-Cretaceous event, the early Cenozoic hyperthermals did not have a long-term or significant impact on the

generic composition of the benthic assemblages. Nonmetric multidimensional scaling (nMDS) ordination based on species composition shows an evolutionary trend with the greatest shifts in

species composition at the K/Pg boundary, the Danian, the PETM, and the late Bartonian (Fig. 4a,d). These shifts are also larger than at the genus-level due to the shorter stratigraphic

ranges of species. Genus-level compositional changes are reflected in just three clusters: the Maastrichtian, early Danian, and remaining Paleocene-Eocene samples (Fig. 4b). The close

resemblance of the Bray-Curtis and Kulczynski dissimilarity measures between centroids suggests that the disparities between intervals are driven by changes in the presence/absence of

genera, rather than changes in the relative abundances (Fig. 4e). The change in composition between the Cretaceous and Paleocene occurs in two steps: the first corresponding to the mass

extinction event and the second reflecting the evolution/immigration of new genera during the post-extinction recovery (Fig. 4b,e), shown by the high origination rates in the mid-Danian

(Fig. 2b). Although clustering close together, a permutational multivariate analysis of variance (PERMANOVA) shows that the Paleocene and Eocene samples are compositionally dissimilar (t = 

3.0, _p_ < 0.01), suggesting that the PETM resulted in a change in benthic community composition (Fig. 4e). This dissimilarity value is, however, not much higher than the Danian and

Selandian values (t = 2.7, p < 0.01; Fig. 4e). In addition, the t-values of the Paleocene and Eocene samples are consistent with high similarities, as demonstrated by nMDS (Fig. 4b). The

nMDS analysis of the functional diversity data did not resolve any clusters of time bins (Fig. 4c). The functional dissimilarity between time intervals does show similar trends, albeit

noisier: the highest dissimilarity occurs between the Cretaceous and Danian samples. The similarity between the Kulczynski index compared to the Bray-Curtis index suggests that changes in

the presence/absence of modes of life within a sample are causing the differences, rather than the relative abundance of modes of life. The PERMANOVA results show that for the early Danian,

late Danian, PETM, and Bartonian, the t-values are highest suggesting functional composition shifts at these times. Unlike changes in taxonomic composition, these shifts do not cause

evolutionary changes but temporary changes to a different functional state. Like changes in the taxonomic composition, the early Cenozoic hyperthermals did not cause dissimilarity values

greater than in the Paleocene (Fig. 4f), and thus, hyperthermals did not lead to significant restructuring in the functional composition of benthic communities. The PERMANOVA shows that the

early Danian samples have a higher dissimilarity (t = 4.7, _p_ < 0.01) when compared to the late Maastrichtian and Danian samples, but they still show an overlap. The higher dissimilarity

in the early Danian samples is due to their higher abundance of deep-infaunal, facultatively motile, deposit-feeders. DISCUSSION The sequence stratigraphic framework of the Gulf Coastal

Plain is characterised by a series of sequence boundaries, condensed units, and variable sedimentation rates that would be expected to overprint the evolutionary history of molluscan

communities along the Gulf Coastal Plain23. The expected impacts on the fossil record would, therefore, be that first and last appearances are concentrated along sequence boundaries leading

to pseudo-extinction and origination events, variations in time-averaging with beds that were deposited during low sedimentation rates having an elevated richness, and rapid facies changes

causing punctuated shifts in the assemblage composition24. Despite these expected impacts on the fossil record, the only interval in the early Cenozoic with elevated extinction and

origination rates is the late Danian, which is not associated with a major hiatus25. In contrast, the PETM on the Gulf Coastal Plain does coincide with a sequence boundary, but elevated

turnover rates are not recorded (Fig. 2). Likewise, glauconitic lithologies, which are interpreted to have lower sedimentation rates than other lithologies, and thus increased time

averaging, are only significantly more diverse than mixed carbonate-clastic rocks but not more diverse than carbonate or clastic rock types (see Supplemental Material). Together, this

suggests that at the resolution of our study neither the sedimentology nor the sequence stratigraphic architecture of the Gulf Coastal Plain is overprinting the long-term ecological changes

during the early Cenozoic. Our results indicate that the early Cenozoic hyperthermals did not lead to long-term changes in faunal diversity or composition at the genus-level. In comparison

to genera, species generally have shorter stratigraphic ranges, and therefore exhibit higher rates of turnover throughout the stratigraphic interval of study. Consequently, the samples of

various ages differ in composition, suggesting that the Cenozoic hyperthermals may have impacted benthic mollusc assemblages at the species-level. Regardless, the species-level turnover

rates do not show unusually high values associated with the early Cenozoic hyperthermals. High generic turnover rates are only associated with the late Danian, which could be a consequence

of a late Danian hyperthermal. This result stands in stark contrast to observations of the biotic shifts associated with mass extinctions driven by similar environmental changes (e.g.,

end-Permian, end-Triassic, and the early Toarcian extinction events26). Our analyses show that the end-Cretaceous mass extinction and the appearance of new taxa in the Danian caused the only

significant shift in the generic composition of benthic communities in the Gulf region during this time (Figs. 2–4). Ivany _et al_.16 also showed that the PETM did not lead to lasting

changes in body size, species-level diversity, or the life history of dominant benthic clades. These observations suggest that the environmental changes associated with the early Cenozoic

hyperthermals, in particular the PETM, were not detrimental to shallow marine molluscan faunas, despite evidence for rapidly rising temperatures, ocean acidification, and

deoxygenation2,7,27. Our results highlight a major distinction between different ecosystems and environments during the PETM; while the molluscan communities studied here were largely

unaffected, others have recorded a notable (50%) extinction of deep-sea benthic foraminifera28, a decline in metazoan reef volume and an associated turnover of reef fauna8,29, loss of the

sedimentary mixed layer30, as well as poleward range shifts13 and extinctions of nannoplankton31. The Gulf Coastal Plain may not record short-term composition changes due to the coarse

binning of samples in this study and hiatuses at sequence boundaries, but such potential biases did not diminish the diversity and composition changes associated with the end-Cretaceous mass

extinction. Still, if the early Cenozoic hyperthermals led to an increased extinction intensity or the extinction of dominant genera, they would have also resulted in significant

composition shifts and elevated turnover rates. Because we did not find significant differences in taxonomic and ecological composition among early Cenozoic samples (Figs. 3 and 4), we

conclude that the hyperthermals did not lead to long-term changes in benthic communities. This conclusion does not rule out transient changes in other ecological attributes. An increase in

temperature, for example, may reduce body-size, promote parasitism, and/or alter habitat preferences, as demonstrated by work on modern systems32 and Mesozoic hyperthermal events33.

Nonetheless, if any ecological changes did occur during the hyperthermal events in the Gulf Coastal Plain they must have been reversible and short-lived. A plausible explanation for the

resilience of benthic marine communities to early Cenozoic hyperthermals could be that communities consisted of genera inferred to have warm-water affinities and thus could tolerate high

temperatures34. The faunal assemblages of the Gulf Coastal Plain are dominated by gastropods, bivalves, and scaphopods, which have dominated benthic assemblages during hothouse periods since

the Mesozoic35 and include many taxa with modern representatives that have high-temperature tolerances36. Likewise, an investigation of the latitudinal ranges of the species in the

Paleocene and Eocene observed in this study shows that ~30-35% of the species have warm water affinities and ~45-50% are eurythermal species that were not recorded north of the warm

temperate zone (up to 43 °N) (Fig. 5). These warm water affinities of molluscs may explain why the early Cenozoic hyperthermals caused significant ecological turnover in ecosystems dominated

by other faunal groups with lower thermal limits, such as the replacement of corals by large-foraminifera as reef builders8. Even though the environmental changes in the Gulf Coastal Plain

did not lead to long-term impacts on benthic mollusc communities, other regions that experienced even higher temperatures during the early Cenozoic hyperthermal events may have experienced

significant biotic turnover; for example, PETM temperatures of 40 °C in Tanzania caused a temporary exclusion of planktic organisms37. Previous investigations of molluscan extinctions along

the Gulf Coastal Plain concluded that changes in temperature drove the disappearance of taxa over time34,38,39. These studies show that climatic cooling at the Eocene/Oligocene transition

resulted in a regional extinction, including the loss of endemic taxa that were adapted to the high early Cenozoic temperatures34,38. The different responses of molluscs on the Gulf Coastal

Plain to cooling and warming events highlight the relationships between evolutionary history, ecological structure, and environmental change. The impact of a major environmental state shift

on biodiversity, therefore, may depend on the starting state of the system40. A rapid warming event in a greenhouse world will have less impact than one in an icehouse world, and vice versa.

This interpretation is based on the expectation that species in hothouse worlds would have greater thermal tolerance to high temperature or greater thermal plasticity than species in

icehouse worlds. This phenomenon can be seen in modern populations of the gastropods _Chlorostoma funebralis_, which suffer less heat stress and recover quicker from heat stress in hot than

in cooler settings41,42. Furthermore, given that the species from the Gulf Coastal Plain are mostly warm-temperate (Fig. 5), which tend to have wider thermal niches than tropical species43,

and provided most live close to their thermal optimum and not at their upper thermal limits (a condition that seems to be more common in tropical species) this offers another explanation for

the low impact of hyperthermals. Accordingly, the early Cenozoic hyperthermals do not provide ideal analogues for predicting the consequences of current and projected climate change.

Whereas the Paleocene-Eocene events occurred in a hothouse climate, projected climate warming takes place under icehouse conditions that developed following the Oligocene44. Consequently,

tropical ectotherms are likely more vulnerable to climate change today than they were during the early Cenozoic45. In this context, studies of projected climate change will likely

underestimate long-term consequences for the biosphere if they rely on the early Cenozoic hyperthermals as frames of reference. Furthermore, projected anthropogenic climate warming – like

other events in Earth history such as the end-Permian, and end-Triassic mass extinctions – is not only associated with high magnitude climate warming but also a high rate of temperature

change46. Even though the magnitude and rates of warming for these ancient episodes of severe environmental perturbation are underestimated46, they may be better analogues for assessing the

impacts on present-day marine ecosystems47. Apart from heat stress, the development of hypoxic dead zones is a major driver of ecological change in modern shallow marine ecosystems48, and is

perhaps one of the most significant drivers of extinction during ancient hyperthermal events49. The formation and spread of the Gulf of Mexico dead zone in the recent past has caused the

complete loss of benthic invertebrates from some regions50. In addition, runoff and nutrient-input into the Gulf of Mexico has caused the dead zone to expand and marine diversity to

decline51. Increases in runoff, nutrient-input, phosphorus regeneration, and dead zone development also affected shallow marine environments during the PETM14,52, including those along the

Gulf Coastal Plain53. Such environmental deterioration possibly also occurred during the other early Cenozoic hyperthermals. Although no lasting impact of any such environmental changes is

recorded in the benthic communities studied here, the constraints on temporal resolution in our dataset do not rule out the possibility of expanded, but transient, dead zones during the

hyperthermal events. Even if temporary dead zones occurred, benthic communities could rapidly recolonise habitats without significant turnover. Ocean acidification is another major threat to

modern marine communities2, and geochemical studies provide support for shoaling of the calcite compensation depth (CCD) and surface ocean acidification during the PETM2,54,55. The shoaling

of the CCD is thought to have driven the extinction of deep-sea benthic foraminifera during the PETM28, and may have even reached unexpectedly shallow depths15. Surface ocean acidification

may have produced some detrimental effects on benthic fauna (especially for planktotrophic molluscan larvae), but the lack of composition change along the Gulf Coastal Plain, the

diversification of pteropods12,56, and minimal changes in boron/calcium ratios of planktic foraminifera12 suggest that it was not a significant stressor. It is possible that the habitats

studied here were too shallow to have been affected by the shoaling CCD during the early Cenozoic and/or that surface ocean acidification may not have developed to lethal levels. The fauna

investigated in this study may have also possessed adaptions for coping with ocean acidification. Given that molluscs regulate pH and carbonate chemistry at their calcification sites and

possess organic shell coatings57,58,59, the lack of major changes in composition indicates that those genera that dominated the assemblages withstood any effects of adverse environmental

deterioration that developed during the early Cenozoic hyperthermals. MATERIALS AND METHODS We compiled a database of benthic molluscan fossil abundances from the US Gulf Coastal Plain,

which were assigned to 16 time intervals that correspond to formations that could be correlated lithostratigraphically (see Figs. 1 and S1, S3). Faunal data come from our own collections and

the literature. Our collections comprise 142 samples collected from 129 localities along the Gulf Coastal Plain (Dataset S1). Samples were taken at the bed-level and fossils were collected

exhaustively until no further specimens were found. Except for careful cleaning and the removal of matrix, no other preparation was necessary. In total, 111 240 mollusc specimens were

identified to the genus-level, and 93 267 were identified to the species-level. Quantitative faunal data from the literature comes from sources that reported original specimen

counts3,4,5,6,7,8,9,10,11,12,13, increasing the database to 958 samples and 316 767 mollusc specimens (this data is also available from the Paleobiology Database https://paleobiodb.org, and

related collection numbers are available in Dataset S1). Taxonomic identifications limited to the family-level were excluded. Samples with <50 specimens were excluded from the analysis,

which reduced the dataset to 608 samples. Each genus was assigned to a mode of life within a modified version of the Bambach ecospace model after Sessa _et al_.60, which is a combination of

tiering, motility, and feeding (see Supplemental Material). The stratigraphic framework for the timing of the early Cenozoic hyperthermals is poorly understood for the Gulf Coastal Plain

with only the PETM identified from just below the Tuscahoma/Hatchetigbee Formation boundary53. Correlation of the hyperthermals on the Gulf Coastal Plain is, therefore, based on the

nannoplankton zonation (Fig. S3). In our dataset, 32% of the species are informally described (i.e., sp. or spp.), which means that the species variation between samples is likely to be

underrepresented in this study. The sample of greatest richness contains 372 species from 212 genera. Recognising that taxonomic richness varies with the number of identified specimens, we

normalised the data with respect to sampling effort by applying the Shareholder Quorum Sampling (SQS) method61 (SQS; quota = 0.6) using SQS version 3.261 in the statistical programming

environment R62 to calculate subsampled diversity for each sample. Samples with only a single taxon were excluded when calculating subsampled diversity using the SQS function. For

investigating functional diversity, diversity metrics were calculated on the functional groups that were a combination of tiering, motility, and feeding and the abundances of multiple genera

within a functional group were tallied for each sample. For each investigated time bin the taxa were characterized as either: confined to interval (FL), only bottom boundary crossed (bL),

only top boundary crossed (Ft), or both boundaries crossed (bt) following Foote63. From this data the extinction rate ((_N__bL_ + _N__FL_) / total number of taxa), origination rate ((_N__Ft_

 + _N__FL_)/total number of taxa), and turnover rate ((_N__bL_ + _N__Ft_ + _N__FL_)/total number of taxa) were calculated. For multivariate analyses, absolute abundances were converted to

relative abundances, since sample sizes were variable throughout the dataset. Some samples were strongly dominated by few genera and so the percentage data were square root transformed to

deemphasize the influence of the most abundant genera. Nonmetric multidimensional scaling (nMDS) with the Bray-Curtis similarity matrix was computed to visualise the groupings of samples

based on their relative abundances. The stress criterion was used to evaluate goodness-of-fit for the final nMDS ordination64. Furthermore, because stress values (the difference between the

rank correlation of inter-sample distances and the distances among samples in the ordination space24) >0.3 are considered to signify an unsatisfactory representation of the data, those

samples with individual stress values > 0.3 were removed from the ordinations. For the species-level analysis, 11 samples were removed prior to the nMDS ordination because they had high

stress values (>0.8) and they were distorting the rest of the ordination. A permutational multivariate analysis of variance (PERMANOVA) was carried out to quantify the differences in

composition between time bins. These tests were carried out using the Bray-Curtis coefficient for the relative abundance data, and Kulczynski coefficient for the presence/absence data, which

are considered the most appropriate coefficients for handling relative and presence/absence data, respectively64. The resulting t-values from pairwise tests give an absolute measure of

separation between the groups with greater values indicating a greater difference between community structures65. The significance of the t-statistic was investigated using _p-_values with a

significance level of 0.05. In samples where fewer than 999 permutations could be generated, Monte Carlo _p_-values were employed. PERMANOVA tests were also carried out on a

presence/absence version of the dataset using the Kulczynski measure of dissimilarity to investigate if the difference between groups is due to changes in the relative abundances or the

presence/absence of taxa. Unipartite and bipartite networks were partitioned into modules based on the Bray-Curtis similarities of the samples and incidence of taxa among the samples,

respectively (see Supplemental Material). The unipartite network was partitioned into non-overlapping (mutually exclusive) modules with the walktrap algorithm22 and the statistical

significance of the result was confirmed with randomisation testing, which showed that partitioning of randomly generated networks of corresponding size and degree distribution only rarely

(<1%) found modules with higher modularity scores20; (Supplemental Material). The bipartite network was partitioned into overlapping (non-mutually exclusive) modules with the community

overlap propagation algorithm (COPRA21). To find the best output, COPRA was run 100,000 times on the bipartite network, and the solution with the highest extended modularity score

(calculated from the sample nodes) was recorded20; (Supplemental Material). REFERENCES * Henson, S. A., Beaulieu, C. & Lampitt, R. Observing climate change trends in ocean

biogeochemistry: when and where. _Glob. Change Biol._ 22, 1561–1571 (2016). Article  ADS  Google Scholar  * Hönisch, B. _et al_. The geological record of ocean acidification. _Sci._ 335,

1058–1063 (2012). Article  ADS  CAS  Google Scholar  * Littler, K., Röhl, U., Westerhold, T. & Zachos, J. C. A high-resolution benthic stable-isotope record for the South Atlantic:

Implications for orbital-scale changes in Late Paleocene–Early Eocene climate and carbon cycling. _Earth Planet. Sci. Lett._ 401, 18–30 (2014). Article  ADS  CAS  Google Scholar  * Burke, K.

D. _et al_. Pliocene and Eocene provide best analogs for nearfuture climates. _P.N.A.S_ 115, 13288–13293 (2018). Article  CAS  PubMed  PubMed Central  Google Scholar  * Edgar, K. M. _et

al_. Symbiont ‘bleaching’ in planktic foraminifera during the Middle Eocene Climatic Optimum. _Geol._ 41, 15–18 (2013). Article  ADS  MathSciNet  Google Scholar  * Schulte, P. _et al_. Black

shale formation during the Latest Danian Event and the Paleocene-Eocene Thermal Maximum in central Egypt: Two of a kind? _Palaeogeography, Palaeoclimatology, Palaeoecology_ 371, 9–25

(2013). Article  ADS  Google Scholar  * Norris, R. D., Turner, S. K., Hull, P. M. & Ridgwell, A. Marine ecosystem responses to Cenozoic global change. _Sci._ 341, 492–498 (2013). Article

  ADS  CAS  Google Scholar  * Kiessling, W. & Simpson, C. On the potential for ocean acidification to be a general cause of ancient reef crises. _Glob. Change Biol._ 17, 56–67 (2011).

Article  ADS  Google Scholar  * Scheibner, C. & Speijer, R. P. Decline of coral reefs during late Paleocene to early Eocene global warming. _eEarth Discuss._ 2, 133–150 (2007). Article 

ADS  Google Scholar  * Thomas, E. & Shackleton, N. J. The Paleocene-Eocene benthic foraminiferal extinction and stable isotope anomalies. _Geological Society, London, Special

Publications_ 101, 401–441 (1996). Article  ADS  Google Scholar  * Yamaguchi, T., Norris, R. D. & Bornemann, A. Dwarfing of ostracodes during the Paleocene–Eocene Thermal Maximum at DSDP

Site 401 (Bay of Biscay, North Atlantic) and its implication for changes in organic carbon cycle in deep-sea benthic ecosystem. _Palaeogeography, Palaeoclimatology, Palaeoecology_ 346–347,

130–144 (2012). Article  ADS  Google Scholar  * Janssen, A. W., Sessa, J. A. & Thomas, E. Pteropoda (Mollusca, Gastropoda, Thecosomata) from the Paleocene-Eocene Thermal Maximum (United

States Atlantic Coastal Plain). _Palaeontologica Electronica_ 19.3.47A, 1–26, https://doi.org/10.26879/689 (2016). Article  Google Scholar  * Speijer, R. P., Scheibner, C., Stassen, P. &

Morsi, A.-M. M. Response of marine ecosystems to deep-time global warming: a synthesis of biotic patterns across the Paleocene-Eocene thermal maximum (PETM). _Austrian J. Earth Sci._ 105,

6–16 (2012). Google Scholar  * Self-Trail, J. M. _et al_. Shallow marine response to global climate change during the Paleocene-Eocene Thermal Maximum, Salisbury Embayment, USA.

_Paleoceanography_ 32, 710–728 (2017). Article  ADS  Google Scholar  * Bralower, T. J. _et al_. Evidence for shelf acidification during the onset of the Paleocene-Eocene Thermal Maximum.

_Paleoceanography Paleoclimatology_ 33, 1408–1426 (2018). Article  Google Scholar  * Ivany, L. C. _et al_. Little lasting impact of the Paleocene-Eocene Thermal Maximum on shallow marine

molluscan faunas. _Sci. Adv._ 4, eaat5528, https://doi.org/10.1126/sciadv.aat5528 (2018). Article  ADS  CAS  PubMed  PubMed Central  Google Scholar  * Ward, L. W. Stratigraphy and

characteristic mollusks of the Pamunkey Group (Lower Tertiary) and the Old Church Formation of the Chesapeake Group - Virginia Coastal Plain. _U.S. Geol. Surv. Professional Pap._ 1346, 1–78

(1985). Google Scholar  * Hansen, T. A., Upshaw, B. III, Kauffman, E. G. & Gose, W. Patterns of molluscan extinction and recovery across the Cretaceous-Tertiary boundary in east Texas;

report on new outcrops. _Cretac. Res._ 14, 685–706 (1993). Article  Google Scholar  * Muscente, A. D. _et al_. Quantifying ecological impacts of mass extinctions with network analysis of

fossil communities. _P.N.A.S_ 115, 5217–5222 (2018). Article  CAS  PubMed  PubMed Central  Google Scholar  * Muscente, A. D. _et al_. Ediacaran biozones identified with network analysis

provide evidence for pulsed extinctions of early complex life. _Nat. Commun._ 10, 911, https://doi.org/10.1038/s41467-019-08837-3 (2019). Article  ADS  CAS  PubMed  PubMed Central  Google

Scholar  * Gregory, S. Finding overlapping communities in networks by label propagation. _New Journal of Physics_ 12 (2010). * Pons, P. & Latapy, M. In _Computer and Information Sciences

- ISCIS 2005_ (eds. Yolum, P., Güngör, T., Gürgen, F. & Özturan, C.) (Springer-Verlag, 2005). * Dockery, D. T. Punctuated succession of Paleogene mollusks in the northern Gulf Coastal

Plain. _Palaios_ 1, 582–589 (1986). Article  Google Scholar  * Patzkowsky, M. E. & Holland, S. M. _Stratigraphic Paleobiology: Understanding the distribution of fossil taxa in time and

space_. (University of Chicago Press, 2012). * Dockery, D. T. In _Late Paleocene-early Eocene climatic and biotic events in the marine and terrestrial records_ 296–322 (1989). * Alroy, J.

Accurate and precise estimates of origination and extinction rates. _Paleobiology_ 40, 374–397 (2014). Article  Google Scholar  * Zhou, X., Thomas, E., Rickaby, R. E. M., Winguth, A. M. E.

& Lu, Z. I/Ca evidence for upper ocean deoxygenation during the PETM. _Paleoceanography_ 29, 964–975 (2014). Article  ADS  Google Scholar  * Thomas, E. Cenozoic mass extinctions in the

deep sea: What perturbs the largest habitat on Earth? _GSA Spec. Pap._ 424, 1–23 (2007). Google Scholar  * Zamagni, J., Mutti, M. & Košir, A. Evolution of shallow benthic communities

during the Late Paleocene–earliest Eocene transition in the Northern. Tethys (SW Slovenia). _Facies_ 54, 25–43 (2007). Article  Google Scholar  * Rodríguez-Tovar, F. J., Uchman, A., Alegret,

L. & Molina, E. Impact of the Paleocene–Eocene Thermal Maximum on the macrobenthic community: Ichnological record from the Zumaia section, northern Spain. _Mar. Geol._ 282, 178–187

(2011). Article  ADS  CAS  Google Scholar  * Gibbs, S. J., Bown, P. R., Sessa, J. A., Bralower, T. J. & Wilson, P. A. Nannoplankton extinction and origination across the Paleocene-Eocene

Thermal Maximum. _Sci._ 314, 1770–1773 (2006). Article  ADS  CAS  Google Scholar  * Pörtner, H. Climate change and temperature-dependent biogeography: oxygen limitation of thermal tolerance

in animals. _Naturwissenschaften_ 88, 137–146 (2001). Article  ADS  PubMed  Google Scholar  * Piazza, V., Duarte, L. V., Renaudie, J. & Aberhan, M. Reductions in body size of benthic

macroinvertebrates as a precursor of the early Toarcian (Early Jurassic) extinction event in the Lusitanian Basin, Portugal. _Paleobiology_, 1–21 (2019). * Hansen, T. A. Extinction of late

Eocene to Oligocene molluscs: relationship to shelf area, temperature changes, and impact events. _PALAIOS_ 2, 69–75 (1987). Article  ADS  Google Scholar  * Foster, W. J. & Sebe, K.

Recovery and diversification of marine communities following the late Permian mass extinction event in the western Palaeotethys. _Glob. Planet. Change_ 155, 165–177 (2017). Article  ADS 

Google Scholar  * Song, H. _et al_. Anoxia/high temperature double whammy during the Permian-Triassic marine crisis and its aftermath. _Sci. Rep._ 4, 4132, https://doi.org/10.1038/srep04132

(2014). Article  CAS  PubMed  PubMed Central  Google Scholar  * Aze, T. _et al_. Extreme warming of tropical waters during the Paleocene–Eocene Thermal Maximum. _Geol._ 42, 739–742 (2014).

Article  ADS  CAS  Google Scholar  * Ivany, L. C., Patterson, W. P. & Lohmann, K. C. Cooler winters as a possible cause of mass extinctions at the Eocene/Oligocene boundary. _Nat._ 407,

887–890 (2000). Article  ADS  CAS  Google Scholar  * Garvie, C. Microgastropod population changes from the early Cretaceous to the Recent in the Gulf Coastal Plain of the USA. _Zoosymposia_

1, 295–308 (2008). Article  Google Scholar  * Barnosky, A. D. _et al_. Approaching a state shift in Earth’s biosphere. _Nat._ 486, 52–58 (2012). Article  ADS  CAS  Google Scholar  * Gleason,

L. U. & Burton, R. S. RNA-seq reveals regional differences in transcriptome response to heat stress in the marine snail _Chlorostoma funebralis_. _Mol. Ecol._ 24, 610–627 (2015). *

Gleason, L. U. & Burton, R. S. Phenotypic evidence for local adaptation to heat stress in the marine snail _Chlorostoma_ (formerly _Tegula_) _funebralis_. _J. Exp. Mar. Biol. Ecol._ 448,

360–366 (2013). * Stuart-Smith, R. D., Edgar, G. J., Barrett, N. S., Kininmonth, S. J. & Bates, A. E. Thermal biases and vulnerability to warming in the world’s marine fauna. _Nat._

528, 88–92 (2015). Article  ADS  CAS  Google Scholar  * Zachos, J. C., Dickens, G. R. & Zeebe, R. E. An early Cenozoic perspective on greenhouse warming and carbon-cycle dynamics. _Nat._

451, 279–283 (2008). Article  ADS  CAS  Google Scholar  * Dee, S. G., Torres, M. A., Martindale, R. C., Weiss, A. M. & DeLong, K. L. The future of reef ecosystems in the Gulf of Mexico:

insights from coupled climate model simulations and ancient hot-house reefs. _Frontiers in Marine Science_ 6, https://doi.org/10.3389/fmars.2019.00691 (2019). * Kemp, D. B., Eichenseer, K.

& Kiessling, W. Maximum rates of climate change are systematically underestimated in the geological record. _Nat. Commun._ 6, 8890, https://doi.org/10.1038/ncomms9890 (2015). Article 

ADS  CAS  PubMed  PubMed Central  Google Scholar  * Blois, J. L., Zarnetske, P. L., Fitzpatrick, M. C. & Finnegan, S. Climate change and the past, present, and future of biotic

Interactions. _Sci._ 341, 499–504 (2013). Article  ADS  CAS  Google Scholar  * Diaz, R. J. & Rosenburg, R. Spreading Dead Zones and consequences for marine ecosystems. _Sci._ 321,

926–929 (2008). Article  ADS  CAS  Google Scholar  * Harnik, P. G. _et al_. Extinctions in ancient and modern seas. _Trends Ecol. Evol._ 27, 608–617 (2012). Article  PubMed  Google Scholar 

* Rabalais, N. N., Turner, R. E. & Wiseman, W. J. Gulf of Mexico Hypoxia, AKA “The Dead Zone”. _Annu. Rev. Ecol. Syst._ 33, 235–263 (2002). Article  Google Scholar  * NOAA. _Gulf of

Mexico ‘dead zone’ is the largest ever measured_, https://www.noaa.gov/media-release/gulf-of-mexico-dead-zone-is-largest-ever-measured (2017). * Dickson, A. J. _et al_. The spread of marine

anoxia on the northern Tethys margin during the Paleocene-Eocene Thermal Maximum. _Paleoceanography_ 29, 471–488 (2014). Article  ADS  Google Scholar  * Sluijs, A. _et al_. Warming, euxinia

and sea level rise during the Paleocene–Eocene Thermal Maximum on the Gulf Coastal Plain: implications for ocean oxygenation and nutrient cycling. _Clim. Past._ 10, 1421–1439 (2014). Article

  Google Scholar  * Zachos, J. C. _et al_. Rapid acidifcation of the ocean during the Paleocene-Eocene Thermal maximum. _Nat._ 308, 1611–1615 (2005). CAS  Google Scholar  * Penman, D. E.,

Hönisch, B., Zeebe, R. E., Thomas, E. & Zachos, J. C. Rapid and sustained surface ocean acidification during the Paleocene-Eocene Thermal Maximum. _Paleoceanography_ 29, 357–369 (2014).

Article  ADS  Google Scholar  * Byrne, M. In _Oceanography and Marine Biology: An Annual Review_ Vol. 49 (eds R.N Gibsen, R.J.A Atkinson, & J.D.M Gordon) 42 (CRC Press, 2011). * Ries, J.

B. A physicochemical framework for interpreting the biological calcification response to CO2-induced ocean acidification. _Geochimica et. Cosmochimica Acta_ 75, 4053–4064 (2011). Article 

ADS  CAS  Google Scholar  * Ries, J. B., Cohen, A. L. & McCorkle, D. C. Marine calcifiers exhibit mixed responses to CO2-induced ocean acidification. _Geol._ 37, 1131–1134 (2009).

Article  ADS  CAS  Google Scholar  * Lebrato, M. _et al_. Benthic marine calcifiers coexist with CaCO3-undersaturated seawater worldwide. _Glob. Biogeochemical Cycles_ 30, 1038–1053 (2016).

Article  ADS  CAS  Google Scholar  * Sessa, J. A., Bralower, T. J., Patzkowsky, M. E., Handley, J. C. & Ivany, L.C. Environmental and biological controls on the diversity and ecology of

Late Cretaceous through early Paleogene marine ecosystems in the U.S. Gulf Coastal Plain. _Paleobiology_ 38, 218–239 (2012). * Alroy, J. In _Quantitative Methods in Paleobiology_ Vol. 16

(eds John Alroy & Gene Hunt) 55–80 (The Paleontological Society Papers, 2010). * R: A language and environment for statistical computing (R Foundation for Statistical Computing, Vienna,

Austria, 2019). * Foote, M. Origination and extinction components of taxonomic diversity: general problems. _Paleobiology_ 26, 74–102 (2000). Article  Google Scholar  * Clarke, K. R. &

Warwick, R. M. _Change in marine communities: An approach to statistical analysis and interpretation (2nd Edition)_. (PRIMER-E, 2002). * Anderson, M. J. & Walsh, D. C. I. Permanova,

Anosim, and the Mantel test in the face of heterogeneous dispersions: What null hypothesis are you testing? _Ecol. Monogr._ 83, 557–574 (2013). Article  Google Scholar  Download references

ACKNOWLEDGEMENTS We are grateful to the Non-vertebrate Paleontology Laboratory staff at the University of Texas at Austin (Jackson School Museum of Earth History), in particular Liath

Appleton, Angie Thompson, Lisa Boucher, and the late Ann Molineux; without their help this study would not have been possible. We would like to thank Chris Lowery for discussions about the

early Cenozoic succession along the Gulf Coastal Plain as well as Matthew Clapham and two anonymous reviewers for their comments. This project was funded by Geo.X grant SO_087_GeoX and is

associated with the DFG Research Unit TERSANE (FOR 2332: Temperature-related stressors as a unifying principle in ancient extinctions). A.D. Muscente thanks NSF EAR grant #1660005 for

postdoctoral funding. AUTHOR INFORMATION AUTHORS AND AFFILIATIONS * Museum für Naturkunde, Leibniz Institute for Evolution and Biodiversity, Invalidenstraße 43, Berlin, 10115, Germany

William J. Foster & Martin Aberhan * University of Potsdam, Institute for Geosciences, Karl-Liebknecht Straße 24-25, Potsdam-Golm, 14476, Germany William J. Foster * University College

Dublin, School of Earth Sciences, Belfield, Dublin, 4, Ireland William J. Foster & John W. Counts * Non-Vertebrate Paleontology Laboratory, Texas Natural Science Center, The University

of Texas at Austin, 10100 Burnet Road, Austin, Texas, 78758, USA Christopher L. Garvie * The University of Texas at Austin, Department of Geological Sciences, 2275 Speedway, Austin, Texas,

78712, USA Anna M. Weiss, A. D. Muscente & Rowan C. Martindale * Cornell College, Department of Geology, Mount Vernon, Iowa, 600 First Street SW, 52314, USA A. D. Muscente Authors *

William J. Foster View author publications You can also search for this author inPubMed Google Scholar * Christopher L. Garvie View author publications You can also search for this author

inPubMed Google Scholar * Anna M. Weiss View author publications You can also search for this author inPubMed Google Scholar * A. D. Muscente View author publications You can also search for

this author inPubMed Google Scholar * Martin Aberhan View author publications You can also search for this author inPubMed Google Scholar * John W. Counts View author publications You can

also search for this author inPubMed Google Scholar * Rowan C. Martindale View author publications You can also search for this author inPubMed Google Scholar CONTRIBUTIONS W.J.F. and R.C.M.

conceived the study. C.G. collected the data, with modifications by W.J.F. and A.M.W. W.J.F., A.D.M. and J.W.C. analysed the data. W.J.F., R.C.M., M.A., C.G., A.D.M., J.W.C. and A.M.W.

contributed to writing the manuscript. CORRESPONDING AUTHOR Correspondence to William J. Foster. ETHICS DECLARATIONS COMPETING INTERESTS The authors declare no competing interests.

ADDITIONAL INFORMATION PUBLISHER’S NOTE Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations. SUPPLEMENTARY INFORMATION

SUPPLEMENTAL DATA. SUPPLEMENTARY MATERIAL. RIGHTS AND PERMISSIONS OPEN ACCESS This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use,

sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative

Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons license, unless indicated

otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds

the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/. Reprints and

permissions ABOUT THIS ARTICLE CITE THIS ARTICLE Foster, W.J., Garvie, C.L., Weiss, A.M. _et al._ Resilience of marine invertebrate communities during the early Cenozoic hyperthermals. _Sci

Rep_ 10, 2176 (2020). https://doi.org/10.1038/s41598-020-58986-5 Download citation * Received: 07 June 2019 * Accepted: 16 January 2020 * Published: 07 February 2020 * DOI:

https://doi.org/10.1038/s41598-020-58986-5 SHARE THIS ARTICLE Anyone you share the following link with will be able to read this content: Get shareable link Sorry, a shareable link is not

currently available for this article. Copy to clipboard Provided by the Springer Nature SharedIt content-sharing initiative