Biogeochemical significance of pelagic ecosystem function: an end-Cretaceous case study

Pelagic ecosystem function is integral to global biogeochemical cycling, and plays a major role in modulating atmospheric CO2 concentrations (pCO2). Uncertainty as to the effects of human activities on marine ecosystem function hinders projection of future atmospheric pCO2. To this end, events in the geological past can provide informative case studies in the response of ecosystem function to environmental and ecological changes. Around the Cretaceous–Palaeogene (K–Pg) boundary, two such events occurred: Deccan large igneous province (LIP) eruptions and massive bolide impact at the Yucatan Peninsula. Both perturbed the environment, but only the impact coincided with marine mass extinction. As such, we use these events to directly contrast the response of marine biogeochemical cycling to environmental perturbation with and without changes in global species richness. We measure this biogeochemical response using records of deep-sea carbonate preservation. We find that Late Cretaceous Deccan volcanism prompted transient deep-sea carbonate dissolution of a larger magnitude and timescale than predicted by geochemical models. Even so, the effect of volcanism on carbonate preservation was slight compared with bolide impact. Empirical records and geochemical models support a pronounced increase in carbonate saturation state for more than 500 000 years following the mass extinction of pelagic carbonate producers at the K–Pg boundary. These examples highlight the importance of pelagic ecosystems in moderating climate and ocean chemistry.


Introduction
Atmospheric CO 2 concentrations ( pCO 2 ) are regulated by a complex, interconnected system of sources and sinks, both abiotic and biotic [1][2][3]. Biological activity in the surface oceans plays a major role in this via the 'biological carbon pump', whereby pelagic organisms take up carbon in the surface ocean, die and sink, sequestering carbon in the deep ocean. In addition, pelagic calcifying organisms (such as coccolithophores and planktonic foraminifera) export CaCO 3 to the deep oceans, sequestering weathering products from land in sediments (the 'alkalinity pump'). This balances alkalinity fluxes, provides a dissolvable carbonate reservoir that buffers the ocean from potentially harmful pH change, and helps to maintain largely equable climates [4,5]. Together, planktonic foraminifera and coccolithophores account for the vast majority of the pelagic carbonate flux [6], which in turn accounts for almost half of total marine carbonate production [7]. As such, pelagic organisms play an important role in biogeochemical cycling and climate regulation.
Human activities (examples among many include injection of CO 2 , overfishing, oxygen depletion and habitat destruction) threaten the function of the pelagic ecosystem [8][9][10], adding uncertainty to the projection of pCO 2 and climate over the coming centuries [3,11]. In part, this uncertainty stems from a lack of available ecological datasets across the spatial and temporal scales that would be relevant in constraining model predictions [12]. Most existing ecological time series are too short to discern trends beyond decadal variation in the climate system and relatively few studies have addressed the link between biodiversity and ecosystem function on geological timescales (see [13] in this issue for an exception). The microfossil record can be a useful resource in addressing these knowledge gaps [14], and placing constraints on the response of the pelagic ecosystem to environmental perturbations and their effect on biogeochemical cycles. During an approximately 1 million year (Myr) interval surrounding the Cretaceous -Palaeogene (K -Pg) boundary 66.04 Myr ago, two very different disturbances are recorded in the marine fossil record. These events provide case studies on the interplay between environmental change, biodiversity and ecosystem function under similar background conditions. Our study thereby begins to address a gap in our current understanding of the relationship between biodiversity and ecosystem function [13,15].
The onset of vast flood basalt volcanism (the Deccan Large Igneous Province, LIP) in the latest Cretaceous resulted in the release of 15 000-35 000 Gt CO 2 and 6400-17 000 Gt SO 2 over a relatively long (more than 100 000 year) timescale [16,17]. In contrast, the impact of an approximately 10 km wide bolide at Chicxulub at the K-Pg boundary [18] led to instantaneous release of SO x , NO x and CO 2 ( [19] and references within), and rapid and transient (probably less than 5 year) acidification of the surface ocean [20,21]. Besides the very different timescales of these environmental perturbations, a critical difference is that the Chicxulub impact coincides with a major mass extinction and Late Cretaceous Deccan volcanism does not [19]. Species loss in the open ocean following the bolide impact, while variable between groups [22], was particularly high in the calcareous plankton (approx. 95% and 90% in planktonic foraminifera and calcareous nannofossils, respectively [23,24]). In contrast, during Late Cretaceous Deccan trap volcanism, biotic disturbance in the open ocean was largely limited to changes in biogeographic ranges [25,26]. Together, these events allow us to contrast the impact of environmental changes on ecosystem function with and without associated loss of pelagic biodiversity. Here, we use carbonate preservation indices to gain a fuller understanding of changes in biogeochemical ecosystem function across this interval, combining new and previously published records of carbonate preservation from geographically disparate deepsea sites with new insights from ocean carbon cycle modelling.

