Oceanic productivity after the Cretaceous/Paleogene impact: Where do we stand? The view from the deep

The Cretaceous/Paleogene (K/Pg) boundary, formerly known as the Cretaceous–Tertiary (K/T) boundary, has attracted the attention of thousands of researchers focused on the study of Earth’s latest mass extinction and its possible causes and patterns of extinction (see Schulte et al., 2010, and Hull et al., 2020, for extensive reviews). As with many other large-scale events, the search for a triggering mechanism for the extinctions and associated paleoenvironmental changes was a long quest, and discoveries commonly came by serendipity. In the 1970s, Walter Alvarez and William Lowrie were determining the geomagnetic polarity sequence of the Upper Cretaceous–lower Cenozoic pelagic limestones at Gubbio, Italy, to provide a time scale for seafloor-spreading magnetic anomalies (Lowrie and Alvarez, 1975, 1981). The biostratigraphy was based on planktic foraminifera studied by Isabella Premoli Silva and Hanspeter Luterbacher (Luterbacher and Premoli Silva, 1962, 1964), and on calcareous nannofossils (Premoli Silva et al., 1976; Monechi, 1979). Isabella showed Walter the thin clay layer at the Cretaceous–Paleogene boundary, and he took samples, hoping to find a way to determine how much time the clay represented. Back in Berkeley, working with his father Luis Alvarez, Frank Asaro, and Helen Michel, he measured the iridium content of the clay. The Earth’s crust is depleted by orders of magnitude in this siderophilic element relative to the core and mantle of the Earth (i.e., the planetary interior) and relative to some differentiated and undifferentiated meteorites such as carbonaceous chondrites (e.g., Goderis et al., 2013, 2021). The Ir in sedimentary layers of the Earth’s crust mostly comes from the rain of micrometeorites from outer space (e.g., Campos, 1997). The Alvarez team assumed that, if this rain were constant, the Ir content in the clay layer would allow them to estimate sedimentation rates. Unexpectedly, they found that the Ir concentration in the K/Pg boundary clay was much higher than could be explained by the normal fall of micrometeorites, and they proposed that the Ir anomaly resulted from the impact of a huge asteroid (~10 km in diameter) on Earth at the K/Pg boundary (Alvarez et al., 1980). Additionally, they postulated that effects of the impact (associated climatic changes, specifically the effect of a large input of dust in the atmosphere) were the cause of the extinction of the biota at the K/Pg boundary. At about the same time, Smit and Hertogen (1980) independently published similar conclusions based on their study of the Caravaca section in Southern Spain. The hypothesis that the impact of an extraterrestrial body caused the K/Pg extinctions was not new (e.g., De Laubenfels, 1956; Russell, 1977, 1979), but the Ir anomaly was the first physical evidence to support such a hypothesis. Nowadays, it is widely acknowledged that the asteroid impact (Alvarez et al., 1980; Smit and Hertogen, 1980) occurred on the Yucatan peninsula (e.g., Hildebrand et al., 1991; Smit et al., 1992, 1996; Morgan et al., 2016) and was the main cause of the catastrophic mass extinctions (e.g., Schulte et al., 2010; Hull et al., 2020; Bralower et al., 2020a). The globally recognized anomalous concentrations of Ir, and the subsequently discovered impact ejecta, including shocked minerals, Ni-rich spinels, and microspherules in K/Pg boundary sediments reflect that impact (e.g., Alvarez et al., 1980, 1992; Smit and Hertogen, 1980; Alegret et al., 2005; Molina et al., 2009; Goderis et al., 2013; Schulte et al., 2010).

The impact event and its environmental and biological consequences left a global geochemical imprint in oceanic records, namely the collapse of the oceanic bottom to surface gradient in carbon isotope values, as measured in benthic foraminiferal and planktic tests and bulk (nannoplankton) calcite (e.g., Hsü et al., 1982; Hsü and McKenzie, 1985; Zachos et al., 1989). This vertical gradient represents the difference in the carbon isotope signal of benthic and planktic fossil shells, and supposedly dissolved inorganic carbon (DIC) in the deep and surface waters of the ocean, where these organisms calcify. The collapse in gradient at the K/Pg boundary has been traditionally explained by a severe drop in oceanic primary and/or export productivity lasting from hundreds of thousands to a few million years after the impact (e.g., Hsü et al., 1982; Hsü and McKenzie, 1985; Zachos and Arthur, 1986; Arthur et al., 1987; Zachos et al., 1989; D’Hondt et al., 1998). As we describe below, these carbon isotope-based hypotheses are not in full agreement with various paleontological and organic geochemical records (e.g., Alegret et al., 2012; Sepúlveda et al., 2019; Bralower et al., 2020b).

Extinctions on land (e.g., non-avian dinosaurs, flying reptiles, and 30–80% of plants) and in the oceans (an estimated 47–53% extinction of genera, extrapolated to ~76% extinction of species [Jablonski, 1994; Alroy, 2008], e.g., rudists, ammonites, mosasaurs, and >90% of the species of planktonic foraminifera and coccolithophores [e.g., Lowery et al., 2020]) were selective, and there was an abrupt, global biological turnover (see Schulte et al., 2010, for a review).

The primary cause of this mass extinction was debated for decades (and to some extent, still is; Keller et al., 2020), at least in part due to the occurrence of the rapid (nearly instantaneous) asteroid impact at some time during the long-term time interval of flood basalt volcanism in the Indian Deccan Traps, although there is no direct observation of the extinction in sediments interlayered with Deccan flows (Fendley et al., 2021). In addition, and despite detailed radiometric age determination, there is no agreement on the exact timing of the extinction relative to the eruption history of the Deccan Trap lava flows as derived from numerical dating (compare Sprain et al., 2019; Schoene et al., 2019, 2021; Burgess, 2019; Basu et al., 2020), but the K/Pg extinction did occur in the ~1 m.y. interval (65.4–66.4 Ma) of eruption of the largest flows in the Deccan Traps (Burgess, 2019; Sprain et al., 2019; Schoene et al., 2021). It was speculated that the impact at Chicxulub triggered an episode of intensive magmatism in the Deccan region roughly at the other side of the Earth (e.g., Richards et al., 2015), a modification of an earlier hypothesis that an impact triggered eruptions at the impact location (Rampino, 1987; Negi et al., 1993).

Numerous studies pointed out the coincidence in timing and the cause-and-effect relationship between the impact and extinctions (see Schulte et al., 2010, and Hull et al., 2020, for reviews). Carbon cycle modeling, combined with global paleotemperature records, helped to disentangle the relative importance of the asteroid impact and volcanism, and led to the conclusion that the impact alone was the primary driver of the mass extinction (including that of the dinosaurs; Chiarenza et al., 2020) and associated carbon-cycle change (Hull et al., 2020). Modeling and records also documented that major outgassing began distinctly before the impact and persisted long after it (Hull et al., 2020), in agreement with the ages of Deccan flows (e.g., Schoene et al., 2021). There was a volcanic emission-linked period of global warming before the end of the Cretaceous (e.g., Li and Keller, 1998; Barnet et al., 2018), but there were no significant extinctions above background in this period (review in Hull et al., 2020). Deccan Trap gas emissions, however, may have affected recovery from the extinctions, especially by pelagic calcifiers (Hull et al., 2020).

