Articles | Volume 18, issue 9
Research article
01 Sep 2022
Research article |  | 01 Sep 2022

Clumped-isotope-derived climate trends leading up to the end-Cretaceous mass extinction in northwestern Europe

Heidi E. O'Hora, Sierra V. Petersen, Johan Vellekoop, Matthew M. Jones, and Serena R. Scholz

Paleotemperature reconstructions of the end-Cretaceous interval document local and global climate trends, some driven by greenhouse gas emissions from Deccan Traps volcanism and associated feedbacks. Here, we present a new clumped-isotope-based paleotemperature record derived from fossil bivalves from the Maastrichtian type region in southeastern Netherlands and northeastern Belgium. Clumped isotope data document a mean temperature of 20.4±3.8C, consistent with other Maastrichtian temperature estimates, and an average seawater δ18O value of 0.2±0.8 ‰ VSMOW for the region during the latest Cretaceous (67.1–66.0 Ma). A notable temperature increase at ∼66.4 Ma is interpreted to be a regional manifestation of the globally defined Late Maastrichtian Warming Event, linking Deccan Traps volcanic CO2 emissions to climate change in the Maastricht region. Fluctuating seawater δ18O values coinciding with temperature changes suggest alternating influences of warm, salty southern-sourced waters and cooler, fresher northern-sourced waters from the Arctic Ocean. This new paleotemperature record contributes to the understanding of regional and global climate response to large-scale volcanism and ocean circulation changes leading up to a catastrophic mass extinction.

1 Introduction

During the late Maastrichtian, greenhouse gas emissions from the emplacement of the vast Deccan Traps large igneous province (LIP) on the Indian subcontinent resulted in global warming (see Hull et al., 2020, and references therein). This warming event (termed the Late Maastrichtian Warming Event or LMWE; Woelders et al., 2018) has been observed in multiple locations and has been dated to approximately coincide with the onset of major Deccan volcanism (66.413±0.067 Ma; Sprain et al., 2019). Anomalous mercury concentrations in sediments (e.g., Font et al., 2016; Sial et al., 2016; Percival et al., 2018; Zhao et al., 2021) and fossil shells (Meyer et al., 2019) link the warming to a volcanic source, while modeling of the climate impacts of hypothesized volcanic CO2 and Hg emissions results in warming consistent with observations (Tobin et al., 2017; Fendley et al., 2019; Nava et al., 2021).

Prior to the LMWE, the global climate was in a “cool greenhouse” state (Scotese, 2021), with climate models (Miller et al., 2005; Ladant and Donnadieu, 2016), sea-ice-indicating dinoflagellate cysts (Bowman et al., 2013), and coastal Antarctic temperatures near the freezing point (Petersen et al., 2016a), indicating the possible existence of land ice on the Antarctic continent. The LMWE may have caused such an ice cap to melt (Bowman et al., 2013; Petersen et al., 2016a) and, subsequently, global sea level to rise. This is potentially visible in the global sea level curve of Haq (2014), but the resolution of this sea level record is too coarse to confidently link cycles to the short-term LMWE.

After the Late Maastrichtian Warming Event, the climate cooled again gradually, potentially due to silicate weathering feedbacks (Petersen et al., 2016a; Tobin et al., 2017), before being thrown into upheaval by the arrival of the Chicxulub meteorite impact, defining the Cretaceous–Paleogene (K–Pg) boundary (Molina et al., 2006). Ejecta from the impact covered the planet and resulted in a short-lived “impact winter” event (Vellekoop et al., 2014, 2016). Taken together, Deccan Traps volcanism (and associated feedbacks) and the Chicxulub bolide impact create a complex pattern of paleoenvironmental disturbances through the end Cretaceous, ultimately leading to the rapid extinction of ∼70 % of species on Earth (Schulte et al., 2010).

Previous studies have reconstructed Maastrichtian marine temperatures at sites around the globe using several different proxies, including δ18O of foraminifera, mollusks, tooth enamel, and bulk carbonate material (e.g., Pucéat et al., 2007; Tobin et al., 2012; Hull et al., 2020), as well as TEX86 measurements (e.g., Vellekoop et al., 2016; Woelders et al., 2017), Mg/Ca ratios (e.g., Woelders et al., 2018), and clumped isotope thermometry (Dennis et al., 2013; Tobin et al., 2014; Petersen et al., 2016a, b, Meyer et al., 2018, 2019; Zhang et al., 2018; Tagliavento et al., 2019; Table S1 in the Supplement). Taken together, these records clearly show a late Maastrichtian warming (summarized by Hull et al., 2020), although the magnitude of warming differs from site to site. At any given site, ocean temperatures across the LMWE may be influenced by the primary CO2-driven global warming, which may have manifested in amplified warming at high latitudes similar to modern Arctic amplification (Petersen et al., 2016a), and also by local changes in ocean circulation, sea level, upwelling, and/or runoff (changes which may themselves be driven indirectly by global climate change, e.g., sea level rise due to melting of the Antarctic ice caps during the LMWE). Therefore, it is important to consider the local paleogeographic context when interpreting paleoclimate trends during this tumultuous period.

Here, we produce a new climate record for a site in northwestern Europe using the clumped and stable isotopic composition of bivalve fossils to better constrain the magnitude of climate change and possible ocean circulation fluctuations during this key geologic transition interval. We focus on the Maastrichtian type area (Jagt and Jagt-Yazykova, 2012) in the Belgium–Netherlands border region near Maastricht, Netherlands. The reconstructed paleogeography of this locality – a shallow epicontinental sea surrounded by several low-lying landmasses (Fig. 1a) – makes this site sensitive not only to global climate changes, but also to local variations in sea level and ocean circulation. We interpret our data in terms of previously described and modeled estimates for local, regional, and global climate changes during this time. Our study adds to a growing body of literature documenting the impacts of volcanism and ocean circulation on climate, especially leading up to the K–Pg extinction.

Figure 1Paleogeography and sample localities used in this study. (a) Late Cretaceous paleogeographic reconstruction of Europe and surrounding regions based on Vellekoop et al. (2022) with the study area (yellow dot) and hypothesized cold- and warm-water currents (blue and red arrows) marked. (b) A map of modern national borders showing the position of the study area. (c) A zoomed-in map of the Maastricht region, showing the location of the ENCI quarry, CBR Romontbos quarry, and Geulhemmerberg subterranean galleries relative to each other and to the city center of Maastricht.

2 Locality and sample collection

Shell specimens were collected from three locations in the Liège–Maastricht region of northeastern Belgium and southeastern Netherlands: the Maastrichtian type section at ENCI quarry, the CBR Romontbos quarry, and the K–Pg boundary (KPB) section at the Geulhemmerberg subterranean galleries (GSGs; Jagt and Jagt-Yazykova, 2012; Fig. 1c). The ENCI quarry is located ∼3.3 km south of the city of Maastricht in South Limburg, Netherlands (Fig. 1c). Maastrichtian-aged strata at this location include the Gulpen Formation (∼67.1–78.9 Ma) overlain by the Maastricht Formation (∼66.0–66.8 Ma; Vellekoop et al., 2022.), the latter of which was sampled for this study of the uppermost Maastrichtian stage. At this location, the stratigraphically highest exposure of the Maastricht Formation in the Meerssen Member predates the KPB by less than 100 kyr (Keutgen, 2018). The CBR Romontbos quarry is located ∼4 km to the southwest of ENCI quarry, near Eben Emael, Liège, Belgium (Fig. 1c). The Maastricht Formation is also exposed at this location, and samples from here were correlated to align with the Maastricht Formation at ENCI quarry (Felder, 1975). We also include a specimen from the nearby Geulhemmerberg subterranean galleries, located ∼7 km east-northeast of ENCI quarry (Fig. 1c). The exposed section in the GSGs lies stratigraphically above the ENCI quarry and preserves the uppermost section of the Meerssen Member, including the KPB itself (Felder and Bosch, 1998; Keutgen, 2018). Between the top of the ENCI section and the short KPB interval recorded in the GSGs, a time gap exists of no more than several tens of thousands of years (Vellekoop et al., 2022).

At ENCI quarry and in the surrounding region, the Maastricht Formation contains multiple conspicuous fossiliferous horizons (Felder, 1975; J. J. P. Zijlstra, 1988, unpublished data, 1997), including nine horizons in the Maastricht Formation at ENCI quarry (Fig. 2). Fossil shell specimens were collected from five of these horizons (Lichtenberg, ENCI, Romontbos, Lava, and Laumont), as well as from the top of the Emaël Member, the basal and middle intervals of the Nekum Member, and the basal interval of the Meerssen Member (Fig. 2). Ages for all stratigraphic horizons, lithostratigraphic boundaries, and specimens were determined by updating the Keutgen (2018) age model using known biostratigraphic markers and bulk carbon isotope data (Vellekoop et al., 2022; Table S2). Most horizon and unit ages have an uncertainty of ±0.05 Ma, with the exception of the Lichtenberg Horizon (uncertainty of ±0.2 Ma), Valkenburg Member, and St. Pieter Horizon (both ±0.1 Ma; Table S2). Two samples were collected from the middle Nekum Member at the CBR Romontbos quarry from an interval correlated with samples above the middle Nekum Member from ENCI but below the Kanne Horizon.

Figure 2Composite stratigraphic section through the Maastricht Formation at ENCI quarry, Romontbos quarry, and Geulhemmerberg subterranean galleries. Samples collected in the Romontbos quarry and GSGs were placed on the ENCI quarry stratigraphic framework using the age model from Vellekoop et al. (2022; Table S2). From each collection horizon or zone, the number of specimens collected, the number that passed diagenesis screening, the number of well-preserved specimens successfully analyzed for clumped isotopes, and the total number of specimens attempted for clumped isotopes (including some poorly preserved specimens) are listed as follows (n= passed/collected, N= good clumped/attempted). Out of 66 shells sample, 54 shells passed screening, 33 were attempted, and 27 well-preserved shells were successful for Δ47, plus one sample of matrix material from the Meerssen Member.


The K–Pg boundary is positively identified in the Maastricht region and in particular in the GSGs through several lines of evidence (Smit and Brinkhuis, 1996; Vellekoop et al., 2020). The typical planktonic foraminiferal “disaster” assemblage (Smit and Zachariasse, 1996), the presence of the earliest Paleocene dinocyst marker taxon Senoniasphaera inornata (Brinkhuis and Schiøler, 1996; Herngreen et al., 1998), and 87Sr/86Sr analyses of well-preserved foraminifera from the clay layers of the Geulhemmerberg underground galleries (Vonhof and Smit, 1996) have all demonstrated that the K–Pg boundary occurs at the base of IVf-7 of the Meerssen Member at the hardground known as the Berg en Terblijt Horizon (Fig. 2). Directly above this horizon, a sequence of interbedded layers of shell hash and clays represents a period of unique deposition interpreted during a series of high-energy events (storms or seismites) in the 100–300 years immediately following the Chicxulub meteorite impact (Roep and Smit, 1996; Smit and Brinkhuis, 1996; Vellekoop et al., 2020). The one oyster sample collected from these shell hash layers was inferred to have been alive either within this short interval following the KPB or to be reworked from just underneath the Berg en Terblijt Horizon and was thus assigned an age of 66.0 Ma (Clyde et al., 2016; Sprain et al., 2018) with an uncertainty of ±0.01 Ma (Vellekoop et al., 2022), equivalent to the KPB at the scale of our time series (Table S2).