Methods: carbonate preservation
Change in deep-ocean carbonate saturation state (V CaCO3 ) is an indicator of broader carbon cycle disturbance that can be readily discerned in the geological record using records of deep-sea carbonate preservation [27]. New and previously published records of a number of different CaCO 3 preservational indices are compiled here from globally distributed deep-sea drill core sediments over a 3.7 Myr interval surrounding the K -Pg boundary, 66.04 million years ago (Ma). Each CaCO 3 preservation metric has associated strengths and limitations, which we discuss at length in the electronic supplementary material. Where possible, our new records of deep-sea preservation use counts of planktonic foraminiferal fragmentation (as in [28,29]). This metric relies on the observation that with decreasing deep-ocean V CaCO3 microfossils progressively dissolve and fragment [30] (see electronic supplementary material, figure S1). New fragmentation data were generated from Shatsky Rise in the Pacific (ocean drilling programme (ODP) site 1209) and Walvis Ridge in the South Atlantic (ODP site 1267). Meaningful fragmentation counts from the Newfoundland Sediment Drifts site in the North Atlantic (International Ocean Drilling Programme (IODP) site U1403) were not attainable owing to extensive dissolution prior to the K -Pg boundary (see electronic supplementary material, figure S2). At this site, weight per cent (wt.%) coarse fraction (greater than 38 mm) was used as a carbonate preservation indicator (though important caveats to this production-sensitive metric are discussed in the electronic supplementary material). Sediment samples were dried and weighed before being disaggregated in de-ionized water on an orbital shaker and washed through a 63 mm (Walvis Ridge, site 1267 and Shatsky Rise, site 1209) or 38 mm (Newfoundland, site U1403) sieve with de-ionized water. Both the greater-than-63 mm/greater-than-38 mm coarse fraction and the fine fraction were then dried at approximately 458C and the coarse fraction weighed to calculate wt.% coarse fraction. For Walvis Ridge (site 1267) and Shatsky Rise (site 1209), the relative abundance of 'complete' tests (i.e. whole tests that show no signs of any breakage or dissolution of chambers) was counted from a representative split (200-400 fossils) of the greater-than-125 mm size fraction. Full details of the methods used to construct age models for each site (including the construction of new age models for previously published data) are given in the electronic supplementary material.