Since the proposal of the impact hypothesis in 1980, the scientific community has thoroughly documented the relation between the impact of a large asteroid on the Yucatan peninsula, the widely distributed impact-derived material, global geochemical anomalies, and mass extinctions (see, e.g., reviews by Kring, 2007; Schulte et al., 2010), and the crater has been drilled twice (Urrutia-Fucugauchi et al., 2004; Morgan et al., 2016). However, the exact killing mechanism(s) linked to the various environmental effects of the impact on different time scales have been, and still are, the topic of an active, lively, highly interdisciplinary and interesting debate. There are many, not all mutually exclusive, proposed mechanisms of extinction (see, e.g., Kring, 2007, for a discussion), including global darkness due to the injection of dust and volatiles (participating in gas-to-particle reactions) into the atmosphere (e.g., Alvarez et al., 1980; Morgan et al., 1997), resulting in a severe decline in photosynthesis, and thus the potential collapse of oceanic primary productivity (e.g., Hsü et al., 1982; Hsü and McKenzie, 1985; Vellekoop et al., 2017; Gibbs et al., 2020). The darkness and cold (“impact winter”) could have persisted for several decades (e.g., Brugger et al., 2017; Vellekoop et al., 2014), possibly in part caused by soot from the fires (Tabor et al., 2020; Lyons et al., 2020). However, oceanographic changes induced by the short impact winter may have triggered a prolonged, 2 k.y. cooling phase in the western Tethys (Galeotti et al., 2004). Other mechanisms include ozone destruction (see references in Kring, 2007) and severe global heat stress immediately after the collision due to an infrared radiation pulse caused by transit through the atmosphere of recondensed ejecta (e.g., Morgan et al., 2013), which could have triggered global wildfires (e.g., Wolbach et al., 1990), though their extent is debated (e.g., Belcher et al., 2005). Oceanic environments could have been affected by acidification of the ocean due to CO₂ emissions caused by the vaporization of limestone in the impact area or, more likely, by the vaporization of sulfate (Alegret et al., 2012; Ohno et al., 2014; Tyrrell et al., 2015; Henehan et al., 2019).

In attempting to explain global extinctions, however, one should keep in mind that not all environmental effects of the impact were global, and therefore regional differences in extinction patterns can be expected. As an example, regionally, terrestrial biota may have been affected by tsunamis running far inland, e.g., along the U.S. Western Interior Seaway up to what is now North Dakota (DePalma et al., 2019). Mass wasting along continental margins may have affected bottom-dwelling biota regionally, e.g., in the Gulf of Mexico and western Atlantic (Smit et al., 1992; Alegret et al., 2002a; Bralower et al., 2002; Norris and Firth, 2002; Sanford et al., 2016). In addition, the effects of extinctions themselves in disrupting ecosystems, e.g., through the extinction of keystone species, which promoted further extinctions, have not been thoroughly explored (e.g., Plotnick and McKinney, 1993; Hull, 2015). The many, in part, interdependent, environmental, and ecological effects of the impact of a large asteroid acted at very different time scales, from instantaneous (expanding impact fireball, fires, and mass wasting) to hundreds of thousands or even a few million years (e.g., collapse of the marine carbon isotope vertical gradient) (see references in Hull et al., 2020).

Here, we briefly report on some aspects of the history of this debate, particularly focusing on extinctions in the marine realm and the evolution of oceanic primary and export productivity and its global geochemical signals. In the years directly after the publication of the asteroid impact hypothesis (Alvarez et al., 1980; Smit and Hertogen, 1980), marine micropaleontological studies focused on the extinction patterns of the most severely affected groups, the pelagic calcifiers (calcareous nannoplankton and planktic foraminifera; e.g., Premoli Silva et al., 1976; Romein and Smit, 1981; Smit, 1982; Pospichal, 1994; Lipps, 1997; Molina et al., 1998; Bown, 2005a, 2005b; Jiang et al., 2010) and their recovery. More recently, however, detailed studies of groups that were less significantly affected by the K/Pg extinction, such as the siliceous-walled diatoms and radiolarians (Harwood, 1988; Hollis et al., 2003), have significantly contributed to this debate. Possibly by chance—like Alvarez´s discovery of the Ir anomaly—the study of deep-sea dwellers (benthic foraminifera) has led to major advances in understanding the effect of the impact on primary as well as export productivity. Studies aimed at understanding what happened in the deep oceans serendipitously led to implications for observed changes in surface ocean pelagic communities and support for the hypothesis that impact-derived acidification of the ocean was an important mechanism for selective extinction of calcifiers (e.g., Alegret et al., 2012).

The deep sea is the largest and one of the most stable habitats on Earth, and deep-sea dwellers are adapted to the extreme, hostile conditions at the seafloor (darkness, cold, lack of food due to arrival of only 1–3% of organic matter produced at the surface), they have long species lives, and provide an excellent record of past environments (e.g., references in Alegret et al., 2021). Among deep-sea dwellers, benthic foraminifera are the most abundant meiobenthos (Gooday, 2003), and it has long been known that they did not suffer significant extinction after the end-Cretaceous impact (e.g., Beckmann, 1960; Emiliani et al., 1981; Douglas and Woodruff, 1982; Miller, 1982; Thomas, 1990a, 1990b; Culver, 2003; Alegret and Thomas, 2005). Their tests thus provide a continuous biotic and isotopic record across the K/Pg transition (e.g., Alegret et al., 2012). Here, we contribute new quantitative data on benthic foraminiferal assemblages and stable isotope data to a compilation of published data and outline multidisciplinary results that help us draw a route map for future K/Pg studies and data re-evaluation specifically related to marine productivity.

To trace the environmental consequences of the K/Pg impact event in deep-sea environments, we look into benthic foraminifera and their assemblage turnover across the boundary at multiple geographic and bathymetric locations. To document changes in these different paleogeographical settings, we compiled published and unpublished quantitative studies that use similar taxonomic concepts and methodologies (size fraction, number of picked foraminifera, and type of quantitative studies). The selection of comparable studies is essential for the integration of faunal data because of the many serious problems in integrating the varying taxonomic concepts used by different investigators (e.g., Arreguín-Rodríguez et al., 2018; Alegret et al., 2021). After careful examination of data sets, we integrated information from 17 oceanic drilling sites and land sections (Fig. 1) in the Pacific and Indian Oceans, the Southern Ocean, the Atlantic Ocean (including sites near the Chicxulub impact crater in the Gulf of Mexico), and the Tethys Ocean (Table 1). Benthic foraminiferal data sets from all sites are based on study of the >63 µm size fraction, with the exception of Nomura (1991), which is based on the >149 µm size fraction. We included this study because it is the only one available from the Indian Ocean. We are cautious in including these data, however, because they differ from all other data sets.

Paleodepth of the sites ranges from outer shelf through upper bathyal at the three Tunisian sections, to upper abyssal (Table 1). A wide range of depositional settings is therefore represented, including open-ocean settings (e.g., Pacific and SE Atlantic sites), Flysch deposits (sections from the Basque-Cantabrian Basin, NE Atlantic), and complex sequences with K/Pg units affected by impact-generated tsunamis, seismic-induced slumping, and other mass wasting processes in the Gulf of Mexico (e.g., Mexican sections; Alegret et al., 2002a, 2002b; Sanford et al., 2016).