Fossil shell specimens collected from the Maastricht Formation at the ENCI quarry, the CBR Romontbos quarry, and the GSGs included 66 bivalve specimens representing 11 identified taxa of mostly oysters, with a few scallops and clams (from most to least common: Acutostrea uncinella, Agerostrea ungulata, Pycnodonte vesicularis, Gryphaeostrea canaliculate, Entolium membranaceum, Neithea regularis, Pinna sp., Pseudoptera sp., Ceratostreon sp., and Rastellum sp. and indet. Ostreida). Six of the samples were only identified as some type of bivalve or oyster. The majority (40 of 66) are either Acutostrea uncinella or Agerostrea ungulata. All identified specimens were sessile epifaunal suspension feeders and are assumed to have avoided significant lateral postmortem transport. In these units, all biogenic aragonite has dissolved, leaving behind hollowed molds. All taxa sampled here originally formed calcite shells, which continue to be preserved as calcite today (Fig. 3). Water depths for this location are estimated to be less than 20 m at all times and often shallower (Hart et al., 2016), so water temperatures reconstructed from bottom-dwelling mollusks are reasonably interpreted to approximate mean annual surface air temperatures.

Figure 3SEM images of nine shell specimens representing varying levels of calcite preservation. Scale bars represent 1 µm. Images were adjusted for brightness in Adobe Photoshop 2022. Trace element concentrations (in ppm) and SEM preservation indices quantifying the level and prevalence of dissolution and secondary growth (SG) on a scale of 0–3 are shown below each image for comparison. Samples in (a) through (e) were interpreted. Samples in (f) through (i) were not interpreted further due to too much evidence of secondary growth and/or trace elements exceeding thresholds (see Fig. 4).


One sample of carbonate-rich matrix material from the Meerssen Member was analyzed for comparison to fossil shells. The carbonate in this specimen represents a combination of calcitic microfossils and possibly secondary, diagenetic calcite. As such, comparison between its isotopic and trace element composition and that of fossil shells can help assess the potential for diagenetic alteration of shell samples.

3 Methods

3.1 Sample preparation

After removal of matrix material through deionized water rinses and scrubbing with a toothbrush, fossil shell specimens were sampled by either hand-drilling powder from larger fragments or manually crushing smaller fragments in their entirety in a mortar and pestle. The one sample of matrix material from the Meerssen Member was scraped off a fossil and powdered in a mortar and pestle to achieve roughly uniform grain size. The sample powders were divided for trace element analysis and isotopic analysis. Fragments of some shell specimens were mounted for scanning electron microscopy (SEM) observation to assess their preservation state (Fig. 3).

Taxa-specific clumped isotope vital effects in bivalves are considered uncommon based on studies of modern bivalve taxa (Eagle et al., 2013; Henkes et al., 2013; Came et al., 2014), although one recent study documented apparent vital effects in the juvenile portion only of the oyster Magellana gigas (Huyghe et al., 2022). In shells sampled for this study, we estimate that juvenile shell material made up a volumetrically small fraction of sampled shell areas. Combined with the fact that our selected species differ from the taxa demonstrating vital effects (Huyghe et al., 2022), we interpret well-preserved specimens as faithful recorders of paleoenvironmental temperatures at the time of life.

3.2 Preservation assessment (SEM and trace elements)

The isotopic composition of fossil shells is susceptible to post-depositional alteration via recrystallization, which can skew interpreted formation temperatures (Brand and Morrison, 1987). Elevated levels of trace elements (Mn and Fe especially) in calcitic fossil material may indicate diagenesis and/or recrystallization of the calcite (Land, 1967; Möller and Kubanek, 1976; Brand and Morrison, 1987; de Winter et al., 2018). Additionally, the prevalence of diagenetic calcite overgrowth on a sample can be identified easily under scanning electron microscopy (SEM; Fig. 3), as diagenetic calcite does not show the same microstructures as molluskan calcite. All but one of the fossil shells collected for this study were screened for recrystallization using a combination of SEM imaging and trace element detection methods (Figs. 3 and 4, Table S3).

Figure 4Iron and manganese concentrations of all samples used for diagenesis screening. Samples with Fe concentrations above 2500 ppm or Mn concentrations above 100 ppm (or both) were not considered for further interpretation. These thresholds were chosen based on prior studies (see text), then adjusted to lower, more conservative levels to screen out most samples showing secondary growth while keeping all samples showing good preservation or solely dissolution under SEM (see Fig. 4).


All but one of the 66 shells (as well as one matrix sample) were analyzed in the Michigan Elemental Analysis Laboratory at the University of Michigan for Mn and Fe concentrations (Fig. 3) using a Thermo Scientific Quadrupole-ICP-MS (iCAP Q) equipped with an Elemental Scientific PrepFAST 2 automation system. KED mode (He gas) was utilized for all isotopes. All samples were analyzed in 2 % HNO3 (w/v) (Optima grade), and sensitivity drift was corrected for with standard-sample bracketing methods (Table S4).

A subset of 21 specimens representing most horizons were photographed under SEM using a JEOL JSM-7800FLV field-emission SEM in the Robert B. Mitchell Electron Microbeam Analysis Lab at the University of Michigan Ann Arbor (Fig. 3). SEM images were used to assign samples a score of 0–3 for both prevalence of dissolution (Diss) and secondary growth (SG; Table S3). Samples that scored a 2 or 3 for SG were deemed “bad” and not considered for isotopic interpretation as the addition of diagenetic calcite with a (presumably) different isotopic composition than the base shell would skew environmental interpretations. Samples with some evidence of dissolution (score of 2 or 3) but no secondary growth were deemed “okay” and were still considered further because the material that remained likely still represented inhabited environmental conditions. Samples with scores of 0 or 1 for both SG and dissolution were deemed “good” and were fully considered.

3.3 Stable and clumped isotope analysis

Powders drilled or crushed from all shell specimens and one matrix sample were analyzed for oxygen (δ18O) and carbon (δ13C) stable isotope compositions (Figs. S1–S3 in the Supplement, Table 1). Powders from 39 shells were analyzed using a Thermo-Finnigan Kiel IV automated carbonate device coupled to a MAT 253 gas source dual-inlet mass spectrometer in the University of Michigan Stable Isotope Laboratory. Data were corrected for acid fractionation and normalized to the VPDB (Vienna Pee Dee Belemnite) scale using a single-step transfer function between two anchor points, the NBS-18 and NBS-19 standards. Precision was maintained at less than 0.1 ‰ for both δ18O and δ13C (Tables 1 and S4).

Table 1All isotope results. Sample names were created following a convention of three to four capital letters relating to the horizon, then three to four letters relating to the taxa, and then one capital letter to separate multiple shells of the same taxa from that horizon. The age model is from Vellekoop et al. (2022; see Table S2). The preservation index (P index) column combines the results of trace element analysis (good trace – GT, elevated iron – Fe, elevated manganese – Mn, or both – FeMn; or NA if not measured) and SEM imagery (GOOD, OK, BAD, or NA if not imaged) in the format of TRACE/SEM. Any sample that failed ≥1 criteria was not considered further. Errors for Δ47, temperature, and δ18Osw are external 1 SE (see text).

Download XLSX

The remaining 27 samples plus seven of the above 39 were analyzed for δ18O and δ13C in conjunction with clumped isotope analysis in the University of Michigan Stable and Clumped Isotopes for Paleoclimatology and Paleoceanography (SCIPP) Laboratory (Figs. S1–S3 in the Supplement, Table 1). During clumped isotope analysis, one measurement simultaneously acquires Δ47, which is a function of carbonate formation temperature (or past seawater temperature in the case of fossil shells), δ13C, and δ18O of CO2 gas evolved from CaCO3 (Eiler, 2011). Samples are typically replicated a minimum of three times, and replicates are averaged to calculate a sample mean. Here, after combining all replicates for each sample, individual replicates were deemed to be outliers if they fell more than 4 standard deviations outside the mean of the remaining three to four replicates for a given sample. A total of 108 replicates representing 34 unique samples were analyzed over a span of 7 years using two different machines and procedures (see below). Only samples with at least three acceptable replicates were considered (n=27 shells + one matrix sample prior to preservation screening; Table 1).

Between 2016 and 2020, a first set of fossil shell powders were measured for clumped isotopic composition (Δ47) in the University of Michigan Stable Isotope Lab using a Thermo-Finnigan MAT 253 dual-inlet isotope ratio mass spectrometer and a manual sample preparation device described in detail by Defliese et al. (2015), with updates to Porapak temperature outlined by Petersen et al. (2016c). Briefly, 3.5–5 mg of powder is reacted at 75 C in a common acid bath, and evolved CO2 is cryogenically separated from water and other contaminants through two −80C ethanol–liquid nitrogen traps and one PorapakQ and silver wool trap. Clean CO2 is measured for five acquisitions of 12 cycles at a major beam (m/z 44) voltage of 16 V. Raw voltages of masses 44 to 49 were converted to bulk isotopic values (δ45 to δ49) and clumped isotopic values (Δ47, Δ48, and Δ49) using IUPAC (Brand) parameters with an R code presented in Petersen et al. (2019). Raw Δ47 values were placed in a gas-based absolute reference frame following methods introduced by Dennis et al. (2011) using theoretical equilibrium values for heated (1000 C) and water-equilibrated (25 C) standard gases (0.0266 and 0.9198, respectively) and the 75 C acid fractionation factor of +0.072 ‰ for Δ47 from Petersen et al. (2019). Boundaries between each 1- to 3-week-long correction window were selected to account for drift in equilibrium gas line slopes and to maintain consistent values of in-house carbonate standards through time. δ18O of carbonate (δ18Ocarb) was calculated from δ18O of evolved CO2 using an empirical acid fractionation factor defined using in-house carbonate standards independently calibrated to NBS-18 and NBS-19 (α=1.007888 for calcite at 75 C; Table S4). Long-term standard deviations of δ13C, δ18O, and Δ47 on the MAT253 were consistent across time at approximately 0.1 ‰, 0.3 ‰, and 0.020 ‰, respectively, based on in-house carbonate standards run many times over the measurement interval (Table S4).