Methods: carbon cycle modelling
The geochemical box model Long-term Ocean Sediment CArbon Reservoir (LOSCAR) v. 2.0.4 [31] was employed to simulate the impacts of volcanic degassing and calcifier extinction on the global carbon cycle, with some modifications. Importantly, to better account for the very different [Ca 2þ ] and [Mg 2þ ] in the K-Pg ocean [32], updated carbonate chemical equilibrium constants from the MyAMI model [33] were substituted into the model, using a [Ca 2þ ] of 42 and [Mg 2þ ] of 20 mmol kg 21 . All plotted model runs (figure 3 and electronic supplementary material, S5 -11) were initiated at a steady state pCO 2 of 600 ppm (in agreement with palaeosol carbonate measurements ( [34,35], and references therein)), and assume a climate sensitivity of 38C per doubling of pCO 2 . This climate sensitivity is in the middle of the range (2.2-4.88C) of observed climate sensitivity over the past 65 Myr [36]. However, a range of other starting atmospheric CO 2 concentrations (400-1000 ppm) and climate sensitivities (0 -58C per doubling) were also explored, with results listed in electronic supplementary material, table S1 (see also electronic supplementary material, Discussion). Our primary experiments (figure 3) also assume a stronger-than-modern silicate weathering feedback to account for a greater abundance of exposed fresh Deccan basalt at low latitudes (see electronic supplementary material for more details), although model runs at a range of feedback rstb.royalsocietypublishing.org Phil. Trans. R. Soc. B 371: 20150510 strengths were also tested (see electronic supplementary material, table S1 and figure S10).
For simulations of Deccan degassing, minimum and maximum emission scenarios (total CO 2 ¼ 4090 or 9500 Gt C; total SO 2 ¼ 3200 or 8500 Gt S [16]) were partitioned into two discrete pulses, in accordance with the proposed second and third stages of volcanism from [37], and the eruptive volumes of [38]. 86.5% of degassing was input over an approximately 140 Kyr interval beginning at the C30n/C29r magnetochron reversal, approximately 360 Kyr before the K -Pg boundary. This corresponds to an observed interval of decreasing seawater 187 Os/ 188 Os (which indicates elevated basalt weathering) [39]. The remaining 13.5% of the volcanic emissions was then released in models at the end of magnetochron C29r in the Danian (250 kyr after the K -Pg boundary). Most other estimates of CO 2 release for the Deccan traps [38,[40][41][42] fall within the range of emissions tested here [16]. To better discern the effects of each gas, scenarios for CO 2 and SO 2 release were also tested in isolation (figure 3). As in [21], SO 2 release and rain-out were simulated by reducing alkalinity in the surface ocean box (see electronic supplementary material for more details). A wide range of possible timescales and modes of degassing were also tested (see electronic supplementary material, table S1).
For simulations of the biogeochemical consequences of the K -Pg mass extinction, a range of carbonate flux reductions were tested, ranging from 10% up to 75%. We tested two types of scenarios: (i) reductions in CaCO 3 flux with no change in the organic carbon flux (i.e. a change in the C CaCO3 : C org flux ratio); and (ii) reductions in overall efficiency of the biological carbon and alkalinity pump (i.e. reducing both C CaCO3 and C org fluxes). For each simulation, changes in fluxes were imposed for 200 Kyr following the K -Pg boundary and then tapered back to pre-event values over a further 200 Kyr to simulate the gradual recovery of early Palaeocene pelagic ecosystems. For further details and discussion about modelling approaches, see electronic supplementary material.

Deccan volcanism, global warming and carbonate dissolution
The main phase of Deccan volcanism [37] is recorded in deepsea sediments by a global decline in 187 Os/ 188 Os [39] just after the C30n/C29r magnetochron reversal [17] at 66.398 Ma (figure 1a). The onset of volcanism and associated release of CO 2 coincides with evidence for a transient warming event (figure 1b) in both geochemical [43][44][45][46] and palaeoecological data [25,52,53]. Our data show a pronounced increase in deep-sea carbonate dissolution in several ocean basins at this time, in response to this volcanism (figure 1c-h). Dissolution is particularly pronounced in the Southern Ocean (ODP site 690 [39,47,48] and figure 1c) and North Atlantic (IODP site U1403 [49] and figure 1d), with wt.% carbonate falling by approximately 20% and approximately 40%, respectively. This result is consistent with enhanced dissolution in high-latitude sediments closest to sites of deep water formation [54], where the impact of increased CO 2 emissions will first be felt. At lower latitudes, increased foraminiferal fragmentation seen at Walvis Ridge (ODP site 1267, figure 1e (this study); DSDP site 527, figure 1f [50]) and Shatsky Rise (ODP site 1209, figure 1g (this study)) is indicative of a shoaling of the lysocline (i.e. the depth at which substantial carbonate dissolution occurs). Reduced planktonic foraminiferal preservation elsewhere on Shatsky Rise (DSDP site 577 [51]; figure 1h), and selective preservation of dissolution-resistant coccolithophores in the Indian Ocean [47] (electronic supplementary material, figure S4) corroborate this observation. In all cases, records of increased dissolution return to roughly pre-event values before the K -Pg boundary (figure 1), restricting the main degassing phase of Deccan volcanism to a distinct less than 200 Kyr interval beginning at the onset of magnetochron C29r, around 66.398 Ma. This supports previous inferences for only transient ocean acidification based on Ir accumulation [39], and suggests Deccan degassing played no direct role in K-Pg mass extinction. New absolute age constraints for the Deccan eruptions [17] have been cited as evidence of a Deccan role in the K-Pg extinction through ocean acidification [55]. Our data (and modelling below and electronic supplementary material, figures S5 and S6) suggest that even these new timescales for eruption are still long enough for surface ocean carbonate saturation to be maintained via carbonate compensation and silicate weathering (see also [56]).

Bolide impact and mass extinction at the K -Pg
In the aftermath of the K-Pg, sediment records from the Pacific and Atlantic ( figure 2) show a pronounced rise in wt.% coarse fraction as a result of both decreased calcareous plankton production and enhanced foraminiferal preservation. Simultaneously, fragmentation of planktonic foraminifera at both Walvis Ridge [59] and Shatsky Rise [28] declines (figure 2), even to essentially no fragmentation at Shatsky Rise. Because some foraminiferal fragmentation is normally expected during sinking and sedimentation even above the lysocline [62,63] (see electronic supplementary material, figure S1), this lack of discernible fragmentation at Shatsky Rise indicates very high [CO 3 22 ] throughout the water column. Rapid and pronounced deepening of the lysocline owing to this enhanced [CO 3 22 ] is evidenced at the Newfoundland Sediment Drift in the North Atlantic (IODP site U1403), where we observed a step-change across the K-Pg from Maastrichtian sediments barren of any planktonic foraminifera to post-boundary sediments in which Danian planktonic foraminiferal species are excellently preserved (electronic supplementary material, figure S2) up until around magnetochron C28r ( figure 2). Similarly, in the South Pacific (IODP site U1370, 5076 m depth) [64], the only carbonate preserved over the last 75 Myr is in the immediate aftermath of the K-Pg boundary, within nannofossil zones NP1 and NP2 [64]. Elsewhere, at the Ontong-Java plateau (ODP site 803, 3410 m depth), carbonate is preserved for a brief interval (less than 1 m, within biozone NP1) around the K-Pg boundary, but is absent above and below [65]. Additional lines of evidence for a rise in oceanic V CaCO3 are also discussed in the electronic supplement material ( §6d). These lines of sedimentological evidence all support the predictions of earlier work [66,67] that reduced pelagic carbonate production owing to extinction of calcareous plankton following the K-Pg bolide impact [68,69] profoundly impaired the marine alkalinity pump (a key pelagic ecosystem function) and prompted a rstb.royalsocietypublishing.org Phil. Trans. R. Soc. B 371: 20150510 period of alkalinity build-up, deepening of the lysocline, and ocean pH rise.

Comparison with carbonate system models
Using the LOSCAR carbon cycle model [31], we attempt to reproduce observed patterns of environmental change and deep-sea carbonate preservation. For pre-boundary volcanism, we find that only high-end Deccan CO 2 emission scenarios can produce the widely observed late Maastrichtian warming of approximately 2-38C at mid-range climate sensitivity (38C/CO 2 doubling), for an initial atmospheric pCO 2 of 600 ppm [34] and an eruptive duration of 140 Kyr. Moreover, with this forcing, only high-strength silicate weathering feedbacks (see electronic supplementary material, Discussion) could draw down CO 2 and temperature within only several hundred thousand years, consistent with observations (figure 1b). For lower CO 2 emission scenarios, either highend late Cretaceous climate sensitivity (greater than 38C per CO 2 doubling) or lower initial pCO 2 are required to produce observed warming (see electronic supplementary material, table S1).
In terms of carbonate cycle perturbation, LOSCAR predicts at most only fleeting reductions (less than 40 figure S5). For even the largest estimates of SO 2 and CO 2 release, LOSCAR suggests that eruptive timescales of less than 100 Kyr are required to produce pronounced lysocline shoaling (electronic supplementary material, figure S6), and even then this shoaling would be briefer (less than 50 Kyr) than indicated in the sedimentary record (approx. 150 -200 Kyr, figure 1c-h). Instead, the dominant long-term signal predicted for Deccan CO 2 release under any modelled emissions scenario is elevated weathering fluxes, a rise in oceanic carbonate saturation, and a deepening of the lysocline (figure 3). We observe little evidence for this enhanced preservation (or 'carbonate overshoot') following the initial dissolution pulse of Deccan volcanism (figure 1c-h). The brevity of dissolution relative to preservation records and the existence of a pronounced carbonate overshoot are consistent in all modelled scenarios, despite different timescales for release, total emissions, starting pCO 2 , equilibrium constants and weathering feedbacks (see electronic supplementary material, Discussion and table S1).
There are multiple possible explanations for this mismatch between empirical observations and model predictions (discussed in depth in the electronic supplementary material), including an overestimation of the duration of the Cretaceous portion of magnetochron C29r (as suggested by recent U-Pb dating; [17]), changes in circulation or productivity, elevated CaCO 3 deposition in shelf settings (see electronic supplementary material, figure S9) or the influence of processes not accounted for in LOSCAR. Another possible explanation is that Deccan-induced warming resulted in a more stratified ocean with more oligotrophic surface waters [26,52,53]. In the modern ocean, oligotrophy favours ecosystems more heavily dominated by coccolithophore production when compared with siliceous and organic-walled primary producers [70]. If this was similar in the Cretaceous ocean, and Deccan warming did indeed result in enhanced stratification and more oligotrophic oceans, it is possible that CaCO 3 production and export rose. A modelled increase in CaCO 3 : C org ratio of 30% during simulated warming succeeds in extending the timescales of deep-ocean carbonate dissolution to approximate agreement with sedimentary records, amplifying atmospheric CO 2 rise, and dampening subsequent carbonate saturation increase (electronic supplementary material, figure S8). This emphasizes the potential importance of accounting for biotic, ecological feedbacks when considering the ocean's response to greenhouse gas forcings.
We also simulate the effects of an extinction of pelagic carbonate producers at the K -Pg boundary (figure 3 and electronic supplementary material, figures S5 and S11). Although the bolide impact [20,21] and a brief reduction in photosynthetic carbon uptake [59] could have induced acidification of surface waters and released CO 2 from the oceans on timescales of less than 10 Kyr [67] (figure 3), the more significant long-term impact on the carbon cycle comes about from the major extinction in both main groups of pelagic calcifiers. This extinction, and loss of abundance, caused changes in carbonate saturation state that persisted for more than 1 Myr (figure 2). We demonstrate that even a conservative 30% reduction of CaCO 3 export flux results in a deepening of the Atlantic carbonate compensation depth (CCD) by 2 km, an increase in surface V Calcite from 6.6 to 10 (electronic supplementary material, figure S11) and a drop in atmospheric CO 2 of approximately 100 ppm (figure 3), consistent with modelled findings of earlier studies [66,67]. This elevation of ocean alkalinity in response to mass extinction could provide a mechanism for low atmospheric pCO 2 estimated for the early Danian [71]. Our modelling suggests a 30% drop in CaCO 3 export would also lower the modelled Pacific CCD  to approximately 4200 m-enough to bring the CCD below the South Pacific Gyre IODP site U1370 but above site U1365, consistent with sedimentary observations [64]. Deepsea sediment cores, though, suggest a much greater reduction in pelagic CaCO 3 production and delivery ( [59] and references within). In our model runs, a reduction in CaCO 3 production by more than 50% would produce sufficiently high supersaturation to initiate abiotic precipitation of CaCO 3 in surface waters (electronic supplementary material, figure S10)-a process which today is not known beyond tropical shelf settings such as the Bahamas or Persian Gulf. While there is perhaps some evidence for this [72], it is also possible that other synchronous changes may have occurred to avoid critical supersaturation. Increased burial of carbonate on shelves to compensate for less deep-ocean burial could have played a role [67] (see also electronic supplementary material, figure S9), although evidence for such an increase is, at best, scant (see electronic supplementary material, Discussion).

Volcanism, impact and the carbon cycle: implications for biodiversity and ecosystem function
Environmental forcing imposed by Deccan emplacement and K -Pg bolide impact produced very different recorded changes in ecosystem function (figures 1 and 2), primarily as a result of very different patterns of ecological response. There is little evidence for loss of species or population abundance in the open ocean plankton during Late Cretaceous Deccan volcanism. The approximately 2-38C warming associated with Deccan CO 2 release resulted in range expansions [25,26,52], dwarfing of some planktonic foraminiferal species [73] and regional assemblage changes [74], but there was no elevation in extinction rates of functionally important marine calcifier species (planktonic foraminifera and coccolithophores) at this time [69,75]. This retention of biodiversity and redundancy among calcifiers, we suggest, was probably important in maintaining the resilience of the pelagic ecosystem (and its associated biogeochemical functions) [15,76]. Consequently, the marine carbonate cycle, coupled with global silicate weathering feedbacks [1], could assimilate Deccan-derived CO 2 over these timescales without drastic, long-lasting effects on surface ocean V CaCO3 (only very modest lysocline shallowing and some reduction in surface ocean pH-see electronic supplementary material, figure  S5-are indicated). This role of pelagic calcifiers in mitigating the impact of CO 2 emissions is underscored by considering similar volcanic episodes before the evolution of pelagic calcifiers [4,5]. Two of these earlier episodes, the end-Triassic Central Atlantic Magmatic Province [77] and Permo-Triassic Siberian Traps [78] volcanism, had profound environmental impacts and resulted in two of the largest mass extinctions in the history of life [5]. The more profound and long-lasting perturbation of surface ocean carbonate saturation we observe over the K -Pg transition arises from mass extinction following the Chicxulub bolide impact. The near-complete loss of the clades responsible for the vast majority of pelagic carbonate cycling ( planktonic foraminifera and coccolithophores [68,69]) resulted in a build-up of alkalinity in the Earth's ocean (as evidenced by improved deep-ocean carbonate preservation; figure 2). This, in turn, may have also drawn down atmospheric CO 2 and prompted climatic changes ( [66], this study). While some evidence suggests export of organic carbon to the deep ocean had largely recovered within a few hundred thousand years [79,80], carbonate preservation (figure 2) suggests recovery of full pre-event biogeochemical function in pelagic ecosystems took more than a million years, coinciding with restoration of micro-and nannofossil biodiversity [81]-an example of the close link between  Figure 3. LOSCAR [31] simulations of Deccan degassing release and K-Pg reduction in CaCO 3 flux. Simulations of atmospheric pCO 2 (a) and deep-water calcite saturation state (b) response to simulated perturbations. Maximum and minimum CO 2 and SO 2 efflux estimates for Deccan eruptions are from ref. [16], partitioned into a main eruptive phase ('phase II') beginning at the onset of C29r (86.5%) and a later one (13.5%) at the close of C29r in the Danian ('phase III') [39]. CaCO 3 flux changes impact the CaCO 3 : C org ratio assuming no change in organic flux. Reduction in biological pump efficiency reduces the efficiency of the biological pump in using the parametrized nutrient pool (see electronic supplementary material, Discussion for more details). Model outputs are plotted against simulated model years, relative to the simulated K-Pg boundary. rstb.royalsocietypublishing.org Phil. Trans. R. Soc. B 371: 20150510 biosphere and geosphere dynamics in the aftermath of mass extinction (see also [82]). Current global change is altering pelagic ecosystems, but the extent of this alteration in biodiversity [9] and ecosystem structure [83] and its ultimate biogeochemical significance, remains unclear [2,3]. In the case of the K -Pg boundary, the extinction was particularly selective against pelagic calcifiers, and post-extinction ecosystems lacked both the diversity and abundance of pre-extinction oceans. Although it is the decline in calcifier abundance that directly accounts for the decline in ecosystem function, it remains an open question how important standing richness, within and across calcifier clades, is in determining calcifier abundance across the event. It is noteworthy in this context that post-extinction biogeochemical function (and by inference the abundance of calcifiers) recovers long in advance of the full recovery of pre-event levels of calcifier diversity (figure 2 and also [58]). This observation suggests that while functional redundancy among latest Cretaceous calcareous plankton may have helped to confer resilience on carbonate export [76] in the face of volcanic CO 2 and SO x emissions and global warming, a much lower standing diversity can still support a comparable carbonate alkalinity pump. For the oceans today, it is crucial to determine where tipping points may lie with regards to shifting the abundance of marine organisms, as it is the aggregate effect of many, many billions that account for pelagic ecosystem function. As we show here, pelagic ecosystem change, particularly in pelagic calcifiers, can profoundly influence the long-term evolution of the Earth system. Data accessibility. The datasets supporting this article have been uploaded as part of the electronic supplementary material (tables S2 and S3), and will in due course be made available via www.Pangaea.de.