The relative abundance of the buliminids sensu lato (s.l.) group (Alegret and Thomas, 2013), which includes biserial and triserial elongate genera of the superfamilies Buliminacea, Bolivinacea, Loxostomatacea, Turrilinacea, Fursenkoinacea, Pleurostomellacea, and Stilostomellacea (Sen Gupta, 1999), was used as a proxy for trophic conditions at the seafloor. In the modern oceans, this group tolerates low oxygen concentrations (e.g., Sen Gupta and Machain-Castillo, 1993; Bernhard et al., 1997), but the major influence on their abundance may be an abundant and continuous food supply, which commonly occurs at lower oxygen values (e.g., Fontanier et al., 2002). Additionally, we document changes in the relative abundance of infaunal (living deeper within the sediment) and epifaunal (living at the seafloor or in the uppermost centimeters of the sediment) habitat-related morphogroups in calcareous and agglutinated taxa (Corliss, 1985; Jones and Charnock, 1985; Corliss and Chen, 1988). The relative abundance of these morphotypes has been commonly used as a proxy for environmental conditions at the seafloor, with infaunal taxa indicating lower oxygenation and/or higher trophic conditions (e.g., Jorissen et al., 1995, 2007). Interpretation of morphogroups is somewhat problematic, however, even for living foraminifera (Buzas et al., 1993), and even more so for fossil assemblages, mostly due to the lack of modern analogs (e.g., Hayward et al., 2012).

For all sites and sections of our compilation, we quantified the overall drop in diversity of benthic foraminifera across the K/Pg boundary (Fisher-α index; Murray, 2006). Changes in the relative abundance of buliminid taxa and infaunal morphogroups across the boundary were calculated as the difference between the average values during the last 170–200 k.y. of the Cretaceous (late Maastrichtian, planktic foraminiferal Plummerita hantkeninoides Biozone) and the first 50 k.y. of the Paleogene (early Danian, P0/G. cretacea and Pv. eugubina Biozones) (Table 1; Fig. 2). The percentage of species that last occur at the K/Pg transition was calculated relative to the number of species present in the uppermost Cretaceous; species that reappear later in the Danian were considered Lazarus taxa, and only species that show their highest occurrence at the K/Pg boundary globally as far as we can determine from available data were used to calculate extinction rates (Table 2). These include species that had their uppermost appearance immediately at or a few centimeters below the boundary, to account for the Signor-Lipps effect. Species that show sparse, rare occurrences a few centimeters above the boundary were included to account for reworking processes.

Additionally, we document carbon and oxygen stable isotopes of bulk sediment and of benthic foraminifera, diversity (Fisher-α) of benthic foraminiferal assemblages, benthic foraminiferal accumulation rates (BFARs, number of foraminifera per cm² per k.y., a proxy for export productivity to the seafloor; Herguera and Berger, 1991; Jorissen et al., 2007), and the percentage of infaunal morphogroups and the Superfamily Buliminacea across the uppermost Cretaceous and lower Paleogene at oceanic drilling sites along a paleolatitude transect from the North Atlantic Ocean to the Southern Ocean (Sites 1049, 1267, 1262, and 690; Fig. 3) and at two sites in the Pacific Ocean (Sites 465 and 1210; Fig. 4) (Alegret and Thomas, 2005, 2007, 2009, 2013; Alegret et al., 2012).

Our compilation includes new data from Ocean Drilling Program (ODP) Site 1267 (Walvis Ridge, eastern South Atlantic Ocean) and from the Urrutxua section in NE Spain (Basque-Cantabrian Basin in the NE Atlantic), (Supplemental Material Tables S1–S2¹). The K/Pg transition was cored along a depth transect on ODP Leg 208 (Angola Basin, near the Walvis Ridge), with the boundary marked by a sharp transition from Maastrichtian clay-bearing nannofossil ooze to Danian dark reddish to brown, clay-rich nannofossil ooze and clay. We analyzed the benthic foraminiferal turnover at Site 1267 (lower bathyal; present depth 4355 m) and compared it with published results from Site 1262 (upper abyssal paleodepth; present depth 4755 m) (Alegret and Thomas, 2007). The record across the K/Pg boundary at Site 1267 appears to be more complete than at Site 1262. There appears to be a short unconformity at Site 1262, which corresponds to a thin layer with small spherules at Site 1267 (Westerhold et al., 2008; Bralower et al., 2020b). No spherules were observed at Site 1262. Here we show results from the uppermost 2.14 m of the Maastrichtian and the lowermost 6.77 m of the Danian at Site 1267 (Table S1) using the cyclostratigraphic age model in Hull et al. (2020). For comparison with published data from Sites 1262, 1210, 465, and 690 (Alegret et al., 2012), we shifted age estimates in that publication by 0.5 m.y. to match the K/Pg boundary age (66.00 Ma in the Geologic Time Scale 2020; Speijer et al., 2020).

Benthic foraminiferal data are also presented from the uppermost Maastrichtian (planktic foraminiferal Plummerita hantkeninoides Biozone) and the lowermost Danian (P0/G. cretacea and Pv. eugubina Biozones) from section Urrutxua, a hemipelagic section in the Basque-Cantabrian Basin (Table S2). During the Late Cretaceous–early Paleogene, this was a deep basin surrounded on its north, south, and east sides by shallow marine carbonate platforms flanked by continental alluvial plains (Alegret et al., 2004a). The K/Pg transition at Urrutxua is similar to that at the better known, nearby Zumaia section (Baceta and Pujalte, 2006): the upper Maastrichtian gray and reddish marls and marlstones contain intercalations of thin, mixed carbonate-siliciclastic turbidites, and the K-Pg clay level varies in thickness between 1 cm and 8 cm due to inter-bed sliding (as inferred from the presence of calcite slickensides). The lower part of the clay level consists of a 0.5–0.8-mm-thick yellowish silty bed with partially dissolved foraminiferal tests, which is followed upward by a reddish to brownish siltstone with microspherules with Ni-rich spinel crystals. An iridium anomaly is recorded between this level and the overlying, thinly laminated dark gray siltstones (Rocchia et al., 1996; Baceta and Pujalte, 2006).

Samples were processed following standard micropaleontological procedures, using H₂O with soap (samples from ODP Site 1267) and H₂O₂ (samples from Urrutxua) to disaggregate the sediment. Then samples were washed with running water over a sieve with a 63 µm net wire (Urrutxua). Those from ODP Site 1267 were run over two sieves (38 µm and 63 µm). Quantitative studies were based on representative splits of ~300 benthic foraminifera per sample, using the >63 µm size fraction of the washed residue. The Fisher-α diversity index and the percentages of infaunal morphogroups and buliminids s.l. were calculated for both locations (Fig. 2; Tables S1–S2). To compare results from ODP Site 1267 with those from the other ocean drilling sites in our compilation (Figs. 3–4), we calculated BFARs (Table S1), diversity of the assemblages, and the percentage of infaunal morphogroups and the Superfamily Buliminacea. We analyzed mono-specific, carbon and oxygen stable isotope records in benthic foraminiferal tests, generated at the University of California, Santa Cruz, after Alegret et al. (2012). For Site 1267, bulk sediment δ¹⁸O and δ¹⁸C and benthic foraminiferal δ¹⁸O were published in Hull et al. (2020), and benthic foraminiferal (Nuttallides truempyi) δ¹³C records are presented here (Table S3; see ^{footnote 1}).