In 2021 and 2022, additional powders were analyzed in the University of Michigan SCIPP lab using a Nu Perspective isotope dual-inlet ratio mass spectrometer connected to a NuCarb automated sample preparation device (Jones et al., 2022). This device reacts sample powders in individual vials under vacuum at 70 C by slowly injecting 150 µL of 103 % phosphoric acid. Vials are isolated for the first 5 min of the reaction to reduce excess bubbling, then opened to a variable temperature trap held at -160 for an additional 15 min, continuously collecting evolved CO2. Next, this variable-temperature trap is elevated to −60C, releasing CO2 to pass through a static trap filled with Porapak Q material held at -30. CO2 is frozen on the far size in another variable-temperature trap set at -160 for 800 s. CO2 is warmed up to the gas phase, and the yield is recorded with a transducer. The Nu Perspective and NuCarb can handle both larger (3–6 mg) and smaller (350–450 µg) samples. Following purification, depending on the transducer reading, the automated sequence directs CO2 from a larger sample to equilibrate into the bellows or freezes CO2 from a smaller sample into a cold finger between the bellows and the change-over block. All sample replicates in this study were 3–5 mg in size and were therefore analyzed in “bellows mode”. Previous studies using a Nu Perspective and NuCarb for clumped isotope analysis described “cold-finger mode” (Mackey et al., 2020). In bellows mode, the bellows are compressed until a desired beam strength of 80 nA is achieved on the major (m/z 44) ion beam. Gas is analyzed for four blocks of 20 reference sample cycles, with 20 s of integration on each half-cycle. Bellows are continually adjusted between each cycle to maintain the initial beam strength at all times.

As above, raw beam intensities of masses 44–49 were converted to bulk isotopic values (δ45 to δ49) and clumped isotopic values (Δ47, Δ48, and Δ49) using IUPAC (Brand) parameters (Brand et al., 2010) with an R code presented in Petersen et al. (2019), adjusted slightly to accommodate differing data formats from the Nu instrument. Isotopic values were converted into the Intercarb Carbon Dioxide Equilibrium Scale (I-CDES25) absolute reference frame using Δ47 values for four ETH standards defined by the Intercarb project (Meckler et al., 2014; Bernasconi et al., 2021) and a Δ47 acid fractionation factor of +0.066 ‰ for 70 C from Petersen et al. (2019). First, a single slope was fitted through ETH 1 and ETH 2 (adjusted by 0.0033 ‰) in δ47 vs. Δ47 space, then all four ETH standards were used for the empirical transfer function step. Some heated gases (1000 and 200 C) and equilibrated gases (25 C) were analyzed for cross-comparison but were not used in reference frame calculations. Carbonate δ18O and δ13C were each calculated via a single-step transfer function defined by the four ETH standards using values published in Bernasconi et al. (2018). Although use of ETH standards does not require an acid fractionation factor when applied to an entire dataset, combining data from gas-based and carbonate-based reference frames requires normalization to the same acid digestion temperature. Long-term standard deviations of δ13C, δ18O, and Δ47 on this instrument were consistent across time at better than 0.05 ‰, 0.1 ‰, and 0.018 ‰, respectively (Table S4).

Clumped isotope (Δ47) values for between two and five in-house carbonate standards were tracked to ensure consistency and comparability of isotope values through time and across machines. Stable isotope (δ13C and δ18O) values for in-house carbonate standards were determined relative to NBS-18 and NBS-19 on the same Kiel IV + MAT253 used for bulk sample δ13C and δ18O (Table S4). Δ47 values for these standards were defined based on comparison to gas standards over multiple years (2015 and 2018–2019 for standards CM, OO, and CORS) or comparison to ETH standards over a 6-month measurement session in 2020 (for all five; Table S4). Close agreement in Δ47, δ18O, and δ13C between long-term values and session-specific values on all three instruments confirms that data from different measurement sessions and instruments can be directly compared (Table S4). Comparison of replicates of the same sample powder run on each machine supports this.

Finally, acid-corrected Δ47 values (CDES25) from both measurement periods were converted to temperature using the composite calibration of synthetic carbonates from Petersen et al. (2019) [Temp(C) =(0.0383×106)/(Δ47-0.258)-273.15]. This equation was chosen in favor of the ETH-standard-based calibration of Anderson et al. (2021) due to the use of a gas-based reference frame for the early measurements and the fact that calibration data acquired on the MAT 253 by Defliese et al. (2015) and Winkelstern et al. (2016) are included in the Petersen et al. (2019) composite calibration equation. Use of the Anderson et al. (2021) equation results in temperatures 1.5–2.5 C cooler and δ18Osw values proportionally lower but does not change any of the conclusions in this study. Carbonate δ18O values and Δ47-derived temperatures were used to calculate seawater δ18O (δ18Osw) for each sample via the calcite–water relationship of Kim and O'Neil (1997). All δ18Osw values are reported on the VSMOW scale (Table 1).

Long-term reproducibility (1SD) of Δ47 in carbonate standards (0.020 ‰ for MAT253 and 0.018 ‰ for Nu) was used to calculate a long-term standard error in Δ47 for each sample. For temperature and δ18Osw values, a long-term standard error was calculated using the long-term standard error in Δ47 and the pseudo-linear relationship between internal 1 SE on Δ47 and either temperature or δ18Osw in all samples (for this study: long-term 1 SE in temperature =306.7× (long-term 1 SE in Δ47)+0.1778, and long-term 1 SE in δ18Osw=66.33× (long-term 1 SE in Δ47)+0.0461). The external 1 SE for all three parameters was selected as the larger of either the internal or the long-term standard error. All replicate-level clumped isotope and stable isotope data have been submitted to the permanent data archive EarthChem as part of the ClumpDB database.

4 Results

4.1 Preservation and sample screening

Isotopic signatures may be affected by both bulk chemical alteration (dissolution and recrystallization replacing original shell material and/or overgrowth of secondary material changing the bulk isotopic composition) and solid-state bond reordering (changes in the number of heavy isotope clumps in the crystal lattice, overprinting the original recorded temperature without any mass transfer). Overall, SEM preservation ratings and trace element concentrations were low in most samples (Figs. 3 and 4, Table S3), ruling out significant chemical alteration, and a mild thermal history suggests reordering is not a concern.

In the 21 shells imaged for SEM, 6 displayed evidence of secondary growth (SG =2 or 3) and were not considered further. One sample showed evidence of dissolution but no secondary growth and was deemed “okay” and included in interpretations. A total of 14 samples showed little to no evidence of secondary growth (SG =0 or 1) and were deemed “good” and included. In the 65 shell specimens tested for trace elements, Mn concentration ranges from 17–210 ppm, while Fe concentration ranges from 96–7151 ppm. The matrix sample was measured to have 60 ppm Mn and 664 ppm of Fe. Threshold values of 100 and 2050 ppm were selected for Mn and Fe, respectively, to exclude most samples that failed SEM screening, while including all samples that exhibited only original material under SEM (Figs. 3 and 4). These thresholds are more conservative than those proposed by other studies (Brand and Morrison, 1987; Voigt et al., 2003; Ullmann et al., 2013, de Winter et al., 2018) but, in total, still only result in the flagging of 9 specimens out of the total 65 (Fig. 4, Table S3).

Bond reordering without exchange of mass can occur if samples are heated above 100–150 C (Henkes et al., 2014; Winkelstern and Lohmann, 2016), such as during burial. Heating breaks and reforms bonds between atoms in the carbonate lattice and changes the clumped isotopic composition – and therefore interpreted temperature – without affecting the bulk δ18O or δ13C composition appreciably. Regional geologic evidence suggests a relatively shallow burial history for these strata. In the southern Netherlands, the Roer Valley Graben contains several basin-wide unconformities dated to the Late Cretaceous and middle Paleocene that signify two basin inversions, leaving the Maastrichtian and Danian strata relatively undisturbed (Luijendijk et al., 2011, and references therein). Outstandingly well-preserved mosasaur skeletons and hollow aragonitic molds in the Maastricht Formation further indicate a lack of burial or significant diagenetic fluid flow during this time (Dortangs et al., 2002; Jagt et al., 2016). Overall, these burial conditions are not expected to result in bond reordering in any samples.

Isotopic evidence supports the idea that selected trace element and SEM thresholds are conservative and that overall diagenetic alteration in this region is mild to absent. Despite being excluded based on trace element and/or SEM results, the isotopic values of excluded samples were often not substantially different than other shells that passed diagenesis screening (e.g., only one of four measured for clumped isotopes had a higher-than-normal Δ47-derived temperature; Fig. 5). Additionally, the matrix sample taken from the Meerssen Member showed trace element concentrations below thresholds (Fig. 4) and a temperature of 19.8±3.7C (Fig. 5), very similar to temperatures from fossil shells, indicating that the bulk matrix – which is more likely to be affected by recrystallization or reordering than a dense calcite shell – has also not undergone noticeable isotopic or geochemical alteration. This also explains how three samples failed SEM screening (Fig. 3) but showed low trace element concentrations (Fig. 4). We infer that in these samples, although some secondary growth was seen under SEM (Fig. 3), it must have been volumetrically insignificant in the total powdered sample, not present in the portion of the shell that was powdered (which differed from the portion imaged under SEM), and/or did not have elevated trace element concentrations.

Figure 5Maastrichtian Δ47-derived paleotemperatures for all samples. Samples shown in red (n=4) did not pass screening thresholds and were not included in further interpretation. Error bars represent the larger of either the internal or long-term 1 SE (see text). Onset of the main Deccan Traps eruptions (grey shading) at 66.413±0.067 Ma (Sprain et al., 2019). Member names and horizon abbreviations are labeled across the top of the figure (see Fig. 2). * Schiepersberg Member.


In total, 12 of 65 shell samples failed one or both screening techniques, with 9 failing trace element thresholds, 6 failing SEM criteria, and 3 of the 6 or 9 failing both. The one sample that was not measured was assumed to pass screening based on the majority of samples passing screening. An additional three samples failed stable isotope measurements for analytical reasons, and two more failed to achieve the three good replicates necessary for Δ47 but were still used for δ18O and δ13C. In the six samples that were measured on both the Kiel and one of the clumped isotope machines, clumped-isotope-derived δ18O and δ13C values were used (but only differed by 0.1 ‰–0.2 ‰ from Kiel-derived values). After removing poorly preserved and insufficiently analyzed samples, the dataset includes 51 well-preserved shell samples analyzed for δ18O and δ13C and 23 well-preserved shell samples with acceptable Δ47 analyses (Table 1). All further interpretation excludes discussion of removed samples.

4.2 Stable isotopes

The 51 shell specimens exhibit a range of carbonate δ18O (δ18Ocarb) values between −2.3 ‰ and −0.5 ‰ VPDB, with 2 samples showing higher values of 0.0 ‰ and 0.7 ‰ VPDB (Fig. S1 in the Supplement). The average δ18Ocarb value throughout the section is -1.2±0.5 (1 SD) VPDB. In contrast, the carbonate fraction of bulk matrix in the one sample from the Meerssen Member measured here and in a high-resolution record from the Lichtenberg Horizon to the top of the Nekum Member (Vellekoop et al., 2022) falls anywhere from ∼0 ‰ to ∼1 ‰ higher than molluskan carbonate (Fig. S1 in the Supplement), with a mean value of -0.6±0.4 ‰ (1 SD) VPDB and a range of −1.7 ‰ to 0.0 ‰ VPDB through the section.