Deep-sea benthic foraminifera underwent temporary changes in assemblage composition but were not significantly affected by extinction across the K/Pg boundary (Beckmann, 1960; Emiliani et al., 1981; Douglas and Woodruff, 1982; Miller, 1982), unlike many other groups of marine organisms. In fact, Cushman (1946) placed the Cretaceous–Paleogene boundary as defined by foraminiferal extinction at the end of the Paleocene, when deep-sea benthic foraminifera did suffer a major extinction (e.g., Thomas, 1990b). He located the major boundary there because total foraminifera (the sum of planktic + benthic) show a much larger species turnover at the end of the Paleocene than at the end of the Cretaceous, because benthic species are so much more numerous than planktic species (though with much lower numbers of specimens). The extensive review by Culver (2003) and further high-resolution, quantitative assemblage studies provided detailed data, though not at high time-resolution for all sites, which supported the earlier observations that out of the many hundreds to a few thousands of benthic species living at any time period in the deep sea (Murray, 2007), relatively few benthic foraminiferal species suffered extinction across the boundary. In spite of these well-documented low extinction rates, at some sites not above background extinction rates (Thomas, 1990b), late Maastrichtian assemblages were different from those in the early Danian, as seen in a global drop in diversity and changes in the relative abundance of species and habitat-related morphogroups (Figs. 2–3) (e.g., Kuhnt, 1990; Thomas, 1989, 1990a, 1990b; Widmark and Malmgren, 1992; Kuhnt and Kaminski, 1993; Coccioni and Galeotti, 1994; Alegret et al., 2001, 2005; Alegret and Thomas, 2005, 2007, 2009, 2013; Peryt et al., 2002; Coccioni and Marsili, 2007).

Decreased benthic foraminiferal diversity indicates that the K/Pg impact caused temporary environmental stress at the seafloor globally, from the shelf to abyssal depths, and at all latitudes, though the effects—and thus possibly the type of stress—varied by location (Fig. 2, Table 1). Paleodistance from the impactor determined the sedimentary complexity of the of K/Pg boundary units, with expanded ejecta units, tsunami and mass wasting deposits, and unconformities clearly recognizable at sites closer to the Chicxulub impact crater (e.g., Alegret, 2003; Molina et al., 2009; Schulte et al., 2010; Sanford et al., 2016). The global drop in benthic foraminiferal diversity, however, shows poor correlation with paleodistance from the impact site (R² value = 0.186), paleolatitude (R² = 0.041), and paleodepth (R² = 0.172) (Fig. S1; see ^{footnote 1}). This is as expected, because many benthic foraminiferal species are cosmopolitan and can repopulate regions of the oceans where they were wiped out due to their motile propagules (e.g., Alve and Goldstein, 2010).

Benthic foraminifera at shallower settings had been suggested to suffer more severe extinction than those living at greater paleodepths (Thomas, 1990b; Kaiho, 1992; Widmark, 1997; Coccioni and Galeotti, 1998), but Culver (2003) concluded that paleodepth from shelf to abyss did not significantly affect the extinction rates of benthic foraminifera.

The evaluation of extinction rates, however, is problematic because estimates by various authors are derived using different criteria, and thus are difficult to compare between sites. These include: time resolution of the study, number of species that disappear at the boundary, their percentage relative to the total number of species within the latest Cretaceous assemblage, relative abundance within the assemblage, and considerations regarding reworking in the lowermost Danian and abundant Lazarus taxa that disappear at the boundary and reappear in the lower Danian, etc. (e.g., Thomas, 1990b). Deep-sea benthic foraminiferal assemblages (like these of other deep-sea biota) typically are species-rich and have distributions with relatively few common species and many rare species (e.g., Douglas and Woodruff, 1982; Murray, 2007). Thus, it is difficult to statistically define whether species became extinct exactly at the boundary, as observed stratigraphic ranges are not continuous for rare species (e.g., Thomas, 1985, 1990a).

We used normalized taxonomic concepts and methodology to calculate the number and percentage of species that last occurred at the K/Pg transition (Table S4; see ^{footnote 1}). Those that were reported to reappear in younger sediments at any location were considered to be Lazarus taxa. Local extinction rates range from 1.6% at ODP Site 690 in the Southern Ocean to 8.9% at La Lajilla, Mexico (Table 2). Thomas (1990b) estimated 14% extinction rates in the Tasman Sea (Lord Howe Rise DSDP Site 208) based on Webb (1973). Coccioni and Marsili (2007) documented the extinction of 3% of benthic species across the K/Pg boundary at Elles, Tunisia, but we arrived at a 5.9% extinction rate for that location by using our methods with their data. This number is similar to extinction rates at other Tunisian sections (6.1% both at El Kef and Aïn Settara) (Table 2; Table S4). Using our methods to standardize estimates of extinction, shallower sites in general do not appear to have been more affected than deeper ones, and the R² values between rate of extinction and water depth are not significant (0.169; Fig. S2; see ^{footnote 1}). Our compilation, however, includes three neritic to outer bathyal (< 500 m) sites only, and more comparable data sets are needed from shallow settings. Similarly, the last occurrences of species (i.e., including the temporary disappearance of Lazarus taxa) at the K/Pg boundary show no significant correlation with paleodepth, distance from the crater, or paleolatitude (Fig. S2).

In conclusion, we see no statistically significant correlations between the local severity of species extinction and water depth, distance from the impact crater, and paleolatitude, and neither do we see a correlation between the local decrease in diversity and water depth, distance from the impact crater, and paleolatitude (Fig. S1) at the time resolution of the compiled studies.

To estimate global extinction of benthic foraminifera, we excluded taxa with open nomenclature included in our compilation (Table S4) as well as species that have been reported from younger sediments. As a result, we identified only eight species (out of the hundreds of deep-sea benthic foraminiferal species) that went extinct globally at the K/Pg boundary (Table S5; see ^{footnote 1}).

The relative abundance of infaunal versus epifaunal morphogroups across the K/Pg boundary has been more unequivocally documented than the extinction percentage and allows for more robust comparisons among studies by different authors, because it is mostly based on the morphology of the tests rather than on species-level taxonomy, and thus it is less affected by taxonomic concepts than, e.g., estimates of extinction rates. The relative abundance of infaunal morphogroups decreased after the K/Pg boundary at most localities studied, but showed little change in the SE Atlantic and Southern Ocean and significantly increased in some locations in the Pacific Ocean and in the western Indian Ocean (Fig. 2). At most sites, these changes are associated with a coeval increase or decrease in abundance of the infaunal buliminids (Fig. 2), except for the Indian Ocean site (Nomura, 1991), where infaunal morphogroups increase whereas buliminids decrease across the boundary. In addition, at two sites (Bidart and Coxquihui) the drop in infaunal morphogroups is not correlated to significant changes in the percentage of buliminids.

Changes in habitat-related morphogroups have been traditionally interpreted in terms of oxygenation and/or trophic conditions at the seafloor, with the dominance of epifaunal taxa suggesting more oligotrophic conditions during the earliest Danian in the Tethys (Coccioni et al., 1993; Coccioni and Galeotti, 1998; Peryt et al., 2002; Coccioni and Marsili, 2007; Alegret, 2008), NE Atlantic (Basque-Cantabric Basin; Alegret et al., 2004a; Alegret, 2007), and the Gulf of Mexico and proto-Caribbean (Alegret, 2003; Alegret and Thomas, 2004; Alegret et al., 2001, 2002a, 2002b, 2005). In contrast, a marked increase in infaunal taxa, coupled with significantly increased BFARs and buliminids (Figs. 2 and 4), has been interpreted to indicate an increase in arrival of food to the seafloor at Pacific Sites 465 and 1210 (Alegret and Thomas, 2005, 2009; Alegret et al., 2012). The lack of significant changes in the abundance of morphogroups and BFARs (Fig. 3), combined with the species turnover, point to changes in the character of food supply (e.g., type of food, thus possibly the type of phytoplankton, and temporal variability in flux) rather than in the amount of food supply to the seafloor in the Southern Ocean (Alegret and Thomas, 2013) and SE Atlantic Site 1262 (Alegret and Thomas, 2007).

Low-oxygen conditions, based on the interpretation of changes in the relative abundance of morphogroups and species composition of the assemblages (in part supported by geochemical studies; e.g., Smit, 1990; Martínez-Ruiz et al., 1992), have been reported regionally (e.g., Agost and Caravaca sections in the Western Tethys; Coccioni et al., 1993, Coccioni and Galeotti, 1994; Alegret et al., 2003). At the Tethyan section of El Kef, Speijer and Van der Zwaan (1996) inferred a sudden drop in productivity (in primary productivity at the surface and in arrival of food at the seafloor) combined with reduced bottom-water oxygenation, from the combination of sedimentary lamination and the short-term increase in abundance of the infaunal Praeglobobulimina quadrata at the base of the boundary clay at El Kef (Speijer and Van der Zwaan, 1996; Alegret, 2003), as also seen at Elles (Coccioni and Marsili, 2007), although geochemical data do not unequivocally indicate low oxygenation. The lowermost Danian assemblages just above this interval at El Kef, Aïn Settara (Galeotti and Coccioni, 2002; Alegret et al., 2004b) and Elles (Coccioni and Marsili, 2007), are very low diversity and strongly dominated by Cibicidoides pseudoacutus. This robust species has large pores on one side, and it belongs to the oxic group as defined by Corliss and Chen (1988). It somewhat resembles the living epifaunal Cibicidoides wuellerstorfi in morphology and thus possibly in habitat (i.e., epifaunal or on objects sticking out above the sediment). Thus, there are no benthic faunal indications of low oxygen levels during this period, and the type of environmental stress causing the low diversity during this interval in the Tunisian sections remains problematic.

The occurrence of several subsequent peaks in the abundance of different species of benthic foraminiferal taxa that might have been opportunistic (e.g., Coryphostoma incrassata gigantea and Spiroplectammina spectabilis) at multiple sites during the earliest Danian has been related to the post-extinction proliferation of different photosynthesizers, which may have triggered changes in the nature rather than in the amount of the organic matter supplied to the seafloor (e.g., Alegret, 2007; Alegret and Thomas, 2007, 2013; Alegret et al., 2004a, 2012).

In addition, carbonate supersaturation in the oceans after the extinction of calcareous pelagic calcifiers may have given competitive advantage, at least locally, to species with large and heavily calcified tests (Alegret and Thomas, 2013). Possibly, the peak in C. pseudoacutus in the North African sections could be explained that way, but a more probable indicator species of higher carbonate saturation is Stensioeina beccariiformis, which increased in abundance above the post-K/Pg boundary for a short time at multiple locations globally (Alegret et al., 2021), with an especially notable occurrence at Pacific Site 1210 (Alegret and Thomas, 2009). This common, cosmopolitan species became extinct during the Paleocene Eocene Thermal Maximum benthic foraminiferal extinction, at a time of ocean acidification (e.g., Thomas, 1998; Penman et al., 2014).

Interpretation of benthic morphogroups should be done cautiously, however, because (1) their use has limitations even in the modern environments, and assignments from morphology have been found to be correct for living foraminifera only in ~75% of cases (Buzas et al., 1993); (2) common, now extinct species in Cretaceous and Paleogene faunas are non-analog to living species in morphology (Thomas and Gooday, 1996; Thomas, 2007; Hayward et al., 2012; Arreguín-Rodríguez et al., 2018), and thus there are no direct observations of their ecological affinities; (3) an in-depth analysis at the species level reveals significant abundance peaks of different taxa within the same morphogroup, at constant total abundance of that morphogroup (e.g., Site 690; Alegret and Thomas, 2009); and (4) their record (and interpretation) strongly relies on preservation of the sedimentary record, as shown in examples A and B in the discussion below. In addition to these examples, the sedimentary record of the K/Pg transition in localities around the Gulf of Mexico is incomplete due to the effects of impact-linked sediment disturbances such as seismic instability, tsunamis, submarine landslides, slumps, and other mass wasting processes (e.g., Bralower et al., 1998; Klaus et al., 2000; Soria et al., 2001; Alegret et al., 2002a; Alegret and Thomas, 2005; Sanford et al., 2016), which prevents environmental interpretation over several thousands of years of the post-impact trophic and oxygen conditions at the seafloor.

Our results from the Urrutxua section confirm an apparent decrease in the percentage of infaunal morphogroups across the K/Pg boundary in the Basque-Cantabrian Basin, as recorded in the nearby Zumaia section (Alegret and Ortiz, 2010) and the Loya and Bidart sections in Southwestern France (Alegret, 2007; Alegret et al., 2004a). However, benthic foraminifera inside trace fossils that cross the K/Pg boundary at Bidart indicate that this interpretation is incorrect at a higher time resolution, because in most sections there is no sediment record preserving short-lived assemblages (e.g., Alegret et al., 2015). The K/Pg transition was thought to be complete at Bidart because the lowermost biozones of the Danian are well represented (Haslett, 1994; Apellaniz et al., 1997), but benthic foraminifera inside trace fossils revealed that a probably very thin section of the lowermost Danian clay layer was eroded, as it is only represented by the sediment infilling burrows preserved in the uppermost Cretaceous strata (Alegret et al., 2015). The benthic foraminifera in the burrows mostly consist of buliminids and thus represent a short interval of high productivity in the few thousands of years represented in deposition of the K/Pg clay layer (Alegret et al., 2015) instead of the oligotrophic conditions originally inferred (Alegret et al., 2004a). Whether this earliest Danian “Bulimina interval” was widespread in the Basque-Cantabrian Basin and lost from the stratigraphic record by, e.g., a regional erosional event, we do not know, but so far the specific buliminid taxa observed inside the burrows at Bidart have not been observed in the other sections of this basin, nor even in the Loya section, at a distance of just 19 km from Bidart. Such short intervals of time, possibly less than 1000 years, shortly after a massive impact, may have been unobserved if present only in a few centimeter-length intervals or removed by erosion. Possibly, similar short-term, high-productivity events are reflected by the peaks in the biserial genera Bolivina and Spiroplectammina in a very thin basal interval of the K/Pg clays in Caravaca (Southeastern Spain; Coccioni and Galeotti, 1994). In another example of very localized preservation of sediment deposited very shortly after the extinction, pockets of clayey sediment in the type Maastrichtian area in the Netherlands may have been deposited in the first few hundreds to thousands of years of the Cenozoic (e.g., Henehan et al., 2019, and references therein). In addition, there are surprising records that are argued to reflect conditions in the first hours to one day after the impact in the Chicxulub crater itself, although the precise timing of these events is debatable (Lowery et al., 2018; Gulick et al., 2019; Bralower et al., 2020a; Schaefer et al., 2020).

The completeness of the K/Pg record, at high time resolution, is therefore of vital importance in reconstructing the environmental consequences of the impact on time scales of tens to hundreds to thousands of years or shorter, and sections traditionally considered to be continuous (such as Bidart) due to the presence of all biozones may be incomplete in detail, requiring re-evaluation of paleoenvironmental interpretations. Sections with deposition of impact materials along the eastern U.S. margin, such as ODP Site 1049 (Blake Nose, Northwestern Atlantic; Alegret and Thomas, 2004) and possibly such sites as Bass River (New Jersey, USA, shelf; Olsson et al., 1997; Esmeray-Senlet et al., 2015), may miss the important first few hundreds to thousands of years after the impact so that immediate post-impact events are not recorded. Alternatively, sections thought to be incomplete, such as the type Maastrichtian, may retain a high-resolution record in small areas (Henehan et al., 2019). Ongoing high-resolution studies on the Urrutxua section, using an unprecedented high-resolution sampling in the Basque-Cantabrian Basin, may throw new light on the earliest Danian reconstructions. Overall, whether one correlates biotic or isotopic records, one needs to carefully evaluate evidence for the completeness of sedimentary records at a specific time resolution while working at the limits of resolution possible in many geological sections.