Carbonate δ13C (δ13Ccarb) values range from 1.4 ‰ to 3.4 ‰ VPDB, with one high sample at 4.0 ‰ VPDB (Fig. S2). The average δ13Ccarb value is 2.8±0.6 ‰ (1 SD) VPDB. Bulk matrix δ13C is clearly lower than most molluskan carbonate by anywhere from 1 ‰ to 2 ‰ VPDB (Fig. S2). The average bulk matrix δ13C value is 1.7±0.1 ‰ (1 SD) VPDB with a range of 1.3 ‰ to 1.9 ‰ VPDB through the section.

In most cases, the number of specimens representing each taxon is too small and/or the stratigraphic range is too limited to make reliable interspecies comparisons (e.g., Rastellum sp.). Even considering this, the pectinids Entolium membranaceum and Neithea regularis appear to record carbon isotope compositions ∼1 ‰1.5 ‰ more depleted relative to other species, which mainly represent Ostreida, while oxygen isotopic compositions appear similar (Fig. S3). This inter-order difference may represent differing use of metabolic (respired) carbon vs. environmental carbon between these two orders of bivalvia. Between the two most common taxa (Acutostrea uncinella and Agerostrea ungulate – both oysters), there does not appear to be any systematic offset in oxygen or carbon isotopic compositions (Fig. S3).

4.3 Temperature estimates

A total of 4 of the 27 shell specimens successfully measured for Δ47 were excluded from further interpretation due to evidence of diagenesis (described above). The 23 samples that pass diagenetic screening criteria exhibit paleotemperatures ranging from 9.5–27.2 C, with a mean value of 20.4±3.8C (1 SD; Fig. 5). In all instances in which multiple taxa were sampled from the same horizon, measured temperatures are within error of each other, suggesting no interspecies offset due to vital effects. In the absence of vital effects in any taxa, all species are expected to record the same temperature due to their common benthic habitat. The high agreement between taxa from the same horizon seen here is similar to that seen in other deposits of the same age from the Gulf Coastal Plain (Meyer et al., 2018), although further replication to reduce error bar size would be necessary to rule out small vital effects.

In most horizons with multiple samples, estimated temperatures overlap within error for all samples (Fig. 5). The exception is the Lava Horizon (66.4 Ma), which shows a large spread in temperatures from 9.5±3.8C (external 1 SE) to 22.3±6.9C (external 1 SE; Fig. 5). This may represent seasonal aliasing in samples (e.g., powdering a portion of shell representing wintertime growth vs. summertime growth), although these samples were not treated any differently than other horizons where this is not observed. This may also represent long-term time averaging into a single horizon, where individual shells may have lived thousands of years apart. It is not a species-related offset, as the highest and lowest temperatures both occur in Acutostrea uncinella (Fig. 5). A similar argument can be extended to other horizons where a larger spread is seen, such as the lower Meerssen Member (66.1 Ma) and the Romontbos Horizon (66.55 Ma), although samples in these layers are still within error at the external 1 SE level.

A fluctuating temporal trend is observed in which paleotemperatures rise from a low at 66.4 Ma to a peak in the basal Nekum Member at 66.325 Ma, followed by a decrease by 66.275 Ma that persists through the KPB (Fig. 5). The onset of warming coincides with the onset of Deccan volcanism at 66.413±0.067 Ma (Sprain et al., 2019; Fig. 5), within the error of our age model (Table S1 in the Supplement). We therefore identify this warming event as the named Late Maastrichtian Warming Event (LMWE; Woelders et al., 2018; Hull et al., 2020).

4.4 Seawater δ18O estimates

Calculated δ18Osw values ranged from −0.9 ‰ to +2.1 ‰ VSMOW (Fig. S4). The average δ18Osw value is +0.2±0.8 ‰ (1 SD) VSMOW, which is above the ice-free mean ocean water estimate of −1.0 ‰ VSMOW. This heavier average δ18Osw value is consistent with the presence of some continental ice volume (presumed to be located on the Antarctic continent) during the cool greenhouse Maastrichtian interval (Miller et al., 2005; Petersen et al., 2016a; Ladant and Donnadieu, 2016). Just like with temperature, where multiple taxa were sampled within a given horizon δ18Osw values are within error of each other, although there is fairly large spread in the Lava Horizon (Fig. S4).

Seawater δ18O values show a similar temporal trend as that of temperature (Fig. 6c and d). We observe peak δ18Osw values at the base of the Nekum Member (1.4±0.7 ‰ VSMOW, horizon mean and 1 SD of three samples), and the lowest δ18Osw values occur at the Lava Horizon in the Emaël Member (-1.0±1.8 ‰ VSMOW, horizon mean and 1 SD of three samples; Figs. 6d and S4). This relationship is not a function of calculations of temperature or δ18Osw and, instead, falls naturally out of the δ18Ocarbδ18Osw–temperature relationship. Across the entire dataset, both temperature and δ18Osw show weak correlation with δ18Ocarb, with neither as the primary sole driver (Fig. S5).

Figure 6Isotopic compositions and paleotemperatures of fossil shells compared to inferred relative sea level in the Maastricht Formation. Error bars represent 1 SD for δ18O (a) and δ13C (b) and represent external 1 SE for temperature (c) and δ18Osw (d) (see Sect. 3.3). All samples failing diagenetic screening criteria have been removed. Colored lines pass through horizon means. Age model from Vellekoop et al. (2022; Table S3). The relative sea level curve (e) is based on Schiøler et al. (1997), adjusted to the age model used in this study. Hermatypic coral appearance data (e) are from Leloux (1999). Onset of the main Deccan Traps eruptions (grey dotted line) at 66.413±0.067 Ma (Sprain et al., 2019).


5 Discussion

5.1 Further evidence for a warmer Maastrichtian world and a Late Maastrichtian Warming Event

New data presented here indicate that mean annual ocean temperatures in the Maastricht region (paleolatitude ∼40 N; van Hinsbergen et al., 2015) varied from 9.5±3.9 to 27.2±3.7C between 67.1 and 66.0 Ma, with a mean temperature of 20.4±3.8C (1 SD) over this interval. This is much warmer than the modern-day climate of Maastricht, which has a mean annual air temperature of 9.8 C (Maastricht Climate, 2021). Warmer-than-modern temperatures in the Maastrichtian have been seen before in many studies (Table S1 in the Supplement) and are expected due to elevated atmospheric CO2 levels during the Maastrichtian of ∼1000 ppm on average over the study interval (Zhang et al., 2018; Henehan et al., 2019; Hoenisch, 2021). Reconstructed Maastrichtian marine temperatures from sites located between 30 and 50 N show an average temperature of 19.8±5.2C (1 SD, n=23 studies; Table S1 in the Supplement), in close agreement with our site mean. Outside the marine realm, Maastrichtian-aged fossil plants also suggest that terrestrial temperatures at 40 N were warmer, measuring around 15 C (Golovneva, 2000), and that near 49 N, air temperatures varied from 10–18 C during this time interval (Wilf et al., 2003). General circulation models using Late Cretaceous paleogeography and elevated atmospheric pCO2 levels also display a similarly warmer climate globally as a result of CO2-induced warming (Tabor et al., 2016; Ladant et al., 2020). Modeled zonal mean temperatures at a paleolatitude of 40 N are 20.2 C for a  CO2 (560 ppm) Maastrichtian scenario or 21.4 C for a  CO2 (1120 ppm) Maastrichtian scenario (Tabor et al., 2016).

Multiple paleoclimate records document a Late Maastrichtian Warming Event (LMWE) (e.g., Stott and Kennett, 1990; Li and Keller, 1998; Petersen et al., 2016a; Woelders et al., 2017, 2018; Barnet et al., 2019; Gao et al., 2021; Nava et al., 2021) beginning gradually at 66.4 Ma with more rapid warming at 66.3 Ma (Hull et al., 2020; Fig. 7). The timing of this LMWE appears to closely follow the onset of LIP volcanism on the Indian subcontinent based on 40Ar/39Ar ages for the oldest basalts of the Deccan Traps (66.413±0.067 Ma from the Jawhar Formation; Sprain et al., 2019) and shifts in marine osmium isotope ratios towards mantle values beginning ∼400 kyr prior to the KPB (Ravizza and Peucker-Ehrenbrink, 2003; Hull et al., 2020), although others argue for a slightly later onset around 66.311±0.051 Ma (Schoene et al., 2019; Nava et al., 2021). Volcanically derived CO2 is widely considered to be the cause of this warming.

Figure 7Temperature evolution through the Late Maastrichtian Warming Event. Comparison of new temperature data to the compilation of late Maastrichtian temperatures by Hull et al. (2020). Grey and light blue vertical bars indicate onset of Deccan volcanism and the KPB, as in Fig. 5. Vertical error bars on horizon averages represent 1 SD of all temperatures from a given horizon. Horizons with no error had only one sample. Horizontal error bars represent uncertainty in the age model of Vellekoop et al. (2022). The brown dashed line indicates a temperature profile calculated from δ18Ocarb assuming a fixed δ18Osw value of −1 ‰ VSMOW and, similarly to the Hull et al. (2020) composite record, shows peak temperatures delayed relative to the Δ47-based temperature peak.


The Δ47-derived temperature record generated here documents a warming in the Nekum Member beginning from a temperature low in the Lava Horizon at 66.4±0.05 Ma (16.6±5.9C, 1 SD, N=4) and increasing to the basal Nekum at 66.325±0.05 Ma (25.7±1.4C, 1 SD, N=3) over 75 kyr (Fig. 5). The magnitude of this warming is calculated as 9.1 C going from these listed horizons or, more conservatively, 5.8 C comparing the warmest horizon (basal Nekum) to a pre-event baseline averaging the preceding seven horizons (19.9±2.0C, N=14). Although this warming event is only barely statistically distinguishable as differing from the mean of the full record when considering the size and meaning of the 1 SE error bars and the fairly small number of samples representing peak warmth (n=3; Fig. 5), the close agreement between all three points from this horizon and the following additional lines of evidence support the interpretation of this warming as an expression of the LMWE.

Based on a newly updated age model (Vellekoop et al., 2022; Fig. S1 in the Supplement), the warmest Δ47-derived temperatures slightly follow or align with the onset of earliest Deccan eruptions (depending on whether a date of 66.413 or 66.3 Ma is used) and the general age of the globally defined LMWE (Hull et al., 2020; Fig. 7). The dinocyst Palynodinium grallator – a marker for the LMWE at northern midlatitudes (Vellekoop et al., 2018; 2019) – first appears in the Maastricht Formation in the Nekum Member (Schiøler et al., 1997; Vellekoop et al., 2019). Thermophilic hermatypic scleratinian macrofossils also first emerge in high abundances in the Nekum Member (Liebau, 1978; Leloux, 1999; Fig. 6e). Together, these records suggest that the LMWE is preserved in the Nekum Member at ENCI quarry and that the warmer Δ47-derived temperatures found there are likely real.