At Walvis Ridge, the K/Pg boundary was recovered at two sites: at Site 1267 (lower bathyal, 2000 m paleodepth), the percentage of infaunal taxa was slightly higher than at Site 1262 (upper abyssal, 2500–3000 m paleodepth; Alegret and Thomas, 2007; Alegret et al., 2012), which was in deeper waters, and thus was expected to be more oligotrophic because of the exponential decrease in food supply with water depth (e.g., Martin et al., 1987). In a deeper setting, the scarce food is taken up mostly by epifaunal morphogroups, so that abundant infauna cannot be sustained. A sharp drop in diversity of the benthic assemblages occurred at the same level as the collapse in the carbon isotope gradient at the K/Pg boundary (Fig. 3), as observed elsewhere. Species extinction rates across the K/Pg boundary were very low at both sites (< 5%) (Table 2). At Site 1267, in contrast with Site 1262, we see that infaunal morphogroups peaked markedly in abundance, as did buliminids, and BFARs were very high in the lowermost Danian (Fig. 3), which indicates enhanced export productivity for the first thousand years after the impact. As suggested by Westerhold et al. (2008), a short unconformity in the lowermost Danian at Site 1262 is probably the cause of this lack of expression of a short-lived event, as supported by the occurrence of a thin, dark red-to-white layer with spherules at Site 1267 and not at 1262 (Fig. 5).

These examples, combined with the benthic foraminiferal literature, highlight several points. (1) It is important to carefully evaluate sedimentary evidence to determine whether K/Pg boundary records are complete and on what time scale (e.g., decades to centuries or hundreds of thousands of years) for the purpose of reconstructing the rapid paleoenvironmental consequences of the asteroid impact at 66 Ma. Misinterpretations of, e.g., export productivity across the K/Pg boundary interval may, in part, arise from the varying and unknown duration of periods from which no sediment is preserved, on time scales of a thousand years or less, as observed in, e.g., Bidart and at Site 1262, in most records from around the Gulf of Mexico and along the NW American margin, such as Blake Nose (Alegret and Thomas, 2005) and in shelf regions in New Jersey, such as Bass River (Olsson et al., 1997). (2) Comparable data sets are needed that are based on quantitative studies of benthic assemblages in the same size fractions, which follow similar taxonomic concepts (even if using informal taxonomy; Arreguín-Rodríguez et al., 2018) that allow for correlation of faunal changes globally (Alegret et al., 2021) and calculation of extinction rates (which must be global, by definition) using the same criteria. (3) Interpretations are needed that are based not only on benthic foraminiferal morphogroups, of which the interpretation is uncertain even for modern assemblages, and this is even more important for extinct morphogroups. Morphogroups are relatively easy to identify, but their use as the only approach to infer conditions at the seafloor may lead to misinterpretations of the fossil record; thus, they must be accompanied by an integrated study of the assemblages (quantitative approach, identification of taxa at the species level whenever and for as many taxa as possible, and BFAR calculations when sedimentation rate estimates make that possible).

Here we address the status of the discussion on oceanic productivity after the end Cretaceous mass extinction. We consider the discussion of oceanic primary and export productivity as an integral part of the discussion on mass extinctions, their severity, and extinction/survival patterns. The publications by Alvarez et al. (1980) and Smit and Hertogen (1980) were evidently among the more transformative and fertile in Earth Sciences: these papers were not the last word on a topic, but they brought mass extinctions to the attention of a broad audience of scientists as well as the general public. They refocused attention on a long-existing basic debate in Earth Sciences (going back to Charles Darwin, Charles Lyell, and Phillips [1860]) about whether mass extinctions are a real phenomenon or reflect basic incompleteness of the geological (and especially fossil) record (e.g., Newell, 1962). These papers were transformative, because they triggered a wealth of fruitful and highly interdisciplinary research persisting until today and predictably into the future that addresses a wide variety of subjects ranging from the physics and chemistry of large impacts (e.g., Morgan et al., 2016) to climate modeling leading to, e.g., the nuclear winter hypothesis (e.g., Turco et al., 1983) and discussions of long-term global warming (e.g., Hull et al., 2020), astronomic hypotheses including those of potential periodicity of impacts (Sepkoski and Raup, 1986; Rampino and Caldeira, 2015), and even the potential existence of a dwarf companion star of our Sun tentatively called Nemesis (Davis et al., 1984). This large literature includes topics relevant to understanding our present and future Earth, such as anthropogenic pollution, global warming, and oceanic deoxygenation and acidification (e.g., Crutzen, 1987; Prinn and Fegley, 1987; Sigurdsson et al., 1992; compilation by Hönisch et al., 2012). Pre-2000 papers on potential acidification effects of an impact, however, did not address ocean acidification, since at the time this environmental effect of CO₂ emissions was little known. They instead looked into the effects of acid rain linked to anthropogenic sulfate and nitrate emissions, which at the time was a major environmental concern (e.g., Likens and Bormann, 1974).

The papers also triggered extensive research into the geological record of fires due to burning of both biomass and sedimentary organic compounds (“fossil fuels”) (Crutzen, 1987; Wolbach et al., 1990; Durda and Kring, 2004; Belcher et al., 2005; Harvey et al., 2008; Morgan et al., 2013), Earth system and ecosystem modeling (e.g., Chiarenza et al., 2020; Gibbs et al., 2020), as well as the detailed evaluation of the fossil record using novel statistical methods and high-resolution stratigraphy. A direct effect of the publication of Alvarez et al. (1980) was the search for and final discovery of the Chicxulub crater (Hildebrand et al., 1991), with a treasure trove of new data generated by drilling of the crater itself (e.g., Urrutia-Fucugauchi et al., 2004; Morgan et al., 2017). As another outcome, the research into the end-Cretaceous impact refocused attention on non-impact causes of mass extinction, such as the eruption of large igneous provinces (with literature ranging from Rampino and Stothers, 1988, to Keller et al., 2020, and Nava et al., 2021), because there is general agreement that the Chicxulub impact occurred during the time interval of eruption of massive lava flows of the Deccan Traps, although there is no agreement on the exact location of the impact relative to the specific lava flows (e.g., Burgess, 2019; Schoene et al., 2021).

Despite these decades of research, we do not have accurate information on the percentage species extinction, which is maybe not surprising because we do not know the number of species living on Earth today within even an order of magnitude (e.g., Mora et al., 2011). We commonly see a citation of a species extinction rate of 75–76%, which does not always refer to a specific source or explain whether this is an estimate for all species or for marine shelled species only. As for our data on benthic foraminifera, species extinction rates in general are extremely difficult to estimate, in part due to severe inconsistencies in taxonomy, and in part due to a lack of time resolution in the age information about species occurrences in data compilations at the genus level (e.g., Raup and Sepkoski, 1982) which average information over a geological stage (Maastrichtian and Danian: ~6.5 and 4.5 m.y.). A source that documents a 75–76% species extinction rate is Jablonski (1994), who uses Raup’s (1979) rarefaction method to extrapolate from Sepkoski’s compilation of a genus extinction rate of 47% (referring to Sepkoski, 1996, which Jablonski [1994] cited as Sepkoski [1994, in press], which took two years to make it through the publication process). Sepkoski’s (1996) compilation refers to marine vertebrates, invertebrates, foraminifera, and radiolarians, and thus it does not include terrestrial organisms (even dinosaurs). Russell (1977, 1979) postulated a species extinction rate of 75% (likewise extrapolated from genus extinction rates) for the end Cretaceous, which included marine and terrestrial species, plants, and animals. Alroy (2008) provides a more recent genus extinction rate of 53%, but only for marine invertebrates.

The publication of the asteroid impact hypothesis (Alvarez et al., 1980) and that of the now widely recognized big five mass extinctions (Raup and Sepkoski, 1982) occurred within a short time. Despite this wide recognition of five Phanerozoic mass extinctions in the popular literature, in which a sixth mass extinction is recognized as being anthropogenic (e.g., Kolbert, 2014), there is no agreement among scientists that there indeed were five major turnover events in the history of life, with some not recognizing the extinction at the end of the Cretaceous as being of the highest importance in the Phanerozoic (e.g., Rojas et al., 2021).

In short, we do not know the K/Pg species extinction rate with any accuracy, although it is obvious that many species of large animals at the top of the food chain, including dinosaurs on land and ammonites and mosasaurs in the oceans, became extinct, which led to a major restructuring of ecosystems. Published extinction rates, however, generally show data over periods of several millions of years (10–20 m.y. in Russell, 1977) and thus cannot differentiate between long-term extinctions over much of the Maastrichtian, which were possibly linked to CO₂ emissions by the Deccan Traps before the impact, and almost instantaneous extinctions right at the boundary. In addition, rapid extinction of key species could trigger more delayed extinctions due to long-term ecosystem interactions. However, an extensive compilation of biotic data for the Maastrichtian (Hull et al., 2020; supplement therein) shows a lack of significant Late Cretaceous extinctions pre-dating the impact, and thus there is no support for the hypothesis that long-term effects of the Deccan Trap lava flows destabilized global ecosystems before the impact.

Following these decades of research, we also have no unanimity about the primary cause of the mass extinctions or the importance of many potential secondary causes. These include primary cooling due to emissions of sulfate exacerbated by secondary cooling due to soot from fires (Tabor et al., 2020; Lyons et al., 2020) and the reverberation of extinctions through ecosystems, by both bottom-up processes such as extinction of photosynthesizers on land and in the oceans because of darkness, and top-down processes, such as the effect on ecosystems of extinction of top-level predators on land and in the oceans. As highlighted above, part of the apparent discrepancies in patterns of extinction observed at different locations may be due to comparison of records that are incomplete on different time scales, such as potential short-term (days to decades) extreme cooling (e.g., Brugger et al., 2017) followed by much longer-term warming due to CO₂ emissions (Pierazzo et al., 1998; Beerling et al., 2002; Hull et al., 2020). These effects are independent of whether sulfates and CO₂ were emitted due to an impact in gypsum-bearing limestones or by large lava flows. We argue that ocean acidification, a separate process than acid rain (D’Hondt et al., 1994; Alegret et al., 2012; Henehan et al., 2019), may have played an important role in the mass extinctions, at least for shallow and surface water-dwelling calcifying taxa, because of the severe selectivity in extinction of pelagic calcifiers (calcareous nannoplankton and foraminifera) as compared to non-calcifiers (e.g., Lowery et al., 2020, and references therein). Ocean acidification after the K/Pg impact probably was due to the effects of SO₂ rather than CO₂ emission (D’Hondt et al., 1994; Alegret et al., 2012; Tyrrell et al., 2015). Against acidification as an important contributing cause, it has been argued (Lowery et al., 2020) that acidification of surface waters at the K/Pg boundary (Henehan et al., 2019) was similar in magnitude (e.g., several tens of pH units) to that of the Paleocene–Eocene Thermal Maximum (Penman et al., 2014; Gutjahr et al., 2017), at which time pelagic calcifiers did not suffer extinction, but we think that this is not a valid argument. The rate of acidification after the impact was practically instantaneous, whereas that during the Paleocene–Eocene Thermal Maximum probably stretched over several thousands of years (Zeebe et al., 2014, 2016), and the severity of environmental effects of acidification at the surface of the ocean and thus the effect on biota depends strongly on that rate (e.g., Hönisch et al., 2012; Lord et al., 2016). Slower rates will lead to more severe effects for deep-sea biota, as observed during the Paleocene–Eocene Thermal Maximum (e.g., Thomas and Shackleton, 1996), and less severe effects for surface biota, as observed during the K/Pg extinction (Alegret et al., 2012).

As to oceanic productivity, publications in the 1980s argued that post-impact oceans were to a large extent deprived of “life,” the “Strangelove Ocean,” though this term was undoubtedly proposed tongue-in cheek (Hsü and McKenzie, 1985). We speculate that the hypothesis that the oceans were almost lifeless was inspired first by the selectively extreme and rapid extinctions of calcifying pelagic autotrophs, mixotrophs, and heterotrophs (calcareous nannoplankton and planktic foraminifera; compiled in Lowery et al., 2020; see also Gibbs et al., 2020) recognized by Alvarez et al. (1980) as coeval with their evidence for an impact, and second, by the collapse of the surface-deep-sea δ¹³C gradient (e.g., Hsü et al., 1982; Hsü and McKenzie, 1985; Zachos and Arthur, 1986). This collapse was originally interpreted as an almost complete cessation of photosynthesis in the oceans (e.g., Hsü et al., 1982; Hsü and McKenzie, 1985; Zachos and Arthur, 1986; Kump, 1991). Such complete cessation of photosynthesis for the duration of the gradient collapse (a few hundred thousand years) disagrees with the observed lack of significant extinction of deep-sea benthic foraminifera (e.g., Alegret et al., 2012), which depend on surface water photosynthesis for survival through bentho-pelagic coupling (e.g., Gooday, 2003; Jorissen et al., 2007). More than a decade later, in the “Living Ocean Model,” the δ¹³C gradient collapse was re-interpreted (D’Hondt et al., 1998) as due to a decrease in organic flux reaching the seafloor rather than a collapse of primary productivity, with the decrease in flux being due to a combination of the extinction of fecal pellet-producing marine animals as well as a decline in the size of primary producers. These authors explained the “low post-extinction δ¹³C gradient by just slightly increasing (from 90 to 95%) the fraction of total production that was degraded in the upper 200 m of the ocean.” Such a relatively minor decline in the strength of the biological pump would be in line with the lack of extinction of benthic foraminifera, but later authors interpreted the “Living Ocean Model” as a much more severe collapse of the biological pump (e.g., Coxall et al., 2006), which would not agree. Henehan et al. (2019) used Earth system modeling to arrive at an estimate of a 50% reduction of globally averaged export production, but show modeled geographic variability in export production.

These interpretations do not fully address the fact that the δ¹³C record in surface- and deep-water biogenic carbonates cannot be simply explained as the record of δ¹³C in DIC, because of the extreme extinction of pelagic calcifiers that carry the surface δ¹³C signal. The deep-water record is based on benthic foraminifera, which did not suffer extinction, so that one can analyze one species that crosses the K/Pg boundary (see e.g., Alegret et al., 2012; Henehan et al., 2019; data on Site 1267 presented here). However, one cannot do that for the surface record: Maastrichtian planktonic foraminifera were relatively large, with a positive δ¹³C signature due to the presence of symbionts, and were replaced by small, non-symbiont–carrying Danian forms with a lighter isotopic signal (e.g., Coxall et al., 2006; D’Hondt, 2005; Birch et al., 2016; Shaw et al., 2021). The Maastrichtian fine carbonate fraction (thus bulk carbonate) was dominated by photosynthesizing calcareous nannoplankton, but in the Danian the fine carbonate fraction was in many places dominated by fragmented calcareous dinoflagellate cysts, which are isotopically light (Zonneveld, 2004; Zonneveld et al., 2007; Alegret et al., 2012; Minoletti et al., 2014), mixotroph and non-photosynthesizing calcareous nannoplankton (Gibbs et al., 2020), and microbially produced micrite (Minoletti et al., 2005; Bralower et al., 2020b). A significant part of the δ¹³C vertical gradient collapse, though probably not all (Birch et al., 2016), thus was a consequence of the extinction, and does not reflect the δ¹³C of DIC in the surface waters but instead reflects a change in the carriers of the δ¹³C signal.