Interestingly, a δ18Ocarb-based temperature reconstruction assuming constant δ18Osw (red dashed line, Fig. 7) incorrectly shows a slightly later peak in temperature. The true peak, as indicated by the Δ47-derived record, is hidden in δ18Ocarb due to simultaneous increases in temperature and enrichment in δ18Osw (Figs. 5, 6, and S4). The timing of warming in our Δ47-derived record appears to predate that of the composite record of Hull et al. (2020) as well (Fig. 7), although uncertainty in age in both records may not sustain this relationship. If age models hold as shown, the later timing seen in the composite record may reflect the same issue of a hidden LMWE peak, as many of the records included are based on foraminiferal δ18Ocarb and include an assumption of constant δ18Osw. Another Δ47-derived temperature record across this period from Seymour Island, Antarctica (Petersen et al., 2016a), also shows simultaneous increases in temperature and δ18Osw as well as a delayed apparent warming in a δ18Ocarb-based temperature reconstruction, although the absolute timing of the LMWE in that record is apparently later (∼66.24–66.2 Ma). Faulty assumptions of constant sedimentation rate between magnetostratigraphic datums could potentially explain the difference in absolute timing.

Following on the assumed identification of the LMWE in the basal Nekum Member, our temperature record indicates that seawater temperatures in this region returned to baseline levels by the middle Nekum Member within ∼100 kyr (Fig. 5). This is consistent with the timing of silicate weathering feedbacks, as previously suggested (Caldiera and Rampino, 1990; Dessert et al., 2001; Petersen et al., 2016a; Tobin et al., 2017; Hull et al., 2020), further supporting the volcanically derived CO2 as the source of the initial warming event. Others have suggested that the emplacement of large igneous provinces can actually be a long-term CO2 sink and lead to longer-term global climate cooling (Schaller et al., 2012; Johansson et al., 2018). This is possibly visible in the slight cooling between the middle Nekum Member and Meerssen Member (Fig. 5), although the uncertainty in the Meerssen Member is large.

5.2 Ocean circulation, sea level, and/or precipitation changes

Paleogeographic reconstructions (Engelke et al., 2017; Scotese, 2021; Fig. 1a) and the presence of seagrass fossils (Hart et al., 2016; van der Ham et al., 2017) indicate that during the Campanian–Maastrichtian, the study site was located on a large, shallowly submerged carbonate platform. Over the Late Cretaceous, eustatic sea level fluctuated multiple times, including one sea level transgression (KMa5) in the latest Maastrichtian at 66.8 Ma (Haq, 2014). Sea level based on the relative abundance of the palynomorph Paralecaniella in ENCI quarry sediments – the same section as our fossils – indicates local sea level fluctuations that differ from the global eustatic curve (Schiøler et al., 1997; Fig. 6e). Since high relative abundances of this acritarch taxon are suggested to reflect marginal-marine to restricted-marine conditions under high hydrodynamic conditions (Brinkhuis and Schiøler, 1996; Schiøler et al., 1997), observed fluctuations are interpreted as local sea level changes. Alternatively, environmental factors other than sea level might have driven the observed high abundances of Paralecaniella in the type Maastrichtian, such as biological responses to local changes in hydrodynamic conditions and salinity, which would also be reflected in the δ18Osw record presented here. Nevertheless, seagrass and foraminifera indicate that water depths remained less than 20 m at all times and often shallower (Hart et al., 2016).

Although this local sea level record derived from palynological abundances may have some issues at a fine scale related to minimal consideration of sedimentology and geochemistry (Vellekoop et al., 2022), we see a notable correlation between it and many of our isotopic parameters (Fig. 6). Peak δ18Osw values, high temperatures, and, to a lesser degree, higher δ13Ccarb values correlate with periods of maximum sea level, while lower δ18Osw values, cooler temperatures, and lower δ13Ccarb values correspond to periods of relatively lower sea level (Fig. 6). In particular, sea level is interpreted to rise from the Lava Horizon up to the Laumont Horizon, peaking around the base of the Nekum Member and then declining again by the Kanne Horizon (Fig. 6e). This aligns with the temperature and δ18Osw peak in the basal Nekum highlighted above (Fig. 6c and d), coinciding with the onset of Deccan volcanism. A second interval of peak sea level occurs during the Gronsveld Member (Fig. 6e), which also displays temperatures and δ18Osw values slightly higher than preceding and following horizons (Fig. 6c and d). δ13Ccarb values also rise into the Gronsveld Member and again rise from the Lava to Laumont horizons (Fig. 6b).

We propose that the observed correlation between environmental conditions and reconstructed local sea level is the result of fluctuating ocean circulation patterns. Paleogeographic reconstructions suggest the presence of a seaway connecting the Arctic Ocean to this study site during the late Maastrichtian (Engelke et al., 2017; Scotese, 2021; Fig. 1a). Modeling studies show that increased isolation of the Arctic made it fresher during the Maastrichtian (Ladant et al., 2020), suggesting that this boreal water mass was both cooler in temperature and depleted in δ18Osw. As a result, the Maastrichtian type region would have been under the influence of relatively cool, depleted water masses. Our study site was situated on a shallow carbonate platform surrounded by small landmasses, the exact size and position of which are difficult to define (Engelke et al., 2017; Scotese, 2021; Fig. 1a). During periods of low relative sea level, the landmasses of the Rhenohercynian Zone (Rhenish Massif, London–Brabant Massif) would have blocked southern currents originating in warmer regions, which would likely be saltier and have higher δ18Osw values due to evaporative enrichment (e.g., the semi-restricted Western Tethys; Fig. 1a). However, as large parts of these landmasses were relatively flat and low-lying during the Upper Cretaceous (Vandenberghe et al., 2014), relatively minor increases in sea level could potentially breach sills or island barriers, changing both the position of ocean currents and the dominant water mass in a given area. Indeed, lithostratigraphic studies in the region have shown that during various phases in the Campanian–Maastrichtian, the structural highs of the nearby London–Brabant Massif were fully submerged (e.g., Dusar and Lagrou, 2007). As a result, during sea level lows, the Maastrichtian type region would have been predominantly under the influence of colder, depleted water masses sourcing from the north, while during sea level highs, the influence of warmer, more enriched water masses sourcing from the south likely increased (Fig. 1a).

Similar interactions between cooler and warmer water masses have been called upon to explain unexpected distributions in warm-water ammonite fauna in the Coniacian–Santonian of Europe (Remin et al., 2016) and remarkably southerly occurrences of cold-water belemnites in late Cenomanian France (Gale and Christensen, 1996). In this case, the authors suggest that warm waters flowing east to west from the proto-Tethys competed with cooler waters flowing southward between Greenland and Scandinavia (see their Fig. 10). We suggest that a similar mechanism applied to late Maastrichtian paleogeography and water masses (Fig. 1a). Ostracod assemblage changes indicate an increasing contribution of warm Mediterranean (Tethyan) waters in the latest Maastrichtian (Bless, 1988). Faunal abundance differences along the US East Coast indicate the presence of an ancient Gulf Stream-like current as far back in time as the Maastrichtian (Watkins and Self-Trail, 2005), which would have warmed the Maastricht region above its expected temperature for latitude when bathed in water from the south, as the modern-day Gulf Stream does for western Europe.

Fluctuations in δ13C appear to correlate with changes in temperature and δ18Osw (Fig. 6). In this hypothetical model of the alternating influence of southern and northern source waters, this could indicate different δ13C of dissolved inorganic carbon (δ13CDIC) in the two water masses. The exact values of δ13CDIC in either water mass are difficult to define because of overprinting by δ13C vital effects in some taxa (Neithia sp., Entollium sp.; Figs. S2 and S3). Looking at the most common taxa (Acutostrea), variations on the order of 0.5 ‰VPDB – or possibly as much as 1 ‰ – are possible between horizons. This is similar to variability seen in the modern Atlantic Ocean (Eide et al., 2017) and is therefore plausibly driven by alternating water masses of differing δ13CDIC. In this scenario, the colder, northern (Arctic) water mass would be more depleted in δ13CDIC than the southern source water.

This dataset is not capable of distinguishing between the above preferred hypothesis and a temporally variable input of isotopically depleted freshwater. During times of lower sea level, partial basin restriction may have occurred, allowing local freshwater inputs to lower δ18Osw values. In maritime regions such as modern Bermuda, despite lack of elevation or substantial landmass, precipitation and precipitation-derived runoff are depleted in δ18O by 4 ‰–5 ‰ VSMOW relative to seawater (Zhang et al., 2021). Additionally, surface runoff may contain dissolved terrestrial organic matter depleted in δ13C, which could explain the correlation between δ18Osw and δ13C through our study period. However, it is difficult to explain how local freshwater in a subtropical environment surrounded by low-lying landmasses could be substantially colder than local seawater, so we favor the variable water mass hypothesis over a local freshwater contribution.

Yet another possible interpretation of the covariation between fluctuating temperature and δ18Osw trends is via climate-related reorganization of atmospheric circulation, moving broad precipitation bands north and south over the study region. An enhanced hydrological cycle has been hypothesized for the latest Cretaceous (e.g., Woelders et al., 2017) due to increased zonal surface winds and ∼20 % stronger Hadley circulation in the Northern Hemisphere (Bush and Philander, 1997). If, during warmer climates, precipitation bands moved such that the Maastrichtian type region experienced less rainfall, δ18O enrichment of surface waters would coincide with warmer temperatures and vice versa. This is difficult to assess, even in global climate model simulations, as model proficiency in reconstructing paleo-precipitation regimes is difficult to validate and/or calibrate. The rough correlation with δ13C is also difficult to explain in this scenario.

6 Conclusions

This study presents a new clumped-isotope-based paleotemperature time series from the latest Cretaceous in northwestern Europe. The Δ47 analyses from the Maastrichtian type section and surrounding outcrops reveal a mean temperature of 20.4±3.8C and an average δ18Osw of +0.2±0.8 ‰ VSMOW, consistent with a subtropical shallow marine environment. This average temperature aligns well with other paleothermometry studies from similar paleolatitudes and with modeled Maastrichtian paleotemperatures. Increasing Δ47-based temperatures at ∼66.4 Ma, in conjunction with other indicators, are interpreted to reflect a regional manifestation of the global LMWE and directly link pre-KPB warming in northwestern Europe to CO2 emissions from Deccan Traps volcanism. Simultaneous increases in temperature and δ18Osw during the LMWE mask the timing of warming in δ18Ocarb-based temperature reconstructions, which may explain a more delayed apparent timing of the LMWE in other studies. Seen now in two spatially distant records, this deserves further exploration and emphasizes the need to use paleothermometry tools that take into account variations in δ18Osw.

Covariations in temperature, δ18Osw, and (to a lesser degree) δ13C from horizon to horizon may reflect the varying influences of different water masses at different times. We hypothesize that colder waters with lower δ18Osw and δ13CDIC may have originated from the Arctic and reached the study site via a passageway between Greenland and Scandinavia (Fig. 1a) during periods of lower sea level when emergent landmasses to the south acted as barriers. Warmer waters with higher δ18Osw and δ13CDIC may have been sourced from a proto-Tethys to the southeast or proto-Atlantic to the southwest (Fig. 1a) and flushed the study region during intervals of higher sea level.

Overall, the fossils of the Maastricht Formation are effective paleotemperature archives that reveal linkages between local sea level changes, Deccan Traps volcanism, and ocean temperature on the European continent in the lead-up to the end Cretaceous.

Data availability

All data presented in this paper are available in the EarthChem database at (O'Hora et al., 2021). Atmospheric CO2 estimates for the Late Cretaceous can be found in Zenodo (Hoenisch, 2021).

Further data are available in the Supplement as follows.

  • Fig. S1: oxygen isotope values of shell specimens and bulk carbonate

  • Fig. S2: carbon isotope values of shell specimens and bulk carbonate

  • Fig. S3: carbon isotopic composition vs. oxygen isotopic composition in well-preserved shells and bulk matrix

  • Fig. S4: reconstructed seawater oxygen isotopic composition through time

  • Fig. S5: covariation of oxygen isotopic composition of carbonate, water, and temperature

  • Table S1: published paleotemperature data for other studies reconstructing Maastrichtian marine temperatures

  • Table S2: age model of Vellekoop et al. (2022)

  • Table S3: preservation data for all samples, including SEM preservation determinations and trace element concentrations;

  • Table S4: summary of in-house carbonate standards


The supplement related to this article is available online at:

Author contributions

SVP, JV, and HEO'H designed the study. SVP and JV collected the samples. HEO'H, SVP, MMJ, and SRS measured and analyzed the data. HEO'H and SVP prepared the paper, with contributions from all co-authors.

Competing interests

The contact author has declared that none of the authors has any competing interests.


Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.


We thank Ashling Neary (SCIPP lab), Angela Dial (MEAL), Lora Wingate, and Kacey Lohmann (SIL) for assisting in sample preparation and measurements at the University of Michigan. The GSG sample was collected as part of a field trip related to the International Conference on Paleoceanography in 2014, and we thank the organizers and co-participants.

Financial support

This research has been supported by the Geological Society of America (grant no. 12956-20) and the University of Michigan Earth and Environmental Science Department's Turner Research Award to Heidi E. O'Hora, the Research Foundation Flanders (grant no. 12Z6621N) to J. Vellekoop, and the National Science Foundation Ocean Sciences Postdoctoral Research Fellowship (NSF grant no. 1420902) to Sierra V. Petersen.

Review statement

This paper was edited by Denis-Didier Rousseau and reviewed by two anonymous referees.


Anderson, N. T., Kelson, J. R., Kele, S., Daëron, M., Bonifacie, M., Horita, J., Mackey, T. J., John, C. M., Kluge, T., Petschnig, P., Jost, A. B., Huntington, K. W., Bernasconi, S. M., and Bergmann, K. D.: A unified clumped isotope thermometer calibration (0.5–1000 C) using carbonate-based standardization, Geophys. Res. Lett., 48, e2020GL092069,, 2021. 

Barnet, J. S. K., Littler, K., Westerhold, T., Kroon, D., Leng, M. J., Bailey, I., Röhl, U., and Zachos, J. C.: A High-fidelity benthic stable isotope record of Late Cretaceous–Early Eocene climate change and carbon-cycling, Paleoceanogr. Paleoclimatol., 34, 672–691,, 2019. 

Bernasconi, S. M., Mueller, I., Bergmann, K. D., Breitenbach, S. F. M., Fernandez, A., Hodell, D. A., Jaggi, M., Meckler, A.N., Millan, I., and Ziegler, M.,: Reducing uncertainties in carbonate clumped isotope analysis through consistent carbonate-based standardization, Geochem. Geophy. Geosy., 19, 2895–2914,, 2018. 

Bernasconi, S. M., Daëron, M., Bergmann, K. D., Bonifacie, M., Meckler, A. N., Affek, H. P., Anderson, N., Bajnai, D., Barkan, E., Beverly, E., and Blamart, D.: InterCarb: A community effort to improve interlaboratory standardization of the carbonate clumped isotope thermometer using carbonate standards, Geochem. Geophy. Geosy., 22, e2020GC009588,, 2021. 

Bless, M. J. M.: Possible causes for the change in ostracod assemblages at the Maastrichtian-Palaeocene boundary in southern Limburg, The Netherlands, Meded. Werkgr. Tert. Kwart. Geol., 25, 197–211, 1988. 

Bowman, V. C., Riding, J. B., Francis, J. E., Crame, J. A., and Hannah, M. J.: The taxonomy and palaeobiogeography of small chorate dinoflagellate cysts from the Late Cretaceous to Quaternary of Antarctica, Palynology, 37, 151–169,, 2013. 

Brand, U. and Morrison, J. O.: Biogeochemistry of fossil marine invertebrates, Geosci. Can., 14, 85–107, 1987. 

Brand, W. A., Assonov, S. S., and Coplen, T. B.: Correction for the 17O interference in d13C measurements when analyzing CO2 with stable isotope mass spectrometry (IUPAC Technical Report), Pure Appl. Chem., 82, 1719–1733, 2010. 

Brinkhuis, H. and Schiøler, P.: Palynology of the Geulhemmerberg Cretaceous/Tertiary boundary section (Limburg, SE Netherlands), Geol. Mijnbouw, 75, 193–213, 1996. 

Bush, A. B. G. and Philander, S. G. H.: The Late Cretaceous: Simulation with a coupled atmosphere-ocean general circulation model, Paleoceanography, 12, 495–516,, 1997. 

Caldeira, K. and Rampino, M. R.: Carbon dioxide emissions from Deccan volcanism and a K/T boundary greenhouse effect, Geophys. Res. Lett., 17, 1299–1302,, 1990. 

Came, R. E., Brand, U., and Affek, H. P.: Clumped isotope signatures in modern brachiopod carbonate, Chem. Geol., 377, 20–30,, 2014. 

Clyde, W. C., Ramezani, J., Johnson, K. R., Bowring, S. A., and Jones, M. M.: Direct high-precision U–Pb geochronology of the end-Cretaceous extinction and calibration of Paleocene astronomical timescales, Earth Planet. Sc. Lett., 452, 272–280,, 2016. 

de Winter, N. J., Vellekoop, J., Vorsselmans, R., Golreihan, A., Soete, J., Petersen, S. V., Meyer, K. W., Casadio, S., Speijer, R. P., and Claeys, P.: An assessment of latest Cretaceous Pycnodonte vesicularis (Lamarck, 1806) shells as records for palaeoseasonality: a multi-proxy investigation, Clim. Past, 14, 725–749,, 2018. 

Defliese, W. F., Hren, M. T., and Lohmann, K. C.: Compositional and temperature effects of phosphoric acid fractionation on Δ47 analysis and implications for discrepant calibrations, Chem. Geol., 396, 51–60,, 2015. 

Dennis, K. J., Affek, H. P., Passey, B. H., Schrag, D. P., and Eiler, J. M.: Defining an absolute reference frame for `clumped' isotope studies of CO2, Geochim. Cosmochim. Ac., 75, 7117–7131,, 2011. 

Dennis, K. J., Cochran, J. K., Landman, N. H., and Schrag, D. P.: The climate of the Late Cretaceous: New insights from the application of the carbonate clumped isotope thermometer to Western Interior Seaway macrofossil, Earth Planet. Sc. Lett., 362, 51–65,, 2013. 

Dessert, C., Dupré, B., François, L. M., Schott, J., Gaillardet, J., Chakrapani, G., and Bajpai, S.: Erosion of Deccan Traps determined by river geochemistry: impact on the global climate and the 87Sr/86Sr ratio of seawater, Earth Planet. Sc. Lett., 188, 459–474,, 2001. 

Dortangs, R. W., Schulp, A. S., Mulder, E. W., Jagt, J. W., Peeters, H. H., and De Graaf, D. T.: A large new mosasaur from the Upper Cretaceous of The Netherlands, Neth. J. Geosci., 81, 1–8,, 2002. 

Dusar, M. and Lagrou, D.: Cretaceous flooding of the Brabant Massif and the lithostratigraphic characteristics of its chalk cover in northern Belgium, Geol. Belg., 10, 27–38, 2007. 

Eagle, R. A., Eiler, J. M., Tripati, A. K., Ries, J. B., Freitas, P. S., Hiebenthal, C., Wanamaker Jr., A. D., Taviani, M., Elliot, M., Marenssi, S., Nakamura, K., Ramirez, P., and Roy, K.: The influence of temperature and seawater carbonate saturation state on 13C–18O bond ordering in bivalve mollusks, Biogeosciences, 10, 4591–4606,, 2013. 

Eide, M., Olsen, A., Ninnemann, U. S., and Johannessen, T.: A global ocean climatology of preindustrial and modern ocean δ13C, Global Biogeochem. Cy., 31, 515–534,, 2017. 

Eiler, J. M.: Paleoclimate reconstruction using carbonate clumped isotope thermometry, Quaternary Sci. Rev., 30, 3575–3588,, 2011. 

Engelke, J., Linnert, C., Mutterlose, J., and Wilmsen, M.: Early Maastrichtian benthos of the chalk at Kronsmoor, northern Germany: implications for Late Cretaceous environmental change, Palaeobio. Palaeoenv., 97, 703–722,, 2017. 

Felder, W. M.: Lithostratigrafie van het Boven-Krijt en het Dano- Montien in Zuid-Limburg en het aangrenzende gebied, in: Toelichting bij geologische overzichtskaarten van Nederland, edited by: Zafwijn, W. H. and van Staalduinen, C. J., 63–72, Haarlem (Rijks Geologische Dienst), 1975. 

Felder, W. M. and Bosch, P. W.: Geologie van de St. Pietersberg bij Maastricht, Grondboor & Hamer, 52, 53–63, 1998. 

Fendley, I. M., Mittal, T., Sprain, C. J., Marvin-DiPasquale, M., Tobin, T. S., and Renne, P. R.: Constraints on the volume and rate of Deccan Traps flood basalt eruptions using a combination of high-resolution terrestrial mercury records and geochemical box models, Earth Planet. Sc. Lett., 524, 115721,, 2019. 

Font, E., Adatte, T., Sial, A. N., de Lacerda, L. D., Keller, G., and Punekar, J.: Mercury anomaly, Deccan volcanism, and the end-Cretaceous mass extinction, Geology, 44, 171–174,, 2016. 

Gale, A. S. and Christensen, W. K.: Occurrence of the belemnite Actinocamax plenus in the in the Cenomanian of SE France and its significance, B. Geol. Soc. Denmark, 43, 68–77, 1996. 

Gao, Y., Ibarra, D. E., Rugenstein, J. K. C., Chen, J., Kukla, T., Methner, K., Gao, Y., Huang, H., Lin, Z., and Zhang, L.: Terrestrial climate in mid-latitude East Asia from the latest Cretaceous to the earliest Paleogene: A multiproxy record from the Songliao Basin in northeastern China, Earth-Sci. Rev., 216, 103572,, 2021. 

Golovneva, L. B.: The Maastrichtian (Late Cretaceous) climate in the northern hemisphere, Geol. Soc. London Spec. Publ., 181, 43–54,, 2000. 

Haq, B. U.: Cretaceous eustasy revisited, Global Planet. Change, 113, 44–58,, 2014. 

Hart, M. B., FitzPatrick, M. E., and Smart, C. W.: The Cretaceous/Paleogene boundary: Foraminifera, sea grasses, sea level change and sequence stratigraphy, Palaeogeogr. Palaeocl., 441, 420–429,, 2016. 

Henehan, M. J., Ridgwell, A., Thomas, E., Zhang, S., Alegret, L., Schmidt, D. N., Rae, J. W., Witts, J. D., Landman, N. H., Greene, S. E., and Huber, B. T.: Rapid ocean acidification and protracted Earth system recovery followed the end-Cretaceous Chicxulub impact, P. Natl. Acad. Sci. USA, 116, 22500–22504,, 2019. 

Henkes, G. A., Passey, B. H., Wanamaker Jr., A. D., Grossman, E. L., Ambrose Jr., W. G., and Carroll, M. L.: Carbonate clumped isotope compositions of modern marine mollusk and brachiopod shells, Geochim. Cosmochim. Ac., 106, 307–325,, 2013. 

Henkes, G. A., Passey, B. H., Grossman, E. L., Shenton, B. J., Pérez-Huerta, A., and Yancey, T. E.: Temperature limits for preservation of primary calcite clumped isotope paleotemperatures, Geochim. Cosmochim. Ac., 139, 362–382,, 2014. 

Herngreen, G. F. W., Schuurman, H. A. H. M., Verbeek, J. W., Brinkhuis, H., Burnett, J. A., Felder, W. M., and Kedves, M.: Biostratigraphy of Cretaceous/Tertiary boundary strata in the Curfs quarry, Netherlands Institute of Applied Geoscience TNO, 61, 1–58, 1998. 

Hoenisch, B.: Paleo-CO2 data archive (Version 1), Zenodo [data set],, 2021. 

Hull, P. M., Bornemann, A., Penman, D. E., Henehan, M. J., Norris, R. D., Wilson, P. A., Blum, P., Alegret, L., Batenburg, S. J., and Bown, P. R.: On impact and volcanism across the Cretaceous-Paleogene boundary, Science, 367, 266–272,, 2020. 

Huyghe, D., Daëron, M., de Rafelis, M., Blamart, D., Sébilo, M., Paulet, Y. M., and Lartaud, F.: Clumped isotopes in modern marine bivalves, Geochim. Cosmochim, Ac., 316, 41–58,, 2022. 

Jagt, J. W., Donovan, S. K., Fraaije, R., Mulder, E. W., Nieuwenhuis, E., Stroucken, J., van Bakel, B., and van Knippenberg, P.: Remarkable preservation of selected latest Cretaceous macrofossils from the Maastrichtian type area (the Netherlands, Belgium), Foss. Rec., 4, 75–78, 2016. 

Jagt, J. W. M. and Jagt-Yazykova, E. A.: Stratigraphy of the type Maastrichtian – a synthesis, Scr. Geol., 8, 5–32, 2012. 

Johansson, L., Zahirovic, S., and Müller, R. D.: The interplay between the eruption and weathering of large igneous provinces and the deep-time carbon cycle, Geophys. Res. Lett., 45, 5380–5389,, 2018. 

Keutgen, N.: A bioclast-based astronomical timescale for the Maastrichtian in the type area (southeast Netherlands, northeast Belgium) and stratigraphic implications: the legacy of PJ Felder, Neth. J. Geosci., 97, 229–260,, 2018. 

Kim, S. T. and O'Neil, J. R.: Equilibrium and nonequilibrium oxygen isotope effects in synthetic carbonates, Geochim. Cosmochim. Ac., 61, 3461–3475,, 1997. 

Ladant, J. B. and Donnadieu, Y.: Palaeogeographic regulation of glacial events during the Cretaceous supergreenhouse, Nat. Commun., 7, 1–9,, 2016. 

Ladant, J.-B., Poulsen, C. J., Fluteau, F., Tabor, C. R., MacLeod, K. G., Martin, E. E., Haynes, S. J., and Rostami, M. A.: Paleogeographic controls on the evolution of Late Cretaceous ocean circulation, Clim. Past, 16, 973–1006,, 2020. 

Land, L. S.: Diagenesis of skeletal carbonates, J. Sediment Res., 37, 914–930,, 1967. 

Leloux, J.: Numerical distribution of Santonian to Danian corals (Scleractinia, Octocorallia) of southern Limburg, the Netherlands, Geol. Mijnbouw, 78, 191–195,, 1999. 

Li, L. and Keller, G.: Abrupt deep-sea warming at the end of the Cretaceous, Geology, 26, 995–998, 1998. 

Liebau, A.: Paläobathymetrische und paläoklimatische Veränderungen im Mikrofaunenbild der Maastrichter Tuffkreide, Neues Jahrb. Geol. P. A., 157, 233–237, 1978. 

Luijendijk, E., Van Balen, R. T., Ter Voorde, M., and Andriessen, P. A. M.: Reconstructing the Late Cretaceous inversion of the Roer Valley Graben (southern Netherlands) using a new model that integrates burial and provenance history with fission track thermochronology, J. Geophys. Res.-Sol. Ea., 116, p. 19,, 2011. 

Maastricht Climate:, last access: 14 June 2021. 

Mackey, T. J., Jost, A. B., Creveling, J. R., and Bergmann, K. D.: A Decrease to Low Carbonate Clumped Isotope temperatures in Cryogenian Strata, AGU Advances, 1, e2019AV000159,, 2020. 

Meckler, A. N., Ziegler, M., Millaìn, M. I., Breitenbach, S. F. M., and Bernasconi, S. M.: Long-term performance of the Kiel carbonate device with a new correction scheme for clumped isotope measurements, Rapid Commun. Mass Sp., 28, 1705–1715,, 2014. 

Meyer, K. W., Petersen, S. V., Lohmann, K. C., and Winkelstern, I. Z.: Climate of the Late Cretaceous North American Gulf and Atlantic Coasts, Cretaceous Res., 89, 160–173,, 2018. 

Meyer, K. W., Petersen, S. V., Lohmann, K. C., Blum, J. D., Washburn, S. J., Johnson, M. W., Gleason, J. D., Kurz, A. Y., and Winkelstern, I. Z.: Biogenic carbonate mercury and marine temperature records reveal global influence of Late Cretaceous Deccan Traps, Nat. Commun., 10, 1–8,, 2019. 

Miller, K. G., Wright, J. D., and Browning, J. V.: Visions of ice sheets in a greenhouse world, Mar. Geol., 217, 215–231,, 2005. 

Molina, E., Alegret, L., Arenillas, I., Arz, J. A., Gallala, N., Hardenbol, J., Salis, K. von, Steurbaut, E., Vandenberghe, N., and Zaghbib-Turki, D.: The global boundary stratotype section and point for the base of the Danian stage (Paleocene, Paleogene, “Tertiary”, Cenozoic) at El Kef, Tunisia-original definition and revision, Episodes, 29, 263–273, 2006. 

Möller, P. and Kubanek, F.: Role of magnesium in nucleation processes of calcite, aragonite and dolomite, Neues Jb. Miner. Abh., 126, 199–220, 1976. 

Nava, A. H., Black, B. A., Gibson, S. A., Bodnar, R. J., Renne, P. R., and Vanderkluysen, L.: Reconciling early Deccan Traps CO2 outgassing and pre-KPB global climate, P. Natl. Acad. Sci. USA, 118, e2007797118, ,, 2021. 

O'Hora, H. E., Petersen, S. V., Vellekoop, J., Jones, M. M., and Scholz, S. R.: Clumped isotope data from Maastrichtian-aged bivalves from the Maastrichtian type region, Version 1.0, Interdisciplinary Earth Data Alliance (IEDA) [data set],, 2021. 

Percival, L. M., Jenkyns, H. C., Mather, T. A., Dickson, A. J., Batenburg, S. J., Ruhl, M., Hesselbo, S. P., Barclay, R., Jarvis, I., and Robinson, S. A.: Does large igneous province volcanism always perturb the mercury cycle? Comparing the records of Oceanic Anoxic Event 2 and the end-Cretaceous to other Mesozoic events, Am. J. Sci., 318, 799–860,, 2018. 

Petersen, S. V., Dutton, A., and Lohmann, K. C.: End-Cretaceous extinction in Antarctica linked to both Deccan volcanism and meteorite impact via climate change, Nat. Commun., 7, 1–9,, 2016a. 

Petersen, S. V., Tabor, C. R., Lohmann, K. C., Poulsen, C. J., Meyer, K. W., Carpenter, S. J., Erickson, J. M., Matsunaga, K. K., Smith, S. Y., and Sheldon, N. D.: Temperature and salinity of the Late Cretaceous western interior seaway, Geology, 44, 903–906,, 2016b. 

Petersen, S. V., Winkelstern, I. Z., Lohmann, K. C., and Meyer, K. W.: The effects of Porapak™ trap temperature on δ18O, δ13C, and Δ47 values in preparing samples for clumped isotope analysis, Rapid Commun. Mass Sp., 30, 199–208,, 2016c. 

Petersen, S. V., Defliese, W. F., Saenger, C., Daëron, M., Huntington, K. W., John, C. M., Kelson, J. R., Bernasconi, S. M., Colman, A. S., Kluge, T., Olack, G. A., Schauer, A. J., Bajnai, D., Bonifacie, M., Breitenbach, S. F. M., Fiebig, J., Fernandez, A. B., Henkes, G. A., Hodell, D., Katz, A., Kele, S., Lohmann, K.C., Passey, B. H., Peral, M. Y., Petrizzo, D. A., Rosenheim, B. E., Tripati, A., Venturelli, R., Young, E. D., and Winkelstern, I. Z.: Effects of improved 17O correction on interlaboratory agreement in clumped isotope calibrations, estimates of mineral-specific offsets, and temperature dependence of acid digestion fractionation, Geochem. Geophy. Geosy., 20, 3495–3519,, 2019. 

Pucéat, E., Lécuyer, C., Donnadieu, Y., Naveau, P., Cappetta, H., Ramstein, G., Huber, B. T., and Kriwet, J.: Fish tooth δ18O revising Late Cretaceous meridional upper ocean water temperature gradients, Geology, 35, 107–110,, 2007. 

Ravizza, G. and Peucker-Ehrenbrink, B.: Chemostratigraphic evidence of Deccan volcanism from the marine osmium isotope record, Science, 302, 1392–1395,, 2003. 

Remin, Z., Gruszczyński, M., and Marshall, J. D.: Changes in paleo-circulation and the distribution of ammonite faunas at the Coniacian–Santonian transition in central Poland and western Ukraine, Acta Geol. Pol., 66, 107–124,, 2016. 

Roep, T. B. and Smit, J.: Sedimentological aspects of the K/T boundary at Geulhemmerberg, Zuid Limburg, the Netherlands, Geol. Mijnbouw, 75, 119–131, 1996. 

Schaller, M. F., Wright, J. D., Kent, D. V., and Olsen, P. E.: Rapid emplacement of the Central Atlantic Magmatic Province as a net sink for CO2, Earth Planet. Sc. Lett., 323, 27–39,, 2012. 

Schiøler, P., Brinkhuis, H., Roncaglia, L., and Wilson, G. J.: Dinoflagellate biostratigraphy and sequence stratigraphy of the type Maastrichtian (Upper Cretaceous), ENCI Quarry, The Netherlands, Mar. Micropaleontol., 31, 65–95,, 1997. 

Schoene, B., Eddy, M. P., Samperton, K. M., Keller, C. B., Keller, G., Adatte, T., and Khadri, S. F.: U–Pb constraints on pulsed eruption of the Deccan Traps across the end-Cretaceous mass extinction, Science, 363, 862–866,, 2019. 

Schulte, P., Alegret, L., Arenillas, I., Arz, J. A., Barton, P. J., Bown, P. R., Bralower, T. J., Christeson, G. L., Claeys, P., and Cockell, C. S.: The Chicxulub asteroid impact and mass extinction at the Cretaceous-Paleogene boundary, Science, 327, 1214–1218,, 2010. 

Scotese, C. R.: An Atlas of Phanerozoic Paleogeographic Maps: The Seas Come In and the Seas Go Out, Annu. Rev. Earth Pl. Sc., 49, 669–718,, 2021. 

Sial, A. N., Chen, J., Lacerda, L. D., Frei, R., Tewari, V. C., Pandit, M. K., Gaucher, C., Ferreira, V. P., Cirilli, S., and Peralta, S.: Mercury enrichment and Hg isotopes in Cretaceous–Paleogene boundary successions: Links to volcanism and palaeoenvironmental impacts, Cretaceous Res., 66, 60–81,, 2016. 

Smit, J. and Brinkhuis, H.: The Geulhemmerberg Cretaceous/Tertiary boundary section (Maastrichtian type area, SE Netherlands); summary of results and a scenario of events, Geol. Mijnbouw, 75, 283–293, 1996. 

Smit, J. and Zachariasse, W.J.: Planktic foraminifera in the Cretaceous/Tertiary boundary clays of the Geulhemmerberg (Netherlands), Geol. Mijnbouw, 75, 187–191, 1996. 

Sprain, C. J., Renne, P. R., Clemens, W. A., and Wilson, G. P.: Calibration of chron C29r: New high-precision geochronologic and paleomagnetic constraints from the Hell Creek region, Montana, Geol. Soc. Am. Bull., 130, 1615–1644,, 2018. 

Sprain, C. J., Renne, P. R., Vanderkluysen, L., Pande, K., Self, S., and Mittal, T.: The eruptive tempo of Deccan volcanism in relation to the Cretaceous-Paleogene boundary, Science, 363, 866–870,, 2019. 

Stott, L. D. and Kennett, J. P.: The paleoceanographic and paleoclimatic signature of the Cretaceous/Paleogene boundary in the Antarctic: stable isotopic results from ODP Leg 113, in: Proceedings of the Ocean Drilling Program, Scientific Results, College Station, Texas, USA (Ocean Drilling Program), 113, 829–848, 1990. 

Tabor, C. R., Poulsen, C. J., Lunt, D. J., Rosenbloom, N. A., Otto-Bliesner, B. L., Markwick, P. J., Brady, E. C., Farnsworth, A., and Feng, R.: The cause of Late Cretaceous cooling: A multimodel-proxy comparison, Geology, 44, 963–966,, 2016. 

Tagliavento, M., John, C. M., and Stemmerik, L.: Tropical temperature in the Maastrichtian Danish Basin: Data from coccolith Δ47 and δ18O, Geology, 47, 1074–1078,, 2019. 

Tobin, T. S., Ward, P. D., Steig, E. J., Olivero, E. B., Hilburn, I. A., Mitchell, R. N., Diamond, M. R., Raub, T. D., and Kirschvink, J. L.: Extinction patterns, δ18O trends, and magnetostratigraphy from a southern high-latitude Cretaceous–Paleogene section: Links with Deccan volcanism, Palaeogeogr. Palaeocl., 350, 180–188,, 2012. 

Tobin, T. S., Wilson, G. P., Eiler, J. M., and Hartman, J. H.: Environmental change across a terrestrial Cretaceous-Paleogene boundary section in eastern Montana, USA, constrained by carbonate clumped isotope paleothermometry, Geology, 42, 351–354,, 2014. 

Tobin, T. S., Bitz, C. M., and Archer, D.: Modeling climatic effects of carbon dioxide emissions from Deccan Traps volcanic eruptions around the Cretaceous–Paleogene boundary, Palaeogeogr. Palaeocl., 478, 139–148,, 2017. 

Ullmann, C. V., Campbell, H. J., Frei, R., Hesselbo, S. P., Pogge von Strandmann, P. A. E., and Korte, C.: Partial diagenetic overprint of Late Jurassic belemnites from New Zealand: Implications for the preservation potential of δ7Li values in calcite fossils, Geochim. Cosmochim. Ac., 120, 80–96,, 2013. 

Vandenberghe, N., Craen, M. D., and Beerten, K.: Geological framework of the Campine Basin, Studiecentrum voor Kernenergie Centre d’étude de l’énergie Nucléaire, 112, 2014. 

van der Ham, R. W. J. M., van Konijnenburg-van Cittert, J. H. A., Jagt, J. W. M., Indeherberge, L., Meuris, R., Deckers, M. J. M., Renkens, S., and Laffineur, J.: Seagrass stems with attached roots from the type area of the Maastrichtian Stage (NE Belgium, SE Netherlands): Morphology, anatomy, and ecological aspects, Rev. Palaeobot. Palyno., 241, 49–69,, 2017. 

van Hinsbergen, D. J., de Groot, L. V., van Schaik, S. J., Spakman, W., Bijl, P. K., Sluijs, A., Langereis, C. G., and Brinkhuis, H.: A paleolatitude calculator for paleoclimate studies, PLOS One, 10, e0126946,, 2015. 

Vellekoop, J., Sluijs, A., Smit, J., Schouten, S., Weijers, J. W. H., Sinninghe Damste, J. S., and Brinkhuis, H.: Rapid short-term cooling following the Chicxulub impact at the Cretaceous-Paleogene boundary, Proceedings of the National Academy of Sciences, P. Natl. Acad. Sci. USA, 111, 7537–7541,, 2014. 

Vellekoop, J., Esmeray-Senlet, S., Miller, K. G., Browning, J. V., Sluijs, A., Schootbrugge, B. van de, Damsté, J. S. S., and Brinkhuis, H.: Evidence for Cretaceous-Paleogene boundary bolide “impact winter” conditions from New Jersey, USA, Geology, 44, 619–622,, 2016. 

Vellekoop, J., Sluijs, A., and Speijer, R. P.: An acme of the dinoflagellate cyst Palynodinium grallator, Gocht, 1970; a marker for the late Maastrichtian warming event at Northern mid-latitudes?, in: 6th International Geologica Belgica Congress: 12-14 September 2018 – Leuven, Belgium, edited by: Elsen, J., Hulsbosch, N., and Stassen, P., Geologica Belgica Conference Proceedings, 3, 298 p.,, 2018. 

Vellekoop, J., Woelders, L., Sluijs, A., Miller, K. G., and Speijer, R. P.: Phytoplankton community disruption caused by latest Cretaceous global warming, Biogeosciences, 16, 4201–4210,, 2019. 

Vellekoop, J., van Tilborgh, K. H. V., Knippenberg, P. V., Jagt, J. W. M., Stassen, P., Goolaerts, S., and Speijer, R. P.: Type-Maastrichtian gastropod faunas show rapid ecosystem recovery following the Cretaceous–Palaeogene boundary catastrophe, Palaeontology, 63, 349–367,, 2020.  

Vellekoop, J., Kaskes, P., Sinnesael, M., Dehais, T., Huygh, J., Jagt, J., Speijer, R. P., and Claeys, P.: A new age model and chemostratigraphic framework for the Maastrichtian type area (southeastern Netherlands, northeastern Belgium), Newsl. Stratigr., 55, 4, 479–501,, 2022 

Voigt, S., Wilmsen, M., Mortimore, R. N., and Voigt, T.: Cenomanian palaeotemperatures derived from the oxygen isotopic composition of brachiopods and belemnites: evaluation of Cretaceous palaeotemperature proxies, Int. J. Earth Sci., 92, 285–299,, 2003. 

Vonhof, H. B. and Smit, J.: Strontium-isotope stratigraphy of the type Maastrichtian and the Cretaceous/Tertiary boundary in the Maastricht area (SE Netherlands), Geol. Mijnbouw, 75, 275–282, 1996. 

Watkins, D. K. and Self-Trail, J. M.: Calcareous nannofossil evidence for the existence of the Gulf Stream during the late Maastrichtian, Paleoceanography, 20, PA3006,, 2005. 

Wilf, P., Johnson, K. R., and Huber, B. T.: Correlated terrestrial and marine evidence for global climate changes before mass extinction at the Cretaceous-Paleogene boundary, P. Natl. Acad. Sci. USA, 100, 599–604,, 2003. 

Winkelstern, I. Z. and Lohmann, K. C.: Shallow burial alteration of dolomite and limestone clumped isotope geochemistry, Geology, 44, 467–470,, 2016. 

Woelders, L., Vellekoop, J., Kroon, D., Smit, J., Casadío, S., Prámparo, M. B., Dinarès-Turell, J., Peterse, F., Sluijs, A., Lenaerts, J. T. M., and Speijer, R. P.: Latest Cretaceous climatic and environmental change in the South Atlantic region, Paleoceanography, 32, 466–483,, 2017. 

Woelders, L., Vellekoop, J., Weltje, G. J., de Nooijer, L., Reichart, G.-J., Peterse, F., Claeys, P., and Speijer, R. P.: Robust multi-proxy data integration, using late Cretaceous paleotemperature records as a case study, Earth Planet. Sc. Lett., 500, 215–224,, 2018. 

Zhang, J. Z., Petersen, S. V., Winkelstern, I. Z., and Lohmann, K. C.: Seasonally variable aquifer discharge and cooler climate in Bermuda during the Last Interglacial revealed by subannual clumped isotope analysis, Paleoceanogr. Paleoclimatol., 36, e2020PA004145,, 2021. 

Zhang, L., Wang, C., Wignall, P. B., Kluge, T., Wan, X., Wang, Q., and Gao, Y.: Deccan volcanism caused coupled pCO2 and terrestrial temperature rises, and pre-impact extinctions in northern China, Geology, 46, 271–274,, 2018. 

Zhao, M., Ma, M., He, M., Qiu, Y., and Liu, X.: Evaluation of the four potential Cretaceous-Paleogene (K–Pg) boundaries in the Nanxiong Basin based on evidences from volcanic activity and paleoclimatic evolution, Sci. China Earth Sci., 64, 631–641,, 2021. 

Short summary
At the end of the Cretaceous period, massive volcanism in India emitted enough carbon dioxide into the atmosphere to warm the climate globally above an already warm background state. We reconstruct late Cretaceous seawater temperatures much warmer than today using the chemistry of fossil oysters from the modern-day Netherlands and Belgium. Covariations in temperature and water chemistry indicate changing ocean circulation patterns, potentially related to fluctuating sea level in this region.