Information from geochemical productivity proxies (Hull and Norris, 2011), organic biomarkers for bacterial productivity (Sepúlveda et al., 2009, 2019), records of non-calcareous phytoplankton such as dinoflagellates (Brinkhuis et al., 1998; Hildebrand-Habel and Streng, 2003) and diatoms (Harwood, 1988; Renaudie et al., 2018), and micrite-producing bacteria (Bralower et al., 2020b) all indicate that post-extinction oceans, from the direct time after impact through possibly a few to tens of thousands of years, were characterized not by a lifeless state, but by spatially and temporally varying, possibly extreme, and in part prokaryote (cyanobacterial) blooms of primary producers, in what has been called the “Heterogenous Ocean” model (Hull and Norris, 2011; Alegret et al., 2012; Esmeray-Senlet et al., 2015), as also modeled by Henehan et al., 2019 (their fig. S24). The varying taxonomic composition of the primary producers from site to site is evidenced by changes in the morphology of microscopic carbonate grains (Bralower et al., 2020b) and at some sites by very low δ¹³C values that indicate a strong contribution from isotopically light calcium-carbonate dinocysts (e.g., Alegret et al., 2012), but at most sites the composition of the fine fraction is not well documented.

Such blooms might have occurred even within the impact crater itself (Bralower et al., 2020a, 2020b). Blooms may have been driven by multiple factors, including the extinction of competing taxa, which might have left large supplies of nutrients available (Henehan et al., 2019). Alternatively or additionally, excess nutrients may have been supplied by nitrogen-oxide input as the result of atmospheric heating by the impact (e.g., Parkos et al., 2015) or by excess erosion and weathering input from land left bare after destruction of vegetation by extensive fires (e.g., Crutzen, 1987; Wolbach et al., 1990; Durda and Kring, 2004; Morgan et al., 2013; Scott et al., 2014). Plankton blooms could have contributed to the extinction of marine animals, if they included toxin-producing algae such as cyanobacteria and dinoflagellates (Castle and Rodgers, 2009).

In contrast to agglutinated and calcareous benthic foraminifera, calcareous pelagic organisms were selectively hard hit by the extinctions. Few (two or three) planktic foraminiferal species may have survived, and Cenozoic taxa may have evolved from these few coastal-dwelling species (Olsson and Liu, 1993; Olsson et al., 1999). Alternatively, Cenozoic planktic foraminifera may have newly evolved from benthic foraminifera (Brinkhuis and Zachariasse, 1988; Arenillas and Arz, 2017) and/or through regional evolution of planktic from benthic taxa at high latitudes (Huber et al., 2020). The severely affected pelagic calcareous nannofossils (e.g., Romein and Smit, 1981; Bown, 2005a, 2005b) may have survived and reradiated in the Cenozoic due to the existence of mixotrophic species that could survive heterotrophically rather than by photosynthesis (Gibbs et al., 2020).

Our new compilation and additional data presented in this paper show broad support for the hypothesis of a “Heterogeneous Ocean” after the K/Pg boundary impact, in which oceanic primary productivity did not decline globally for thousands of years to even hundreds of thousands of years, as had been interpreted from the carbon isotope records. Locally to regionally, benthic foraminifera may have become temporarily extinct due to a lack of food, similar to local extinctions after volcanic eruptions (Hess and Kuhnt, 1996), but on a global scale there must have been sufficient refugia for benthic foraminifera to repopulate through propagules. Locally to regionally, benthic foraminifera qualitatively indicate an increased food flux to the seafloor (Figs. 3–4), as also indicated by geochemical proxies (Hull and Norris, 2011) and by modeling (Henehan et al., 2019). Different benthic foraminiferal-based qualitative and semiquantitative proxies used to estimate organic flux to the seafloor (BFAR, buliminids, infaunal taxa) broadly agree, but there are major differences in detail, and more high-resolution studies are needed to figure out to what extent such proxies indicate food flux overall, or aspects of the food flux, such as seasonality, nature of the organic matter (more or less refractionary, potential lateral transport of continent-derived organic matter along continental margins), or type of phyto-zooplankton.

We do not see a correlation between the severity of extinction (species extinction percent) in deep-sea benthic foraminifera, and factors such as latitude, water depth, and distance from the impact crater. Such a correlation would, however, be difficult to prove: only a few of the many hundreds of species of benthic foraminifera became extinct, and the percent of species extinction thus depends more on the initial species richness in the area of consideration than on the number of species becoming extinct. For example, more than 90% of nannoplankton species became extinct globally, but extinction rates were lower (~75%) at high southern latitudes (Jiang et al., 2010). These high-latitude assemblages had lower diversity before the extinction and contained the survivor species, which assumedly were adapted to more unfavorable (e.g., low light) environmental conditions, so that a higher regional survivorship is expected.

In contrast, changes in benthic foraminiferal diversity are more robust; likewise, they were not linked to any of these factors. This lack of correlation, however, might be due to the lack of time resolution in local or regional records, such as in the areas around the Gulf of Mexico, where the sediment record is clearly and severely affected by the consequences of the Chicxulub impact. More detailed research is necessary to evaluate whether short term events are preserved in geological records (on time scales of less than 10–20 k.y.) at specific sites. Extinction is forever; thus, it is observed even in highly incomplete records, but the record of a short-term (even instantaneous) environmental change causing the extinction (flash acidification, severe cooling, and darkness) may be missing at specific sites due to the incompleteness of the record at such short time scales.

Walter Alvarez brought mass extinctions to the attention of the general public more than four decades ago. Since then, the scientific community has produced and will continue to produce thousands of multidisciplinary works inspired by his seminal works. An integral part of the discussion about the K/Pg mass extinction is the analysis of oceanic primary and export productivity, which was triggered by the impact of an asteroid on Earth at 66 Ma. Our compilation of published and new benthic foraminiferal data supports the hypothesis of a heterogeneous ocean, in which neither oceanic primary productivity nor export productivity declined severely and globally for thousands to a few millions of years after the K/Pg impact. Proxies for export productivity show strong temporal and geographic variability, possibly through extensive plankton blooms after the impact.

Neither the global decrease in diversity nor species extinction rates of benthic foraminifera based on comparable data sets, with normalized taxonomic concepts and methodology, show significant correlation with paleodepth, distance from the Chicxulub impact crater, or paleolatitude. The observed temporal and geographic variability, however, may be related, at least in part, to geographically and bathymetrically variable incompleteness of the record at high time resolution. Careful evaluation of the completeness of the K/Pg records on short timescales is highly recommended for future studies.

We thank Walter Alvarez for inspiring us and many other scientists, and Philippe Claeys, Christian Koeberl, and Sandro Montanari for their invitation to contribute to this volume. This work was supported by the Spanish Ministry of Economy and Competitiveness and FEDER funds (PID2019-105537RB-I00), the Government of Aragon (Reference Group E33_17R, co-financed with FEDER 2014-2020 “Building Europe from Aragon”), and partially by the National Science Foundation (OCE 1536611). This study includes data from samples provided by the Integrated Ocean Drilling Program (IODP).