Microbial metabolism amplified warming in three Phanerozoic hyperthermal events

Domain: earth_science
Published: 2025-12-12
Editorial outcome: major_revision

Abstract
The ocean’s biological carbon pump (BCP) plays a crucial role in climate regulation by facilitating long-term removal of atmospheric CO2. However, BCP changes during past hyperthermals remain poorly understood. Here we use the vertical carbon isotope (δ13C) gradients between surface and mid-waters to trace the BCP changes during the Permian−Triassic, early Toarcian, and Palaeocene−Eocene hyperthermals, detecting a two- to three-fold increase in the vertical δ13C gradients. We found that enhanced organic matter remineralization, driven by increased microbial respiration under warming, are the primary trigger for increased vertical δ13C gradients. Warming-induced enhancements in microbial metabolism caused a two- to six-fold reduction in the transfer efficiency of particulate organic carbon from the surface to the deep ocean during the three hyperthermals, potentially reducing carbon sequestration in the deeper ocean interior. This implies that the warming-induced enhancement of microbial metabolism may have amplified global warming during ancient hyperthermals via a positive feedback mechanism.

Introduction
Rapid global warming events in the Phanerozoic, referred to as ancient hyperthermals, serve as analogues for predicting future climate changes1,2. These hyperthermal events3–6 are linked to remarkable perturbations in the global carbon cycle, as evidenced by large negative carbon isotope excursions (CIEs), rising atmospheric CO2 levels (pCO2)7–9, and biotic extinctions10. The rapid and substantial temperature rises are typically attributed to the release of massive quantities of CO2 from volcanic activities11,12, or biogenic methane emissions from the seafloor13,14. While these mechanisms emphasise carbon emissions as the primary driver of temperature increases, the role of the ocean’s biological carbon pump (BCP) in natural carbon sequestration and climate regulation is often overlooked.

The ocean’s BCP involves carbon fixation via photosynthesis at the surface, followed by the sinking of particulate organic carbon (POC) undergoing microbial respiration (remineralization) in the ocean interior, thereby cycling carbon and nutrients back to the surface15,16. This downward transport process stores carbon within the ocean interior where carbon can be sequestered on centennial and longer timescales17,18. Transfer efficiency, defined as the fraction of POC transported from the surface ocean (0−100m) to the ocean interior, serves as a widely used proxy for estimating the biological carbon pump sequestration19–21. Consequently, projecting the transfer efficiency of POC within the BCP in response to warming is essential to understand future climate change21,22.

Previous studies on changes in the BCP during ancient hyperthermals have mainly focused on changes in primary productivity or POC export production in surface water in response to rapid temperature rises23–27. These studies suggest either a weakened biological pump due to ocean stagnation caused by intensified stratification23,24 or a much stronger biological pump with high productivity driven by increased nutrients due to warming-induced intensive weathering25–27. However, little is known about the transfer efficiency of POC from the surface to the deep ocean during ancient hyperthermals, which hinders a clear understanding of the feedbacks between oceanic carbon sequestration and climate change across the Phanerozoic.

Temperature is widely recognized as a key factor that influences the rates of microbial photosynthesis and respiration28,29. Previous studies have demonstrated that respiration rates can double or triple with a ~10°C increase in ambient temperature28,29. Such temperature-driven changes in microbial remineralization strongly affect the degradation of sinking organic matter30,31. Together with ocean circulation and water mass ventilation, which regulate circulation timescales and pathways of sinking organic carbon, these processes collectively determine the transfer efficiency of POC21,32,33. Meanwhile, the vertical δ13C gradients of dissolved inorganic carbon (DIC), i.e., the decrease in DIC δ13C values (δ13CDIC) with increasing depth in the upper ocean (>1000m), are mainly controlled by the remineralization of particulate organic matter (POM) that releases isotopically lighter carbon at greater depths32. Therefore, vertical δ13C gradients are important geochemical indicators for tracing the transformation and remineralization of POC33–35. Here, we integrate proxy δ13C and simulated δ13C in an Earth system model to investigate changes in vertical δ13C gradients during three key Phanerozoic hyperthermals4,10: the Permian−Triassic Thermal Maximum (PTTM; ~251.9Ma), the early Toarcian Thermal Maximum in the Early Jurassic (ETTM; ~183.0Ma), and the Palaeocene−Eocene Thermal Maximum (PETM; ~56.0Ma). Utilising data-constrained cGENIE (carbon-centric Grid ENabled Integrated Earth system model), we evaluated the response of POC transfer efficiency to elevated temperatures during ancient hyperthermal events.

Results
The three ancient hyperthermals examined in this study feature prominent and rapid global temperatures rises and well-documented carbon isotopes from globally distributed sites, which allows for the reconstruction of δ13C depth profiles in response to warming. This is in contrast to other ancient hyperthermals with relatively smaller, slower temperature rate increases10 and/or poorly documented carbon isotope data (Fig. 1). The PTTM (also known as the Permian−Triassic mass extinction) is marked by ~10°C increase in sea surface temperatures (SSTs) estimated from low-latitude conodont oxygen isotopes6 and coeval negative CIEs (~4‰) recorded in global marine carbonates (δ13Ccarb)7. During the ETTM (also known as the early Toarcian oceanic anoxic event), oxygen isotopes and Mg/Ca ratios of the belemnite from the middle-latitude European Shelf show a ~7°C increase in SST, coinciding with a ~3‰ negative CIE in marine carbonates3,8. Biomarker proxy data based on TEX86 (TetraEther indeX of tetraethers consisting of 86 carbon atoms) from low-latitude regions indicate a ~4°C increase in SST during the PETM, which is further supported by Mg/Ca and oxygen isotope data5. Moreover, foraminifera shells exhibit a ~2.5‰ negative CIE during the PETM9. The pronounced negative CIEs recorded in marine carbonates and foraminifera shells allow for the reconstruction of δ13C depth profiles in the upper ocean.

Here, we compiled published δ13C data from 174 locations (Supplementary Fig. 1; Supplementary Data 1; Supplementary Table 1) and identified three depth-dependent depositional settings, including shallow shelf (0−80m, average ~40m), deep shelf (80−180m, average ~130m), and basin (180−280m, average ~230m) facies. The compiled δ13C data from the modern ocean demonstrate that δ13C values from continental shelf environments provide a robust proxy of the global ocean’s carbon isotopic composition (Supplementary Fig. 2). These depths correspond to the depths of the three uppermost ocean grid layers in the cGENIE Earth system model36 for direct model-data comparison (a sensitivity test for different scenarios of assumed depth is provided in Supplementary Table 2). Then, we reconstructed δ13C depth profiles during the pre-CIE (pre-warming) and the peak CIE (peak warming) phases but focused on global stacked δ13C depth profiles (Fig. 2 and Supplementary Figs. 3–6), because δ13C data from some oceanic regions are insufficient to reconstruct local gradients. One inherent challenge in reconstructing the vertical δ13C gradients is local carbonate diagenesis, including that induced by meteoric water and burial processes, which can produce anomalously large CIE37. We vetted our dataset by removing δ13C data altered by diagenesis and focused on relative changes to minimise the effects of diagenesis in our analysis (see “Methods”). Moreover, the average δ13C values in surface water (shallow shelf) were set as the reference standard and were adjusted to 0‰, while other proxy data were recalculated as offsets to eliminate isotopic differences between oceanic DIC and carbonate resulting from carbon isotope fractionation and carbonate diagenesis (see “Methods”). Surface water δ13C values were normalised to zero for both proxy records and model scenarios to facilitate consistent intercomparisons. Planktonic foraminiferal shells can record oceanic δ13C when using constrained shell size fractions to minimise physiological effects (see “Methods”), because smaller-size shells typically show lower δ13C and δ18O values due to vital effects34. Although planktonic foraminifera flourished after the Jurassic38 and, therefore, data are only available for the PETM, they provide additional evidence to support the observations obtained from δ13C data (see “Methods”).

Compiled global δ13C records exhibit a pronounced negative carbon isotope excursion at three water depth levels (Supplementary Fig. 7; and Supplementary Table 1), which is commonly attributed to the injection of large amounts of 13C-depleted carbon into the atmosphere-ocean system11,12. Although δ13C variations in the ocean-atmosphere carbon inventory would theoretically induce equivalent isotopic shifts in both surface and deep waters, our data indicate that the CIE magnitude in the deep shelf or basin is greater than that in the shallow shelf (e.g., PTTM: 4.5‰ vs. 3.5‰). This discrepancy mainly attributed to an enhanced vertical δ13C gradient39,40. Our compiled global δ13C records show a two- to three-fold increase in the vertical δ13C gradients between pre-CIE and peak CIE phases (Fig. 2; and Supplementary Table 1). Specifically, δ13C data indicate that the vertical δ13C gradient rose from 0.5±1.6‰/hm (pre-CIE) to 1.7±1.7‰/hm (peak CIE) across the PTTM, from 1.7±1.3‰/hm to 2.6±1.5‰/hm during the ETTM, and from 1.0±0.9‰/hm to 2.8±1.9‰/hm during the PETM. An increase in the vertical δ13C gradient was also recorded in the δ13C record of planktonic foraminifera during the PETM (Supplementary Figs. 8 and 9). Similar trends have also been observed in marine carbonates from depth-dependent depositional settings in South China during the PTTM40,41 as well as in planktonic foraminiferal shells from the warm early Eocene and mid-Miocene34,42, further supporting the global pattern of increasing vertical δ13C gradients during these ancient hyperthermals.

Considerable δ13C variability exists even within the same depositional settings (e.g., 3.8‰ in shallow shelf during the P–Tr interval), exceeding the 2–3‰ variability observed in modern δ13C (Supplementary Fig. 2). This greater variability may reflect mixing of DIC sources and fluctuating paleodepths when carbonate δ13C is used as a proxy. Both factors contribute to uncertainty in estimates of vertical δ13C gradients, potentially making changes in gradients statistically insignificant. Firstly, the δ13C values recorded in deep shelf and basin environments potentially represent a mixture of deep and shallow water DIC sources43. To evaluate potential mixing effects on gradient amplification, we implemented a simple isotopic mixing model (see “Methods”). Sensitivity tests reveal that increased gradient during the warming, as indicated by the δ13C data, show limited dependence on mixing source (Supplementary Fig. 10). The estimates (shown in Fig. 2) without regard to mixed sources provide conservative values, as the magnitude increase when raising the mixing proportion (Supplementary Fig. 10). Furthermore, in order to evaluate the potential uncertainty resulting from palaeodepth fluctuations (shallowing or deepening) at locations, we simulated 2000 different palaeodepth scenarios for each hyperthermal event to assess the vertical gradient responses to depth changes during the event (see “Methods”), as these changes could not be reliably quantified through existing proxy data. Results demonstrate that vertical gradients during the peak CIE stage remained significantly enhanced in most scenarios (Supplementary Figs. 11 and 12), though gradient magnitudes diminished, as more locations exhibited depth changes. Collectively, the amplified vertical δ13C gradients documented across widely distributed sites and multiple geological intervals strongly suggest that local effects, including palaeodepth variations, mixed DIC sources, or diagenetic alterations, do not constitute the dominant control.

To test the mechanistic drivers of the increased vertical δ13C gradients during the three investigated hyperthermals, we conducted a series of experiments using the cGENIE Earth system model (Fig. 2; Supplementary Table 4). We prescribed a range of atmospheric pCO2 values (0.5–9× pre-industrial level, 140–2520 ppmv) and obtained simulated SSTs to match the proxy SST data from the low- to mid-latitudes (Fig. 1a). Previous studies suggest that vertical δ13C gradients are mainly controlled by ocean nutrient conditions41, thermohaline circulation40, and temperature-induced changes in microbial metabolism39; however, a comprehensive assessment of their roles in driving δ13C gradients during hyperthermals through modelling is lacking. Accordingly, we constructed three scenarios for each hyperthermal event—circulation change, eutrophication, and a warming-driven increase in microbial metabolism (hereafter, microbial metabolism). The “circulation change” scenario was simulated solely by temperature-driven changes in physical ocean circulation within variant oceanic phosphate levels (1× the modern [PO4] level) between pre-CIE and peak CIE phases, whereas in the “eutrophication” scenario, increased [PO4] (2× the modern level following previous cGENIE experiments39) was set at the peak CIE. Both the “circulation change” and “eutrophication” scenarios were run using a static biological pump, in which global POM export production and remineralization are almost entirely controlled by nutrient availability (here, [PO4])26. In contrast, the “microbial metabolism” scenario used a constant 1× [PO4] and a dynamic biological pump, controlled by both nutrient availability and temperature-dependent metabolic rates. In particularly, the dynamic model incorporates a temperature-dependent scheme for microbial remineralization of organic matter to regulate POC flux in mid-waters44. In the “microbial metabolism” scenario, we reduced the POM sinking rate by ~22% in the PTTM experiments, following Hülse et al. based on marine animal body-size data39 (see “Methods”) to represent a shallower remineralization depth before the early Mesozoic. This is because modern marine phytoplankton (e.g., coccolithophorids, dinoflagellates, and diatoms) did not flourish until the Late Triassic38. The sinking rates for the ETTM and PETM were set at modern levels following previous modelling experiments22,34.

Our model shows that vertical δ13C gradients increase significantly during the peak CIE in the “microbial metabolism” scenario, matching proxy data better than the other two scenarios (Fig. 2). In this scenario, the differences between proxy and model data at peak CIE phase are −0.3±2.7‰/hm for the PTTM, 0.7‰/hm ±2.4‰/hm for the ETTM, and −0.04±2.9‰/hm for the PETM, while those for the other two scenarios are notably larger (1.0–2.3 ‰/hm) (Supplementary Table 4). This is further supported by root mean square error calculations that show the lowest deviation between proxy and model data under the “microbial metabolism” scenario (Supplementary Table 4). However, the modelled vertical δ13C gradient in the “microbial metabolism” scenario does not fully account for the observed data, possibly due to uncertainties in model boundary conditions (e.g., set temperature values) and/or proxy data (e.g., local diagenesis or mixing of deep and shallow water DIC37). In particular, there are large discrepancies between the proxy data (most from European shelf) and model results during the peak CIE phase of the ETTM (Fig. 2b, e). This mismatch likely arises from the relatively low spatial resolution of the cGENIE model, which does not adequately capture the restricted conditions of the European shelf and therefore fails to reproduce the 13C-depleted DIC at mid-water depths observed in the proxy data.

Our best-fit model scenario (i.e., “microbial metabolism”) supports the proposition that temperature-dependent microbial metabolism has a significant impact on vertical δ13C gradients while nutrients and ocean circulation are of secondary importance (Fig. 2; and Supplementary Fig. 6). Specifically, increased ocean temperatures during hyperthermals may have increased microbial respiration and resulted in a higher remineralization rate of sinking organic matter in mid-waters, generating more 13C-depleted DIC and larger vertical δ13C gradients (Supplementary Fig. 13). Our simulation highlights the impacts of temperature-dependent microbial metabolism on vertical δ13C gradients, which have been overlooked in previous studies (e.g., refs. 40,41), while not excluding the effects of nutrients and ocean circulation33,45. The “eutrophication” scenario demonstrates that doubled nutrient levels stimulate productivity (Supplementary Fig. 14), concurrently increasing the vertical δ13C gradient. However, this effect remains substantially weaker than the gradient increase induced by the “microbial metabolism” scenario (Fig. 2). Only extremely high nutrient levels (e.g., 10× [PO4]41) or regionally limited vertical mixing (e.g., in the modern Black Sea40) could generate large vertical δ13C gradients. Additionally, sensitivity experiments under the “microbial metabolism” scenario reveal that the reduced POC sinking rate in the cGENIE model exhibit limited influence on the gradients (Supplementary Fig. 15), reinforcing the primacy of temperature-dependent microbial metabolism.

Our simulations indicate an increase in POC export production in surface water during the three studied hyperthermals, particularly in the “microbial metabolism” and “eutrophication” scenarios compared to a minimal change in the “circulation change” scenario (Fig. 3a and Supplementary Fig. 14). In the “microbial metabolism” scenario, 122%, 83%, and 40% increases in POC export production are associated with the PTTM, ETTM, and PETM, respectively (Fig. 3a). A doubling of oceanic phosphate concentrations in the “eutrophication” scenario results in a ~100% increase in POC export production (Supplementary Fig. 14). Some model26 and proxy25,27,46 studies also support a rise in primary production or export production during these hyperthermals as a result of increased nutrient input (e.g., phosphate), driven by intensified weathering, similar to our “eutrophication” scenario. However, in the “microbial metabolism” scenario—as the best-fit model scenario supported by vertical δ13C gradients (Fig. 2)—the increased POC export production is primarily attributed to more recycled [PO4] returning to the surface water due to higher remineralization rates. This effect is further enhanced by warming-induced increases in photosynthetic rates, leading to higher productivity (Fig. 4). The “microbial metabolism” scenario incorporates the temperature dependency of the nutrient uptake rate in the model44. In contrast, the static model used in the “circulation change” and “eutrophication” scenarios does not account for temperature dependency (Supplementary Fig. 14). This implies that temperature-driven nutrient recycling plays a crucial role in stimulating POC export production.

Export production is an important metric for evaluating the sequestration capacity of the biological carbon pump15,47. However, quantifying export production alone is insufficient to assess the biological pump’s impact on atmospheric CO2, as its efficiency also depends on the residence time of exported carbon in the ocean17–19. Specifically, carbon in the upper ocean (<1000m) typically returns to the atmosphere on sub-centennial timescales, limiting its contribution to long-term CO2 storage, whereas carbon transported to the deep ocean (>1000m) or buried in sediments can remain sequestered for much longer periods20. Thus, the POC transfer efficiency at ~1000m is critical for predicting ocean carbon sequestration20,21. Accordingly, we focused on the ratio of the POC flux at ~1040m to that at the surface in the model as a proxy for the global mean transfer efficiency of POC21.

Despite the increased POC export production in the surface water, our results indicate that this compensatory flux is insufficient to fully offset the reduction in mid-water POC flux due to warming-induced increases in microbial respiration, leading to a significant decrease in POC transfer efficiency (Fig. 3b, c). Simulations showed a two- to six-fold decrease in POC transfer efficiency during warming, from 9.2% to 1.5% for the PTTM, 13.8% to 4.2% for the ETTM, and 1.9% to 0.9% for the PETM (Fig. 3c). Sensitivity analyses of other mid-water depths (~350m and ~640m) also revealed a significant decrease in transfer efficiency (Supplementary Fig. 16). Along with our sensitivity analyses of ocean temperature, the independence of transfer efficiency from [PO4] (Supplementary Fig. 17) and spatial pattern of higher transfer efficiency in the high latitudes in our simulations (Supplementary Fig. 18) further support the idea that transfer efficiency is inversely correlated with ocean temperature. Our spatial pattern is consistent with the results of data-assimilated biological pump models and modern field observations31, though quantitative discrepancies persist with direct sediment trap measurements20. This discrepancy potentially arises from the sole temperature-dependent remineralization process incorporated in our model, without considering zooplankton migration48 or size-dependent plankton communities49.

The primary triggers for rapid temperature rise during hyperthermals are typically massive volcanic CO2 and thermogenic methane release from large igneous provinces, with possible contributions from felsic volcanisms11,12,50, biogenic methane emissions from the seafloor13,14, and enhanced reverse weathering driven by the collapse of marine siliceous organisms51. Biological carbon pumps play a pivotal role in sequestering carbon emissions and reducing atmospheric pCO2. However, a reduction in POC transfer efficiency due to a warming-induced increase in microbial metabolism—as indicated by our model results (Figs. 3 and 4)—potentially lead to less CO2 sequestration in the ocean interior. This process may create potential positive feedback in amplifying the greenhouse effect during the three studied hyperthermal events, which persisted for hundreds of thousands to millions of years (Fig. 1).

We further perform “multi-factor” experiments varying [PO4] and temperature dependency simultaneously to investigate the interplay between nutrient-dependent productivity and temperature-dependent microbial metabolism on the BCP. Model results reveal that these competing drivers collectively reduce mid-water POC flux during the PTTM, while enhancing fluxes during the ETTM and PETM, despite maintaining consistently reduced POC transfer efficiency across all events (Supplementary Fig. 19). This opposite mid-water POC flux pattern may account for the contrasting OC burial responses observed among different hyperthermals52–55. Some hyperthermal events are associated with increased organic carbon (OC) burial, such as the widespread black shale deposition observed in the Toarcian OAEs52. However, there are also contrasting examples, such as the sudden decline in OC accumulation rates in South China during the PTTM53 and the suppression of global OC burial during mid-Miocene warming54. These contrasting observations of OC burial may imply that competitive factors control OC burial during hyperthermals56. Ocean anoxia, together with higher productivity triggered by increased weathering, has been widely used to explain elevated OC burial during OAEs52. In contrast, warming-induced increases in microbial metabolism led to reduced POC transfer efficiency and lower POC flux within the ocean interior, which likely have a counter effect on OC burial54. However, our simulated mid-water POC fluxes should not be conflated with OC burial, critical preservation processes of OC are not involved in our experiments57. A systematic integration of proxy data (e.g., TOC) with models explicitly incorporating OC burial processes will be essential to elucidate how these warming-induced competitive factors control OC burial across different hyperthermals.

Our findings highlight that temperature-dependent microbial metabolism may significantly influence BCP changes and thus climate regulation. Specifically, a temperature increase of 4–10°C during three ancient hyperthermals could reduce transfer efficiency of BCP by 2–6 times, thereby lowering the flux of POC into the ocean interior and potentially counteracting OC burial (Fig. 3). This feedback likely resulted in reduced CO2 sequestration in the ocean interior, amplifying ancient warming (Fig. 4). Conversely, a cooling-induced decline in microbial metabolism increased the efficiency of BCP over the past 15 Myr34. It should be noted, however, this critical feedback becomes more effective over longer timescales, as the effects of organic carbon burial on carbon cycling typically become apparent over million-year timescales57. Consequently, the amplification effect during the PETM may have been relatively limited due to its shorter duration than those of the PTTM and ETTM (Fig. 1), which restricted the establishment of sustained positive feedbacks. Furthermore, a key limitation of our simulations is the inability to quantify how the 2–6-fold reduction in transfer efficiency decreases the amount of CO2 sequestered in the ocean by the BCP over geological timescales, which should be addressed in future work.

Despite this limitation of our simulation, changes in the BCP during ancient hyperthermals in the Phanerozoic provide a strong analogue for the future BCP and anthropogenic climate change. Projections of POC transfer efficiency in response to 21st century climate change exhibit conflicting results across the Coupled Model Intercomparison Project Phase 6 (CMIP6) ensembles21. These discrepancies reflect diverse physical and biogeochemical processes represented by the CMIP6 model. The cGENIE model incorporates temperature-dependent microbial processes, particularly temperature-dependent remineralization (see “Methods”)21. Based on this, previous study using the cGENIE model project that globally averaged transfer efficiencies at a depth of 1040m will decline from ~11% to 4–9% under various carbon emission scenarios22. With a smaller fraction of organic matter reaching the ocean interior, the food supply and marine ecosystem in mid-waters are likely altered22. Moreover, a reduced amount of carbon transported to the seafloor and buried in sediments potentially creates a positive feedback loop, exacerbating global warming on long timescales, similar to the three ancient hyperthermals studied. However, this positive feedback on carbon sequestration depends on timescales. Model projections suggest a net increase in oceanic carbon storage by 2100 due to short timescales (<100 years), despite a decrease in transfer efficiency21. Furthermore, this study does not account for other factors, including trait-based functionality types of the plankton community49, oxygen-dependent processes58, and mineral ballasting effects59, all of which may influence transfer efficiency. Nevertheless, our integrated approach combining proxy data compilation with mechanistic modelling highlights that temperature-dependent microbial metabolism may play a significant role in modulating the BCP dynamics during both ancient hyperthermal events and future climate scenarios.

Methods
Carbonate carbon isotope compilation
We collected marine carbonate carbon isotope (δ13Ccarb) profiles from various water depth-dependent depositional settings on carbonate platforms or ramps to represent the δ13C in the upper ocean (Supplementary Data 1 and references in Supplementary Information). Although our compiled δ13C data are derived primarily from continental shelf settings and lack pelagic records, modern δ13C data demonstrate that δ13C values from continental shelf environments provide a robust representation of the global ocean (Supplementary Fig. 2). The selected δ13C profiles adhered to the following principles: (1) The δ13C data should predominantly be extracted from carbonate platforms or ramps to represent the upper ocean without drilling cores from the deep sea; (2) some of the extremely large CIE magnitudes, such as −26‰ caused by local authigenic carbonate at the Salisbury Embayment in the USA60, were excluded from our dataset as they are caused by regional factors; and (3) the high-resolution δ13C profiles preserve a clear expression of the negative CIE with steady pre-CIE, and prominent peak values were selected because this feature provides an unambiguous indicator and chronostratigraphic control of periods before and during hyperthermals. Biostratigraphy (conodonts for the PTTM, ammonites for the ETTM, and benthic foraminifer or nannofossils for the PETM) was also used to further ensure the age reliability of the carbon isotope profiles. Within the data compilation, we recognised three depositional settings based on lithology, special deposits (e.g., microbialites and oolites), depth-dependent fossils, and paleodepth calculations (Supplementary Data 1; Supplementary Information). These included shallow shelf (0–80m), deep shelf (80–200m), and basin facies (200–500m), which represent the three distinct water depth levels within the upper ocean. For example, for the Permian–Triassic, shallow shelf deposits are represented by microbialites or oolites, whereas a deep shelf is indicated by marlstone intercalated with mudstone, and a deep basin is indicated by siliceous limestones or rocks with abundant radiolarians (see Supplementary Information). Our dataset comprises a total of 85 δ13C profiles from the PTTM, 56 profiles from the ETTM, and 33 profiles from the PETM (Supplementary Data 1; Supplementary Table 1). The δ13C dataset for the Permian–Triassic interval is based on our previous research7 with additional data updates, incorporating a total of 60 shallow shelf profiles, 17 deep shelf profiles, and 8 basin profiles. The δ13C dataset for the Toarcian encompasses 11 shallow shelf profiles, 31 deep shelf profiles, and 14 basin profiles. Additionally, the δ13C dataset for the PETM includes 10 shallow shelf profiles, 17 deep shelf profiles, and 6 basin profiles.

Estimates of vertical δ13C gradients
To reconstruct the vertical δ13C gradients before and during hyperthermals, the average δ13C values (mean ± 1 s.d.) were calculated for each of the three palaeo-depths during the pre-CIE (pre-warming) and peak CIE (peak warming) phases. The Shapiro‒Wilk test (function “shapiro.test” in R) was performed to demonstrate that the δ13C values were statistically normally distributed. We then assumed the mean depths of the three depositional settings to be ~40m, ~130m, and ~230m, which correspond to the three uppermost ocean levels in the cGENIE Earth system model36. There are multiple reasons for simplifying the water depth of depositional settings, as follows: (1) Our aim was to comprehend the relative change in the vertical δ13C gradient between pre-CIE and peak CIE phases, rather than reconstructing continuous temporal patterns; (2) photosynthesis and organic matter remineralization occur in the surface ocean (<~100m) and mid-water (>100m), respectively36, which correspond to our assumed surface depth of 40m vs. mid-water depths of ~130 and ~230m in the model; and (3) both proxy and modelled δ13C using the same water depth facilitate their direct comparisons. The gradients calculated by the δ13C differences between the shallow shelf (40m) and the deep shelf or basin (130m or 230m) were then calculated to determine the gradients. These gradients were expressed as δ13C change per hundred metres in water depth (mean ± 1 s.d. ‰/hm). We concentrated specifically on the gradients calculated by the deep shelf because the data from the basin were insufficient for reliable gradient estimation (Supplementary Fig. 6). A two-sample t-test (function “t.test” in R) was performed to demonstrate the statistical significance of the difference in δ13C between surface and mid-water (Supplementary Table 1). Additionally, we calculated the vertical δ13C gradients under other two assumed water depth scenarios in depositional settings for comparison purposes (Supplementary Table 2).

The observed δ13C values at the deep shelf or basin may reflect a mixing signal of both deep and shallow water DIC. To assess the uncertainty in the mixing proportion, we randomly sampled 1000 different mixing proportions of shallow water DIC (0~90%). Using a simple isotopic mixing model of two end-members (Supplementary Information), then calculated the actual δ13C values of deep water (δ13C_DIC-deep) and the vertical gradient (Supplementary Fig. 10). Furthermore, in order to assess the impact of depth uncertainty on the estimated vertical gradient, we consider two extreme scenarios: (1) the shallow shelf deepening to a deep shelf (“depth increase” scenario) and (2) the deep shelf becoming shallower to a shallow shelf (“depth decrease” scenario). In the “depth increase” scenario, we randomly sample several shallow-shelf sections and reclassify them as deep-shelf locations. In the “depth decrease” scenarios, we randomly sample several deep-shelf locations and reclassify them as shallow-shelf sections. For each scenario, 1000 random sampling iterations were performed (Supplementary Figs. 11 and 12).

Carbon isotopes of planktonic foraminifers
In addition to carbonate δ13C, the δ13C of planktonic foraminiferal shells can indicate oceanic δ13C when considering shell size fractions and physiological disequilibrium effects61. Meanwhile, planktonic foraminiferal δ18O is a powerful proxy for sea temperature that decreases with depth, indicating paleodepth34,62. Therefore, the paired planktonic foraminiferal δ13C and δ18O profiles can be used to reconstruct the δ13C depth profiles. In our dataset, the shells of 24 profiles used a constrained size fraction (212–355μm) following previous work34, while the remaining 8 samples either lack size information or used a smaller size fraction (180–212μm). Smaller shells typically show lower δ13C and δ18O values due to vital effects34, which could lead to lower estimated δ13C and shallower inferred water depths. Consequently, the use of smaller shell sizes may slightly exaggerate vertical δ13C gradients, although this bias is likely limited given the constrained size fractions applied to most samples. Deep-ocean benthic foraminifera are excluded from this approach, as previous studies have noted that the vertical δ13C gradient between planktonic and deep benthic foraminifera decreased during the PETM63,64. One possible explanation is that the smaller CIE magnitude observed at pelagic sites may result from incomplete preservation due to carbonate dissolution64, rendering them unsuitable for reconstructing the vertical gradient. However, this method is only possible for the PETM and, therefore, was considered a supplementary method in this study. A total of 32 paired planktonic foraminiferal δ13C and δ18O profiles from nine locations in the PETM were collected (Supplementary Data 1 and references in the Supplementary Information). The estimation of paleodepth was based on planktonic foraminiferal δ18O and modelled salinity and temperature, following a previous study34 (Supplementary Information).

Description of the Earth system model
We used ‘cGENIE’ in this study, an Earth system model of intermediate complexity composed of a 3-D dynamic ocean circulation model with a 2-D energy and moisture balance model of the atmosphere, a dynamic ocean-climate model65, and a representation of the biogeochemical cycling of elements and isotopes in the ocean36. Our simulations integrated two types of biological pump schemes. First, a static model with a fixed double-exponential profile of POM flux has been applied in previous simulation works26. Second, an alternative approach is a dynamic model that incorporates temperature-dependent ocean microbial processes44. Specifically, the dynamic model incorporates Arrhenius-type temperature-dependent scheme for the microbial remineralization of organic matter to regulate POC flux in mid-waters44. Nutrient inventories, which are crucial to controlling the biological pump, have been simplified in our simulations by including only phosphate66. The biological pump parameters used in our simulations were based on modelling of middle Miocene cooling34. All the simulations were configured on a 36×36 grid with 16 vertical levels in the ocean. The palaeobathymetry, continental configuration, and wind stress of the three hyperthermal events were primarily derived from previous cGENIE simulations (PTTM67, ETTM68, and PETM12).

Description of the numerical experiments
The δ13C values of CO2 forcings were unchanged (−6.5‰) in the model during pre-CIE and peak CIE, as our simulations did not aim to investigate the impact of the massive light carbon release on the ocean δ13C. The changes in ocean δ13C are exclusively caused by various temperature-induced changes to the biological pump, including “circulation change”, “eutrophication”, “microbial metabolism” and “multi-factor” scenarios. A comparison of the static and dynamic biological pump models and the change in ocean phosphate concentration ([PO4]) between the pre-CIE (pre-warming) and peak CIE (peak warming) phases was designed to distinguish the proposed scenarios (Supplementary Table 4). Four scenarios are presented, as follows: (1) “circulation change” (temperature-driven physical ocean circulation change): 1.0× [PO4] during both pre-CIE and peak CIE with the static model; (2) “eutrophication” (nutrition-driven): 1.0× [PO4] at pre-CIE and 2.0× [PO4] at peak CIE with the static model; (3) “microbial metabolism” (warming-driven increase in microbial metabolism): 1.0× [PO4] at pre-CIE and peak CIE with the dynamic model; (4) “multi-factor”: 1.0× [PO4] at pre-CIE and 2.0× [PO4] at peak CIE with the dynamic model. A total of six experiments were designed for each hyperthermal, which included two pre-CIE experiments and four peak CIE experiments.

We prescribed a range of atmospheric pCO2 values (0.5–9× pre-industrial level, that is, 140–2520 ppmv) chosen to match the increase in SST from low or middle latitudes reconstructed from proxy records (Fig. 1). The modelled SST in the tropical Tethys Ocean for the PTTM rises from 27.3±1.2 °C to 37.0±1.1 °C; the modelled SST for the middle-latitude European Shelf for the ETTM increases from 20.8±0.6 °C to 27.8±0.6 °C; the simulated SST in the tropical Atlantic for the PETM increases from 33.1±0.4 °C to 37.3±0.4 °C. Furthermore, it has been observed that modern marine phytoplankton (e.g., coccolithophorids, dinoflagellates, and diatoms) occurred and flourished after the Mesozoic Marine Revolution38,69–71. Therefore, as for “microbial metabolism” scenario, we set the default sinking rate of the POM during the ETTM and the PETM after the Mesozoic Marine Revolution to be the same as the modern-day value (125 m yr-1). For the PTTM, a relatively slower sinking rate (98 m y-1) was determined based on changes in marine animal body sizes72, as suggested in a previous cGENIE simulation39.

Collectively, we performed three ensembles of 18 main experiments based on the Permian-Triassic, Toarcian, and late Palaeocene continental configurations. To further explore the contributions of environmental (temperature, phosphate, and continental configuration) and biotic (sinking rate) factors to the ocean carbon cycle, a series of sensitivity tests were conducted. We ran 50 model experiments using an ensemble of ocean temperatures (1, 2, 4, 8, and 16× the pre-industrial level of 280 ppmv), phosphate (1.0, 1.5, 2.0, 2.5, and 3.0× the modern [PO4] level of 2.159 μmol P kg−1) and a dynamic or static biological pump model with the Permian–Triassic configuration (Supplementary Fig. 13). An ensemble comprising five experiments was conducted with varying sinking rates (90, 100, 110, 120, and 130 m yr−1) based on the Permian–Triassic continental configuration using the dynamic model (Supplementary Fig. 15). All the simulations were run for 10 kyr to reach equilibrium.

Model and proxy data comparison
As vertical δ13C gradients vary across different latitudes according to observations from the modern ocean (Supplementary Fig. 2d-f). All averages of modelled δ13C values were calculated based on modelled results gathered from low to middle latitudes (48.6° S to 48.6° N), where most of the sections are located (Supplementary Fig. 1). As sampling model results based on latitudinal regions helps avoid issues arising from spatial unevenness. Before any comparison between the proxy and model, δ13C values need to be adjusted. The average proxy and modelled δ13C values at the surface water (40m) are designated as the reference standard and adjusted to 0‰, and then other proxy and model data were shifted accordingly by their corresponding reference standard. The proposed centring method eliminated the negative CIE in proxy data during peak CIE phases resulting from carbon emissions not included in the simulations. This also removed any discrepancy between the modelled δ13C values and proxy δ13C values caused by fractionation between oceanic DIC and carbonate (invariant and independent of temperature) and diagenesis. Root mean square error calculations (function “rmse” within the package “Metrics” in R) were performed to compare the vertical δ13C gradient among the three model scenarios and proxy data and determine the scenario that best fits the proxy data (Supplementary Table 2).

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "Enhanced microbial metabolism during hyperthermal events would have led to diminished Corg sequestration in the ocean interior amplifying global warming, but meanwhile during warm periods, enhanced stratification caused by slowered ocean circulation and promoted marine primary productivity would have collectively led to expansion of OMZ and favored Corg burial. Thus, the carbon cycling and related feedbacks during hyperthermal events remain unclear. The paper attempts to better understand the processes and feedback through data compilation and Earth system modeling for the three representative hyperthermal events during the Phanerozoic. The paper focuses on the role of biological C pump to the global warming, which present two-fold main processes through Earth system modeling, involving elevated rates of Corg remineralization and reduced transfer efficiency of POC from the surface ocean to the deep ocean. The paper also explores relationship between evolution of phytoplankton and BCP and further projects into the current and future global warming. This is a well-written, easy-reading manuscript but bears significant implications for the role of BCP under the global warming conditions, and for readers of broad interests. I therefore recommend publication after minor revisions.\nLine 19: What Earth system modeling used in the present study. Would be good to explain it for more information.\nLine 52: delete relatively; little is little.\nIncreased POC transfer efficiency after the rise of modern phytoplankton, but with evolution of various phytoplankton and various particle size etc would also impact the POC transfer efficiency. Thus, this process may be illustrated for an individual event, i.e., before and after, regarding the microbial evolution, but when compared among events, it might be more complicated.\nFigure1. Explain the regression curves and shades in b\nFigure 2. Authors may squeeze a-f panels vertically to save some space given the large white space in these panels.\nSupplementary Fig. 2. What are the different symbols for? The eight regional areas? If so, please indicate it in the Figure caption. Also, how to recognize proxy and modeled data? Need to indicate it too. I guess the questions work for Suppl. Fig. 3 and Fig. 4, too.\nSupplementary Fig. 5. h-j, should be g-i?—the caption is also a bit confusing what does peak CIE vs pre-CIE mean? Does that mean peak value minus pre-CIE value? Need more explanation.", "reviewer_2": "Yuyang Wu et al. with manuscript entitled “Microbial metabolism amplified warming in three Phanerozoic hyperthermal events”, address a very important topic of highly interdisciplinary interest for both understanding deep time and future climate change. In particular, the authors address if there were potential change in the mean ocean biological carbon pump (BCP) during three major global warming events in Earth history (Permian-Triassic, Toarcian and Palaeocene-Eocene events).\nTo test for a change in BCP) the authors focus on carbon isotopes in carbonates from different interpreted ancient water depths. Importantly, the results are then compared to different scenario testing with a global intermediate complexity earth system mode. Based on the large data compilation and model experiments the authors conclude that there likely was a significant reduction of BCP during the events.\nWhile I am vary sympathetic to the authors approach, method and results, there are three major problems with the implicit assumptions in the interpretation of the observations:\n1) The authors assume they can infer paleowater depth of deposition from sedimentological observations or petrology, while in fact their proxies used are only relative indicators in part related to climate oceanographic changes such as storm wave base, and thermal stratification. To put it another way, a marl deposited at one water depth before an event could very well be deposited at a shallower or deeper depth during the event. This would compromise the approach.\n2) The authors assume the carbonate sampled at a given interpreted depth where also precipitated at that depth. Carbonate mud or particles may be dispersed from shallow to deep depending on climate. As far as I can tell the authors do not present or discuss indicators that can confirm that the d13C signal at depth is not derived from shallow water and if this redistribution process could be related to the climate change during the events. Presumably the carbonate in marls interpreted as their deep shelf is coming from micritic carbonate, was this micrite originally produced in shallow water or locally? It is well known that the carbonate sources in platform carbonate factory have quite variable d13C_carb composition. Also, could it be that platform carbonate mud from the platform was not as efficiently dispersed to greater depth during sea level rise, and therefor indirectly causing an increased gradient in d13C.\n3) The authors assume that any change in the vertical d13C gradient is related mainly to a d13C change at depth. This assumption builds on the zeroing of near surface d13C_carb. I see no arguments that rule out that the change is actually a change in surface d13C_carb while the deep d13C_carb is semi fixed? This would implicate that it it could be a d13C mixed layer&atmosphere change and not a result of change in BCP?\nOn a side note the title is misleading as well as many portions of the text. The authors talk about microbial metabolism or microbial respiration when in fact they can only talk about remineralization rate of sinking organic matter. We can not know if the reminiralizatrion is a result mainly of single cell microbes or multicellular eukaryotes that performed the respiration.\nIn conclusion, considering point 1-3, the manuscrip is too immature for publication at this stage. However I do look forward to see further enhancements to the analysis in the future.\nChristian J. Bjerrum", "reviewer_3": "Point-by-point response to the reviewers’ comments.\nThose changes appear as highlights in red colour in the manuscript. Below are our responses (blue\ncolour) to reviewers’ comments. Line numbers in our response refer to the revised manuscript.\nREVIEWER COMMENTS\nReviewer #1 (Remarks to the Author):\nEnhanced microbial metabolism during hyperthermal events would have led to diminished Corg\nsequestration in the ocean interior amplifying global warming, but meanwhile during warm periods,\nenhanced stratification caused by slowered ocean circulation and promoted marine primary\nproductivity would have collectively led to expansion of OMZ and favored Corg burial. Thus, the\ncarbon cycling and related feedbacks during hyperthermal events remain unclear. The paper\nattempts to better understand the processes and feedback through data compilation and Earth system\nmodeling for the three representative hyperthermal events during the Phanerozoic. The paper\nfocuses on the role of biological C pump to the global warming, which present two-fold main\nprocesses through Earth system modeling, involving elevated rates of Corg remineralization and\nreduced transfer efficiency of POC from the surface ocean to the deep ocean. The paper also\nexplores relationship between evolution of phytoplankton and BCP and further projects into the\ncurrent and future global warming. This is a well-written, easy-reading manuscript but bears\nsignificant implications for the role of BCP under the global warming conditions, and for readers of\nbroad interests. I therefore recommend publication after minor revisions.\nWe appreciate your positive comments for this manuscript.\nLine 19: What Earth system modeling used in the present study. Would be good to explain it for\nmore information.\nThanks for your suggestion. We have pointed out that the Earth system model cGENIE is used in\nthis study: “Simulations with the Earth system model cGENIE…”.\nLine 52: delete relatively; little is little.\nCorrection made\nIncreased POC transfer efficiency after the rise of modern phytoplankton, but with evolution of\nvarious phytoplankton and various particle size etc would also impact the POC transfer efficiency.\nThus, this process may be illustrated for an individual event, i.e., before and after, regarding the\nmicrobial evolution, but when compared among events, it might be more complicated.\nWe appreciate the reviewer’s insightful comment that different phytoplankton types and particle\nsizes indeed influence POC transfer efficiency, similar to the concerns raised by Reviewer #3.\nFollowing Reviewer #3’s suggestions, we have removed related inferences, as our simulation does\nnot provide robust results to assess the impact of different phytoplankton types on transfer efficiency.\nFigure1. Explain the regression curves and shades in b\nWe have explained the regression curves and shades in the caption of figure. “The blue and red lines\nrepresent the LOESS fit curve for δ13C and SST, respectively, while light blue and red shaded area\nrepresent 68% confidence intervals”.\nFigure 2. Authors may squeeze a-f panels vertically to save some space given the large white space\nin these panels.\nThanks. We have reduced the space between panels in the Figure 2.\nSupplementary Fig. 2. What are the different symbols for? The eight regional areas? If so, please\nindicate it in the Figure caption. Also, how to recognize proxy and modeled data? Need to indicate\nit too. I guess the questions work for Suppl. Fig. 3 and Fig. 4, too.\nThe different symbols represent proxy data and three model scenarios. The legends of symbols have\nbeen added in the Supplementary Figures 2 to 4, which are same as those in the Figure 2.\nSupplementary Fig. 5. h-j, should be g-i?—the caption is also a bit confusing what does peak CIE\nvs pre-CIE mean? Does that mean peak value minus pre-CIE value? Need more explanation.\nCorrection made. “peak CIE vs pre-CIE” has been changed to “peak CIE minus pre-CIE”, and we\nadd more explanation.\nReviewer #2 (Remarks to the Author):\nYuyang Wu et al. with manuscript entitled “Microbial metabolism amplified warming in three\nPhanerozoic hyperthermal events”, address a very important topic of highly interdisciplinary\ninterest for both understanding deep time and future climate change. In particular, the authors\naddress if there were potential change in the mean ocean biological carbon pump (BCP) during three\nmajor global warming events in Earth history (Permian-Triassic, Toarcian and Palaeocene-Eocene\nevents).\nTo test for a change in BCP) the authors focus on carbon isotopes in carbonates from different\ninterpreted ancient water depths. Importantly, the results are then compared to different scenario\ntesting with a global intermediate complexity earth system mode. Based on the large data\ncompilation and model experiments the authors conclude that there likely was a significant\nreduction of BCP during the events.\nMany thanks. We feel that is a fair summary in-line with our intention.\nWhile I am vary sympathetic to the authors approach, method and results, there are three major\nproblems with the implicit assumptions in the interpretation of the observations:\n1) The authors assume they can infer paleowater depth of deposition from sedimentological\nobservations or petrology, while in fact their proxies used are only relative indicators in part related\nto climate oceanographic changes such as storm wave base, and thermal stratification. To put it\nanother way, a marl deposited at one water depth before an event could very well be deposited at a\nshallower or deeper depth during the event. This would compromise the approach.\nWe agree that uncertainty of water depth change should be evaluated, especially the shallower or\ndeeper depth during the event due to the limitations of our proxies. The interval we focused from\npre-CIE to peak CIE is usually short, lasting 0.05~0.2 Myr during three hyperthermal events\n(Burgess et al., 2014; Nordt et al., 2022; Cramwinckel et al., 2018). Thus, the depositional setting\nusually keep constant, but the water depth potentially fluctuates and become shallower or deeper.\nIt is difficult to quantify the change of water depth by proxies, but we could assume multiple possible\nscenarios and evaluate the impact of depth change on the estimated vertical gradient. We consider\ntwo extreme scenarios: (1) the shallow shelf deepening to a deep shelf (“depth increase” scenario)\nand (2) the deep shelf becoming shallower to a shallow shelf (“depth decrease” scenario). For\nexample, in the “depth increase” scenario, we randomly sample several shallow-shelf locations and\nreclassify them as deep-shelf locations at the peak CIE stage. The results of 2000 random sampling\niterations for each hyperthermal show that the vertical gradients still larger obviously at peak CIE\nin most scenarios, except for some shallower scenarios from early Toarcian (Supplementary Figs.\n10 and 11). These analyses suggest that even allowing for extreme scenarios of depositional setting\nchange, increased trend of vertical gradient during the warming, as indicated by the δ13C data, is\ncarb\nstill relatively robust.\nWe have added relevant sensitivity tests (Supplementary Fig. 9) and a description of the methods\n(lines 353–365) to demonstrate that our key conclusion remains robust despite variations in water\ndepth. Some sentences are added on lines 140-149: “Furthermore, palaeodepths of locations may\nhave fluctuated (shallowing or deepening) during the event. Despite these changes could not be\nreliably quantified through existing proxy data, we assume thousand deepening or shallowing\nscenarios during the event to assess vertical gradient responses to depth changes during the event\n(see Methods). Results demonstrate that vertical gradients during the peak CIE stage remained\nsignificantly enhanced in most scenarios (Supplementary Figs. 10 and 11), though gradient\nmagnitudes diminished, as more locations exhibited depth changes. Collectively, the amplified\nvertical δ13C gradients documented across widely distributed sites and multiple geological\nDIC\nintervals strongly suggest that local effects, including palaeodepth variations, mixed DIC sources,\nor diagenetic alterations, do not constitute the dominant control”.\nSupplementary Fig. 10. Effects of “depth increase” scenario on vertical carbon isotope gradients. We\nrandomly select some locations of shallow shelf deepening to a deep shelf at peak CIE. The results of\n1000 random sampling iterations show that the vertical gradients still larger obviously at peak CIE\nin most scenarios. Red lines represent LOESS fit through the data. (a) Permian-Triassic; (b) Early\nToarcian; (c) Paleocene-Eocene.\nSupplementary Fig. 11. Effects of “depth decrease” scenario on vertical carbon isotope gradients. We\nrandomly select some locations of deep shelf becoming shallower to a shallow shelf at peak CIE. The\nresults of 1000 random sampling iterations show that the vertical gradients still larger obviously at\npeak CIE in most scenarios. Red lines represent LOESS fit through the data. (a) Permian-Triassic;\n(b) Early Toarcian; (c) Palaeocene-Eocene.\n2) The authors assume the carbonate sampled at a given interpreted depth where also precipitated at\nthat depth. Carbonate mud or particles may be dispersed from shallow to deep depending on climate.\nAs far as I can tell the authors do not present or discuss indicators that can confirm that the d13C\nsignal at depth is not derived from shallow water and if this redistribution process could be related\nto the climate change during the events. Presumably the carbonate in marls interpreted as their deep\nshelf is coming from micritic carbonate, was this micrite originally produced in shallow water or\nlocally? It is well known that the carbonate sources in platform carbonate factory have quite variable\nd13C_carb composition. Also, could it be that platform carbonate mud from the platform was not\nas efficiently dispersed to greater depth during sea level rise, and therefor indirectly causing an\nincreased gradient in d13C.\nWe acknowledge that the potential mixing of deep and shallow water DIC may influence the\nestimates. Specifically, the observed carbonate δ13C at the deep shelf may reflect a mixing signal of\nboth deep and shallow water DIC (Shen et al., 2011). If this is the case, the actual DIC δ13C of deep\nwater is likely lower than the carbonate δ13C observed at the deep shelf, meaning that our previous\nresults may have underestimated the vertical gradient.\nTo assess the uncertainty in the mixing proportion, we randomly sampled 1,000 different\nmixing proportions of shallow water DIC (pr = 0 to 90%). Using a simple isotopic mixing model of\ntwo end-members, we then calculated the actual δ13C of deep water (δ13C ) and the vertical\nDIC-deep\ngradient.\nδ13C = pr*δ13C + (1−pr)* δ13C\ncarb-deep DIC-shallow DIC-deep\nδ13C = (δ13C −pr*δ13C )/(1−pr)\nDIC-deep carb-deep DIC-shallow\nThe δ13C data from both the shallow and deep shelves are available (e.g. δ13C = 1‰;\ncarb carb-shallow\nδ13C = −1‰), and assuming δ13C could represent the δ13C . According to the\ncarb-deep carb-shallow DIC-shallow\nequation, when pr =0.5, δ13C is −3‰ and DIC gradient is −4‰. First, regardless of the mixing\nDIC-deep\nproportion, the gradient during the peak stage is generally larger than that during the pre-CIE stage,\nand increases as the proportion increases (Supplementary Fig. 9). Even if the proportion differs\nbetween the two stages due to climate change or sea level rise, as suggested by the reviewers, our\nanalysis, which includes 1,000 × 1,000 scenarios (Supplementary Table 3), shows that in most cases,\nthe vertical gradient still increases. For example, 91% of the scenarios at the P-Tr boundary indicate\nan increased gradient (Supplementary Table 3). These analyses suggest that increased gradient\nduring the warming, as indicated by the δ13C data, is largely independent of the mixing proportion\ncarb\nbetween shallow and deep water DIC. Moreover, previous estimates assuming no mixing (pr = 0)\nrepresent a conservative estimate of the gradient increase during warming.\nWe have added some brief description of sensitivity analyses in the Methods (lines 353-357) and\ndetails in supplementary materials. Some sentences discussing the mixing proportion are added in\nmanuscript. Lines 134-140: “The δ13C values recorded in deep shelf and basin environments\ncarb\npotentially represent a mixture of deep and shallow water DIC sources44. To evaluate potential\nmixing effects on gradient amplification, we implemented a simple isotopic mixing model (see\nMethods). Sensitivity tests reveal that increased gradient during the warming, as indicated by the\nδ13C data, show limited dependence on mixing source (Supplementary Fig. 9). The estimates\ncarb\n(shown in Fig. 2) without regard to mixed sources provide conservative values, as the magnitude\nincrease when raising the mixing proportion (Supplementary Fig. 9)”.\nSupplementary Fig. 9. Effects of mixing source of shallow-water DIC on vertical carbon isotope\ngradients. Mixing proportions (0~90%;1000 random samples) of shallow-water DIC to carbonate\nsediments on the deep shelf. The vertical δ13C gradients (mean + s.e.m.) are calculated from the δ13C\nDIC\ndifference between surface water (~40 m) and mid-water (~130 m), expressed as the δ13C change per\nDIC\n100 m of water depth (‰/hm). Blue and red points represent estimates shown in Fig. 2 without\nconsidering a mixing source. Sensitivity tests reveal that the trend of increased gradients at the peak CIE\nstage, as indicated by the δ13C data, is largely independent of the mixing proportion during the\ncarb\nPermian-Triassic (a), Early Toarcian (b) and Palaeocene-Eocene Thermal Maximum (c).\nSupplementary Table 3. Evaluation of potential vertical gradient changes during peak CIE under varying\nmixing source conditions. This comprehensive analysis examines 1,000 × 1,000 possible scenarios\n(Supplementary Fig. 9) comparing vertical gradient increases versus decreases at the peak CIE stage,\nunder varying mixing source conditions.\nAge Number of increased Number of decreased Ratio of increased Ratio of decreased\ngradient scenarios gradient scenarios gradient scenarios gradient scenarios\nPTTM 911,307 88,693 91% 9%\nETTM 722,435 277,565 72% 28%\nPETM 888,295 111,705 89% 11%\n3) The authors assume that any change in the vertical d13C gradient is related mainly to a d13C\nchange at depth. This assumption builds on the zeroing of near surface d13C_carb. I see no\narguments that rule out that the change is actually a change in surface d13C_carb while the deep\nd13C_carb is semi fixed? This would implicate that it it could be a d13C mixed layer&atmosphere\nchange and not a result of change in BCP?\nWe appreciate the reviewer's insightful comment regarding the interpretation of vertical δ13C\ngradient changes. It is necessary to analyze the original carbonate δ13C data (Supplementary Fig. 5\nand Supplementary Table 1). The δ13C at different depth all become negative during the event (i.e.\ncarb\nnegative carbon isotope excursion), because of the massive 13C-depleted carbon injected into\natmosphere from large igneous province (Gutjahr et al., 2017; Wu et al., 2023). While atmospheric\nδ13C changes would indeed produce equivalent isotopic shifts in surface and deep waters, our data\nshow that the CIE magnitude at deep shelf (4.5‰) is larger than those from shallow shelf (3.5‰).\nThis depth-dependent signal cannot be explained solely by atmospheric changes, but need\nadditionally mechanisms involved increased vertical δ13C gradient related to biological pump\nchange (Hülse et al., 2021).\nThe method zeroing surface δ13C is to eliminate the negative CIE in proxy data during peak CIE\ncarb\nphases resulting from carbon emissions not included in the simulations, and also removed any\ndiscrepancy between the modelled δ13C values and proxy δ13C values caused by fractionation\nDIC carb\nbetween oceanic DIC and carbonate (invariant and independent of temperature) and diagenesis.\nTo strengthen the argumentation, we have made the following substantive revisions in the lines 116-\n123:“Compiled global δ¹³Ccarb records exhibit a pronounced negative carbon isotope excursion\n(CIE) at three water depth levels (Supplementary Fig. 6; Supplementary Table 1), which is\ncommonly attributed to the injection of large amounts of ¹³C-depleted carbon into the atmosphere-\nocean system11, 12. Although δ¹³C variations in the ocean-atmosphere carbon inventory would\ntheoretically induce equivalent isotopic shifts in both surface and deep waters, our data indicate\nthat the CIE magnitude in the deep shelf or basin is greater than that in the shallow shelf (e.g.,\nPTTM: 4.5‰ vs. 3.5‰). This discrepancy mainly attributed to an enhanced vertical δ¹³C gradient39,\n40”.\nOn a side note the title is misleading as well as many portions of the text. The authors talk about\nmicrobial metabolism or microbial respiration when in fact they can only talk about remineralization\nrate of sinking organic matter. We can not know if the reminiralizatrion is a result mainly of single\ncell microbes or multicellular eukaryotes that performed the respiration.\nThe remineralization scheme in the cGENIE model follows an Arrhenius-type temperature\ndependency (Crichton et al., 2021). This simplified scheme applies to microbes (e.g., bacteria,\narchaea, and single-celled eukaryotes) rather than complex multicellular eukaryotes, which has been\nused to investigate ocean microbial remineralization in deep time, such as the P-Tr mass extinction\nand Cenozoic (Hülse et al., 2021; Boscolo-Galazzo et al., 2021). Therefore, it is clear that our focus\nis on microbial metabolism rather than multicellular eukaryotes.\nWe have added details of this temperature-dependent scheme in the Methods (lines 386-388):\n“Specifically, the dynamic model incorporates Arrhenius-type temperature-dependent scheme for\nthe microbial remineralization of organic matter to regulates POC flux in mid-waters”.\nIn conclusion, considering point 1-3, the manuscrip is too immature for publication at this stage.\nHowever I do look forward to see further enhancements to the analysis in the future.\nChristian J. Bjerrum\nThank you for your constructive comments. These suggestions really help improve our manuscript.\nReferences\nBoscolo-Galazzo, F., Crichton, K.A., Ridgwell, A., Mawbey, E.M., Wade, B.S. and Pearson, P.N.,\n2021. Temperature controls carbon cycling and biological evolution in the ocean twilight zone.\nScience 371, 1148-1152.\nBurgess, S.D., Bowring, S. and Shen, S.. High-precision timeline for Earth's most severe extinction.\nProceedings of the National Academy of Sciences 111, 3316–3321 (2014).\nCramer, B.S., Toggweiler, J.R., Wright, J.D., Katz, M.E. & Miller, K.G. Ocean overturning since\nthe Late Cretaceous: Inferences from a new benthic foraminiferal isotope compilation.\nPaleoceanography 24, PA4216 (2009).\nGutjahr, M. et al. Very large release of mostly volcanic carbon during the Palaeocene–Eocene\nThermal Maximum. Nature 548, 573-577 (2017).\nHülse, D. et al. End-Permian marine extinction due to temperature-driven nutrient recycling and\neuxinia. Nat. Geosci. 14, 862-867 (2021).\nNordt, L., Breecker, D. & White, J. Jurassic greenhouse ice-sheet fluctuations sensitive to\natmospheric CO dynamics. Nat. Geosci. 15, 54-59 (2022).\n2\nShen, B. et al. Carbon, sulfur, and oxygen isotope evidence for a strong depth gradient and oceanic\noxidation after the Ediacaran Hankalchough glaciation. Geochim. Cosmochim. Acta 75, 1357-\n1373 (2011).\nWu, Y.Y. et al. Volcanic CO degassing postdates thermogenic carbon emission during the end-\n2\nPermian mass extinction. Sci. Adv. 9, eabq4082 (2023).\nReviewer #3 (Remarks to the Author):\nWu et al., calculate the vertical gradient of ∂13C observations across prominent hyperthermal events\nto show significant increases that may reflect a change in watercolumn cycling of organic carbon\n(the Biological Pump). They further use the cGENIE Earth system model to help interpret these\nsignals in the context of changes in ocean circulation, nutrient inventories, and changes in particulate\norganic carbon flux attenuation (transfer efficiency) with depth due to temperature dependent\nremineralisation. The model results demonstrate that changing transfer efficiency with warming can\nsufficiently explain the ∂13C gradients. The authors conclude that the biological pump may\ntherefore have had more of a role in regulating the carbon cycle during these events than previously\nthought.\nI think the integration of observations and modelling here is robust and informative. In particular,\nthe model is used to help robustly identify TE from a range of potential drivers as an explanation of\nthe ∂13C gradients. I have a few comments on clarifying this approach below but I think this is as\nrobust as a paleo data-model integration can be (I can’t comment on the collation and treatment of\nobservations directly). This outcome is very interesting, potentially significant, and to my\nknowledge has not been shown before.\nHowever, I think the manuscript goes too far in terms of interpretation and implications beyond the\ncore hyperthermal TE and ∂13C result. The work on sinking rates as a function of community\nstructure is too simple and not robust to support the level of interpretation and discussion in the\nmanuscript. Likewise, the discussion of implications for future changes in the biological pump and\ncarbon sequestration is not well supported by this work. The authors need to be much more careful\nabout inferring changes in carbon sequestration from TE alone; not overinterpreting the model w.r.t.\na more complex reality; and consider the differences in future projection timescales (100 years) to\nequilibrium changes presented here. I think the paleo-perspective on future changes is incredibly\nvaluable in general, but it has to be properly and carefully caveated due to the differences in the\nresolution of timescales, spatial scales, and processes. Unfortunately, I don’t think latter part of the\nmanuscript achieves this.\nThanks for your positivity comments. We greatly appreciate your insights and recognize the\nimportance of carefully framing our interpretations and implications. Below, we address your key\nconcerns and outline the revisions we have made (details see the responds to each comment).\n1)The section of “Increased POC transfer efficiency after the rise of modern phytoplankton” and\nrelative Figure 3d have been deleted.\n2) Some strong statements about implications for changes in carbon sequestration and future BCP\nchange have been weakened or removed.\n3) We have performed another three multi-factor experiments varying both [PO ] and temperature\n4\ndependence to evaluate the combined effects by these competitive factors.\n4) We have revised the texts to point out some limitations in cGENIE model.\nGeneral Comments\n∂13C gradient interpretation: There needs to be some clarification on how vertical ∂13C gradients\nare being interpreted. An increase in production, without a change in TE, can lead to an increased\nvertical ∂13C gradient (if measured at a single deep horizon) because both the surface and deep\nvalues are changing. A change in TE will have a similar effect but also changes the distribution over\na range of depth levels. Effectively you show this in Figure 2 g-i as the difference between blue bars!\nIf only one depth level is used then I’m concerned you cannot discern the difference. I was unsure\nif the 0 permil reference standard (lines 105 – 106) are an attempt to control for the production\nscenario? If so, this needs to be much more clearly linked. This is not a major problem because\nFigure 2 shows the production change to empirically to be a smaller effect than with TE. However,\nthis is not clear in the discussion of interpreting vertical gradients. Additionally, there are citations\nthat explore vertical ∂13C gradient interpretations that could be cited here to help such as Morée et\nal., (2018, 2021).\nIn the “eutrophication” scenario, doubling nutrient levels enhance primary production, which can\nresult in a greater vertical δ13C gradient, as shown by the blue bars in Figure 2g–i. We have added\nexpanded discussion of this interpretation. Lines 194-197: “The “eutrophication” scenario\ndemonstrates that doubled nutrient levels stimulate productivity (Supplementary Fig. 13),\nconcurrently increasing the vertical δ13C gradient. However, this effect remains substantially\nDIC\nweaker than the gradient increase induced by the “metabolism” scenario (Fig. 2)”.\nDifferent depth levels had been applied to calculate the gradient to make the conclusion robust,\nincluding δ13C difference between surface water (~40 m) and mid-water (~130 m) showed in Fig.\n2, and δ13C difference between ~40 m and ~230 m (Supplementary Fig. 5). We have added the\ncitation of Supplementary Fig. 5.\nThe method zeroing surface δ13C was consistently applied to both palaeoproxy data and model\noutputs to enable direct comparison. We have clarified it (Lines 110-111): \"Surface water δ13C\nDIC\nvalues were normalized to zero for both proxy records and model scenarios to facilitate consistent\nintercomparisons.\"\nThank you for recommending these relative literatures, which significantly strengthened our\nmechanistic interpretation. Relevant citations have been incorporated throughout the introduction\nand discussion .\n∂13C errors: there are notable uncertainties with the ∂13C values from cGENIE (Figure 2), mainly\nat depth. The authors have subsampled the latitudinal regions most represented by the observations,\nbut does sampling model ∂13C from the nearest location reduce the uncertainty? I this this\nessentially asks if the sampling across shelf vs. open ocean environments in cGENIE adds to the\npresented uncertainty?\nThank you for your suggestions. In the previous manuscript, we had sampled modeled δ¹³C from\nthe nearest locations corresponding to each study region (Supplementary Figs. 2–5). Based on these\nestimates, we recalculated the average gradients and obtained steeper gradients with larger\nuncertainties (Figure R1). Distinguishing between shelf and open-ocean environments did not\nappear to effectively reduce uncertainty, despite slight differences in values. This is likely due to the\nuneven spatial distribution of sampling locations. Furthermore, the relatively low spatial resolution\nof the cGENIE model limits its ability to accurately resolve exact locations, often leading to\nconfusion between shelf and open-ocean environments. Therefore, we have decided to retain our\nprevious approach of sampling model results based on latitudinal regions, as this helps avoid issues\narising from spatial unevenness.\nIn response to the Reviewer, we have revised the Methods section to clarify this approach. Lines\n438-439:“As sampling model results based on latitudinal regions helps avoid issues arising from\nspatial unevenness”.\nFigure R1. Comparison of proxy and simulated vertical DIC δ13C gradients in the upper\nocean during the three hyperthermals. The average of simulated δ13C (mean ± 1 s.d.) is\nDIC\ncalculated based on model results from nearest locations.\nMulti-factor experiments: much of the later text in the manuscript discusses net effects of changing\nexport and TE (see lines 238 – 241). This manuscript only presents single factor experiments which\ndon’t directly quantify this. Multi-factor experiments, e.g., varying both [PO4] and temperature\ndependence, might be illuminating here. I don’t see a reason why it should be assumed that these\ninteractions are linearly additive and can be assumed from the current outputs.\nWe have performed another three multi-factor experiments varying both [PO ] and temperature\n4\ndependence to evaluate the combined effects by these competitive factors. Model results reveal that\nthese competing drivers collectively reduce mid-water POC flux during the PTTM, while enhancing\nfluxes during the the ETTM and PETM, despite maintaining consistently reduced POC transfer\nefficiency across all events (Supplementary Fig. 18).\nWe have expanded the discussion through new multi-factor experiments. Lines 253-260: “We\nfurther perform “multi-factor” experiments varying [PO ] and temperature dependency\n4\nsimultaneously to investigate the interplay between nutrient-dependent productivity and\ntemperature-dependent microbial metabolism on the BCP. Model results reveal that these competing\ndrivers collectively reduce mid-water POC flux during the PTTM, while enhancing fluxes during\nthe the ETTM and PETM, despite maintaining consistently reduced POC transfer efficiency across\nall events (Supplementary Fig. 18). This inverse mid-water POC flux pattern may account for the\ncontrasting OC burial responses observed among different hyperthermals”.\nThe limitations and future directions are clarified. Lines 269-273: “However, our simulated mid-\nwater POC fluxes should not be conflated with OC burial, critical preservation processes of OC are\nnot involved in our experiments. A systematic integration of proxy data (e.g. TOC) with models\nexplicitly incorporating OC burial processes will be essential to elucidate how these warming-\ninduced competitive factors control OC burial across different hyperthermals”.\nSupplementary Fig. 18. Simulated changes in biological carbon pump for metabolism and multi-\nfactor scenarios (varying both [PO ] and temperature dependence) during the three hyperthermal\n4\nevents. (a) Global mean particulate organic carbon (POC) exported at surface water (40 m). (b)\nGlobal mean POC flux at mid-water (1040 m). (c) POC transfer efficiency (i.e., the ratio of POC at\nthe 1040 m to those at the surface water).\nSinking rates: I really struggled with the discussion of the sinking rates to explain the difference in\npre-CIE TEs for PTTM and ETTM in Figure 3. As I understand from the manuscript, a reduced\nsinking rate was applied in the first instance for the PTTM experiment to reflect the absence of\nmodern phytoplankton (lines 152 - 156). The text around lines 250 – 259 seems to then try to\ninterpret the difference in TE despite similar temperatures and concludes it’s a reduced sinking\nrate…which seems circular! Are the authors trying to demonstrate that a reduced sinking rate is\nrequired to help explain the ∂13C gradient? In which case, this is not coming across and would need\nthe equivalent of Figure 2g but from the sensitivity experiments. In general, the “Increased POC\ntransfer efficiency after the rise of modern phytoplankton” section is predicated on the applied\nchange in sinking rate being representative and there is no discussion of this. The temperature-\ndependence findings emerge from the underlying multiple constraints represented by the model.\nThese are therefore much stronger than an applied change in a parameter that directly leads to a\npredictable outcome. See also the specific comment on lines 296 – 299 below.\nAfter thorough evaluation, we acknowledge the limitations in the original analysis of sinking rates\non transfer efficiency. Consequently, we have removed the section titled “Increased POC transfer\nefficiency after the rise of modern phytoplankton” and associated Figure 3d, based on two critical\nconsiderations: 1) single parameter of sinking rate does not consider the various abundances,\ndiameters, and diversities of different plankton community, which could not generate robust\nmodeling results that support the statement in this section. 2) The deleted section constituted\nsupplementary speculation rather than model-supported findings, and its removal does not impact\nour core conclusions regarding temperature-dependent microbial metabolism.\nInstead, we added brief discussion state that sinking rate has slight influence on the gradient based\non sensitivity experiments (Supplementary Fig. 14). Base on this, we revise the statement and\nfigures as follows:\n1)The section of “Increased POC transfer efficiency after the rise of modern phytoplankton” and\nrelative Figure 3d have been deleted.\n2) We added discussions on the effect of the sinking rate on gradients based on sensitivity\nexperiments. Lines 199-202: “Additionally, sensitivity experiments reveal that POC sinking rate in\nthe cGENIE model exhibit limited influence on the gradients (Supplementary Fig. 14), reinforcing\nthe primacy of temperature- dependent microbial metabolism”.\n3) All strong claims about sinking-rate effects in the “Implications for future climate and the BCP”\nsection have been removed or revised (see detailed responds to the specific comments).\nSpecific Comments\nLines 24 – 25: “leading to diminished organic carbon sequestration in the ocean interior”. The\nmanuscript infers this but does not quantify it so you can explicitly say this.\nThere is only quantified transfer efficiency (TE) in our study and no modelled results about how\nreduced TE exacerbated the greenhouse effect actually. We revised text to clarify TE’s potential\nfeedback role rather than asserting direct greenhouse effects, and delete all strong statements.\nLines 24-25: “potentially leading to diminished organic carbon sequestration in the ocean interior”. Also\nsee the consistent revision about comments of Lines 222 – 227.\nLines 42 – 46: Sequestration occurs across the whole water column §and is somewhat evenly\ncontributed from remineralisation above and below 1000m (Siegel et al., 2023). Siegel et al., (2023)\ndemonstrates that you could derive a metric of remineralization weighted mean first-passage time\nto reflect this in one quantity…but…the first-passage time is not easily quantifiable! Instead transfer\nefficiency is an approximate metric of where remineralisation is occurring across the water-column.\nThis is somewhat pedantic but I think an important point in the context of later discussions of TE\nand carbon storage.\nWe thank the reviewer for highlighting this critical distinction. We have revised the text to clarify\nthat transfer efficiency (TE) serves as an approximate proxy for the efficiency of biological carbon\nsequestration. Relative citations are added (e.g. Siegel et al., 2023).\nLines 41-48: “This downward transport process stores carbon within the ocean interior where\ncarbon can be sequestered on centennial and longer timescales17, 18. Transfer efficiency, defined as\nthe fraction of POC transported from the surface ocean (0−100 m) to the ocean interior, serves as\na widely used proxy for estimating the biological carbon pump sequestration. Despite transfer\nefficiency represents an operational approximation rather than a direct sequestration measurement,\nit provides critical insights into natural carbon sequestration process in the ocean.”.\nLines 60 – 63: the citation Kwon et al. (2009) does not comment on vertical gradients of ∂13C\nThanks, correction made. We have cited another relevant paper (Kroopnick et al., 1985) here, which\ndescribes the distributions of δ13C and their vertical gradients.\nKroopnick, P.M., 1985. The distribution of 13C of ΣCO in the world oceans. Deep Sea Research\n2\nPart A. Oceanographic Research Papers 32, 57-84.\nSupplementary Figures 2 to 4: It’s not clear which symbols are data and which are model.\nThe symbols have been added to Supplementary Figures 2 to 4, consistent with those in Figure 2.\nLines 137 – 138: see general comment above on productivity impact on gradients – is this what is\nmeant by ocean nutrient conditions?\nIn the “eutrophication” scenario, doubling nutrient levels enhances production, leading to a greater\nvertical δ13C gradient. Relevant revisions can be found in the responses to the general comments.\nLine 142: The “stagnation” terms implies a slow down in overturning– is this the actual case in your\nexperiments? The supplementary figure 10 suggests slightly higher export which to me suggests\nmaybe actually greater overturning. Perhaps the name is not quite appropriate and something like\n“continental and climate” might be more representative.\nWhat you think is correct. The slightly greater overturning causes a slightly larger export in the\n“ocean stagnation” scenario. To avoid this misunderstanding, we have replaced “ocean stagnation”\nwith “circulation change” throughout the manuscript.\nLines 147 – 148: you need to be more explicit about what a “static biological pump” is. This could\nbe construed as fixed export production (not responding to a change in nutrient supply) which I\ndon’t think is the case. I think here this means fixed attenuation of POC with depth only.\nWhat you think is correct. The “static biological pump” refers to a fixed attenuation of POC with\ndepth only. We have clarified it as follows: Lines 167-168: “static biological pump characterized\nby a fixed attenuation of the POC depth profile”.\nLines 193 – 195: “aided by warming-induced increases in photosynthetic rates” – is this in cGENIE\nor a wider statement? You need to be clear whether any temperaturedependence applies to export\nproduction in the model (Crichton et al., 2021 for example introduced this, a temperature term in\nthe production equations exists in earlier productivity schemes). In general, the increased POC\nexport with lower TE is well-known (see Kwon et al., 2009 for example). It might be interesting to\nseparate out the impact of vertical redistribution and export on ∂13C gradients though! (I think this\nis beyond this manuscript, so just an interested comment).\nThank you for raising these important points. In our study, the temperature-dependent\nphotosynthetic rate was incorporated into the “metabolism” scenario using a dynamic model,\nwhereas the static model used in the “circulation change” and “eutrophication” scenarios does not\naccount for temperature dependency. Thus, there is only a negligible change in productivity in the\n“circulation change” scenario (Supplementary Fig. 13). These two approaches align with the STND\nand TDEP schemes proposed by Crichton et al. (2021).\nWe have clarified this in the revised manuscript. Lines 217-221: “This effect is further enhanced by\nwarming-induced increases in photosynthetic rates, leading to higher productivity (Fig. 4). The\n“metabolism” scenario incorporates the temperature dependency of the nutrient uptake rate in the\nmodel45. In contrast, the static model used in the “circulation change” and “eutrophication”\nscenarios does not account for temperature dependency (Supplementary Fig. 13)”.\nYour comment on “separating out the impact of vertical redistribution and export on δ13C gradients”\nis truly inspiring. We will explore this in future work.\nLines 207 – 208: This is fairly unsurprising because you have a fixed sinking rate and a\nremineralisation rate that is varying with local temperature in cGENIE. Since temperature is\nmonotonically decreasing with depth and there are no other processes, such as zooplankton\nmigration, you essentially have a power-law curve (Cael & Bisson 2018). In reality this may not be\nthe case if there are processes that are not monotonic with depth I think.\nThank you for raising these important points. We have added some discussions to clarify the\nlimitations of our model about the non-monotonic processes (e.g., zooplankton migration): “without\nconsidering zooplankton migration” (line 242). More relevant revisions can be found in the next\nresponse.\nLines 210 – 216: There are some details here that feel overstated or not fully stated. 1) The\nindependence of global mean [PO4] and TE is totally unsurprising! At least in cGENIE, there is no\nmechanism linking the distribution of remineralisation rates across depth with production. This is\nnot necessarily always the case, e.g., a sinking rate dependent on ecosystem size structure which\nitself is some partial function of total nutrient supply. 2) The latitudinal distribution of TE in cGENIE\nis solely dependent on temperature remineralisation as it’s the only process resolved for POC that\nleads to spatial variability. TE variability simply falls out of the latitudinal distribution of\ntemperature. 3) DeVries & Weber (2017) is one estimate of TE but it doesn’t agree with all\nobservations, such as deep sediment traps (Henson et al., 2012)!\nWe thank for your insightful comments. We have revised the texts to weaken claims and point out\nsome limitations in cGENIE model on Lines 234-242: “Along with our sensitivity analyses of ocean\ntemperature, the independence of transfer efficiency from [PO ] (Supplementary Fig. 16) and\n4\nspatial pattern of higher transfer efficiency in the high latitudes in our simulations (Supplementary\nFig. 17) further support the idea that transfer efficiency is inversely correlated with ocean\ntemperature. Our spatial pattern is consistent with the results of data-assimilated biological pump\nmodels and modern field observations31, though quantitative discrepancies persist with direct\nsediment trap measurements20. This discrepancy potentially arises from the sole temperature-\ndependent remineralization process incorporated in our model, without considering zooplankton\nmigration46 or size-dependent plankton communities47”.\nLines 222 – 227: The claim that reduced TE exacerbated the greenhouse effect and was an important\ncontributor to persistent warming is dependent on how big this effect actually was but is not\nquantified here! Is there any way to quantify this relative to the other carbon cycle changes? If not,\nthese statements need to be toned down towards TE having a potential role.\nThere is only quantified TE in our study and no modelled results about how much reduced TE\nexacerbated the greenhouse effect actually. We revised this part to clarify the TE only having a\npotential role, and delete some strong statements on lines 248-252: “However, a reduction in POC\ntransfer efficiency due to a warming-induced increase in microbial metabolism—as indicated by\nour model results (Figs. 3 and 4)—potentially lead to less CO sequestration in the ocean interior.\n2\nThis process may create a potential positive feedback in amplifying the greenhouse effect during the\nthree studied hyperthermal events, which persisted for hundreds of thousands to millions of years\n(Fig. 1).”\nLines 250 – 259: SST is not really the best metric because remineralisation is driven by in-situ\ntemperature. The mean temperature across the ocean volume might be more representative.\nWe agree that the mean temperature across the ocean volume might be more representative.\nAnyway, we have removed this section. For more details, please refer to our response to the general\ncomments on sinking rates.\nLines 296 – 299: This conclusion cannot be drawn from this study. You applied a change in the\nsinking rate based on the assumption of plankton community changes. That can reproduced by\nsimply taking an empirical flux curve representing the flux of POC with depth (e.g., Cael & Bisson\n2018) and changing the sinking rate with the same assumption. To support this statement, you’d\nneed to show that the plankton community is driving TE independently of a parameter choice you\nmade and which emerges from your model simulations. In other words, you got out what you put\nin, which doesn’t support this conclusion.\nWe agree what you suggest. We have deleted these strong statements about phytoplankton community\nand focused on impact of the temperature dependent remineralization. For more details, please refer\nto our response to the general comments on sinking rates.\nLines 300 – 302: this is far too strong a statement to be supported by this manuscript. The manuscript\nonly considers temperature dependent remineralisation, whereas there are a wide range of other\nprocesses such as oxygen dependence; ballasting (protection via minerals); organic matter lability.\nSee also above comment on plankton community structure and sinking velocities.\nWe have deleted the overstated statements about phytoplankton community and emphasized a wide\nrange of other processes might influence the transfer efficiency on lines 292-295: “Furthermore,\nthis study does not account for other factors, including trait-based functional types of the plankton\ncommunity45, oxygen-dependent processes51, and mineral ballasting effects52, all of which may\naffect transfer efficiency”.\nLines 302 – 308: Again, these statements are too strong and not consistent with the findings of the\nmanuscript. A ∆Carbon storage (by the biological pump) is proportional to ∆TE only in isolation\nand is also dependent on timescale. Export production and the residence time of regenerated carbon\nare changing simultaneously in projected future changes leading to a net increase in carbon storage\nby 2100 despite the general decrease in TE (Wilson et al., 2022). The increase occurs because of\nthe centennial timescale where the decreases in export and TE haven’t had time to propagate through\nthe ocean circulation. This change in carbon storage will differ over longer timescales (>2100) as\nthose competing physical and biological signals are convolved, but, the magnitude and direction is\nan open question which I don’t think this manuscript is able to provide a clear perspective on.\nThank you for highlighting the importance of timescales in carbon sequestration. We agree that our\noriginal statements were overstated. We have revised the text to clarify the importance of timescales\nin carbon storage and to emphasize some model limitations. These revisions better align our\nconclusions with the scope of the manuscript.\nLines 290-298: “However, this positive feedback on carbon sequestration depends on timescales.\nModel projections suggest a net increase in oceanic carbon storage by 2100 due to short timescales\n(<100 years), despite a decrease in transfer efficiency21. Furthermore, this study does not account\nfor other factors, including trait-based functional types of the plankton community45, oxygen-\ndependent processes51, and mineral ballasting effects52, all of which may affect transfer efficiency.\nNevertheless, our integrated approach combining proxy data compilation with mechanistic\nmodeling highlights that temperature- dependent microbial metabolism may play a significant role\nin modulating the BCP dynamics during both ancient hyperthermal events and future climate\nscenarios”.\nPoint-by-point response to the reviewers’ comments.\nThose changes appear as highlights in red colour in the manuscript. Below are our responses (blue\ncolour) to reviewers’ comments. Line numbers in our response refer to the revised manuscript.\nREVIEWER COMMENTS\nReviewer #1 (Remarks to the Author):\nThe authors have addressed all the issues that I mentioned in the first round of review. I have no\nfurther comments on the ms and recommend acceptance, and hope to see the paper in print soon.\nWe appreciate your positive comments for this manuscript.\nReviewer #3 (Remarks to the Author):\nThank you to the authors for addressing my previous comments in detail and engaging in a\nconstructive discussion. I am happy that the authors have addressed the general comments raised on\nthe multi-factor experiments; the rise of modern plankton; and the strength of some statements.\nOverall, the manuscript is more focused and stronger in its findings because of this.\nWe sincerely thank the reviewer for taking time to evaluate our manuscript and positive feedback\non our previous revision.\nThe final future climate section now feels more appropriate as a discussion-level section. The\nauthors have added some clarifications on the complexities of considering future changes based on\nmy previous comments. However, I am still a bit hesitant on this part of the manuscript because by\nfocusing heavily on implications for the future BCP, the authors are predominantly relying on\nfindings in previous papers and are not informing the discussion with the results presented in this\nmanuscript. This is not totally detrimental to the manuscript but is a lost opportunity to stimulate\nfurther discussion in this area!\nThank you for your positive comments on our revisions. We agree that the final section lacks some\ndiscussions in relation to the results presented in this manuscript. In response, we have made detailed\nrevisions as outlined in the following reply, based on your suggestions.\nIn my view, the strength and novelty of the results presented by the authors is that they strongly\nconfirm that a temperature-dependent biological pump is a notable feedback in the carbon cycle and\nEarth System. Previous work by John et al., (2013) and John et al., (2014) showed its impact over\nthe long-term geological background but it’s role as a feedback is an assumption based on this – this\nmanuscript confirms we need to include it to help explain features of transient warming events. The\ndiscussion-level questions that follow in my mind are: what’s the size or significance of this\nfeedback; what timescale is it relevant over; is it sufficiently resolved in Earth System models?\nThe future is a part of this discussion but not the sole focus.\nThank you for your constructive suggestions on the final section. They were very helpful in\nimproving the manuscript. In response, we have added further discussion to emphasize the\nimportance of temperature-dependent microbial metabolism. Specifically, we now summarize the\nmagnitude of change in transfer efficiency resulting from a 4–10 °C temperature increase, and we\nalso clarify the relevant timescales associated with this process. In addition, a comparison of transfer\nefficiency across different Earth system models is available in the previous version (lines 301-305).\nThe title of this section has also been revised to “Implications for BCP change and future climate”\nto better reflect the content and its relevance to the results.\nLines 285-298: “Our findings highlight that temperature-dependent microbial metabolism may\nsignificantly influence BCP change and thus climate regulation. Specifically, a temperature increase\nof 4–10℃ during three ancient hyperthermals could reduce transfer efficiency of BCP by 2–6\ntimes, thereby lowering the flux of POC into the ocean interior and potentially counteracting OC\nburial (Fig. 3). This feedback likely resulted in reduced CO sequestration in the ocean interior,\n2\namplifying ancient warming (Fig. 4). Conversely, a cooling-induced decline in microbial\nmetabolism enhanced the efficiency of BCP over the past 15 Myr34. It should be noted, however, this\ncritical feedback becomes more effective over longer timescales, as organic carbon burial operates\non a million-year timescale 58. Consequently, the amplifying effect during the PETM may have been\nrelatively limited due to its shorter duration compared to the PTTM and ETTM (Fig. 1), which\nrestricted the establishment of sustained positive feedbacks. Furthermore, a key limitation of our\nsimulations is the uncertainty in how the 2–6-fold reduction in transfer efficiency decreases the\namount of CO sequestered in the ocean by the BCP for a geological timescale, which should be\n2\naddressed in future work.”.\nLine 45 – 48: I think the preceding statement about TE being a proxy for sequestration is sufficient\nhere and responds well to points raised in the previous review. Therefore, these following lines don’t\nseem as necessary here.\nMany thanks. We have deleted these lines to avoid redundancy.\nLines 142 – 143: “assume thousand deepening or shallowing scenarios” – just needs clarification as\nthis is ambiguous as written. It could be interpreted as 1000m or 1000 different scenarios.\nThank you for pointing out this ambiguity. We agree that the original phrasing was unclear. To clarify,\nthe sentence refers to 1000 different scenarios of deepening or shallowing for each hyperthermal\nevent. We have revised the sentence in the manuscript accordingly. Lines 140-142: “…we assume\n1000 different scenarios of deepening or shallowing for each hyperthermal event to assess the\nvertical gradient responses to depth changes during the event (see Methods)”.\nLine 166: reference to “ocean stagnation” scenario left over from previous manuscript version.\nCorrection made. The“ocean stagnation” is corrected as “circulation change”\nLines 167 – 168: I think it’s worth being explicit by noting that whilst the remineralisation profile\nis fixed export production is still dynamic in response to surface ocean nutrient concentrations.\nYes, export production of static model remains dynamic and responds to surface ocean nutrient\nconcentrations. We have revised relevant sentence to explicitly state that. Lines 167-168:\n“…whereas export production could still vary in response to surface nutrient concentrations”.\nLines 168 – 169: You should state the metabolic scenario includes both temperature-dependence on\nproduction and remineralisation. This is noted later on in the manuscript but it’s of particular\nrelevance here in my opinion.\nWe agree this. We have added relevant description on lines 169-171: “…higher temperatures raise\nthe microbial metabolic rate, thereby altering photosynthesis and POM remineralization depth\nprofile”.\nLines 170 – 175: Apologies if I missed this previously – is the reduced sinking rate applied in all\nscenarios or just the metabolic scenario? This may be open to interpretation, but I would see this as\na long-term background feature similar to a warmer climate/higher CO2 on top of which the\nhyperthermal occurs. You seem to have demonstrated this is tested for in lines 200 – 202 but again\nit’s worth being explicit about the experiment design.\nWe apologize for our mistake. The reduced sinking rate is only applied in metabolism scenario. We\nhave revised relative text in discussion and methods to clarify this.\nLines 170-171: “In the “microbial metabolism” scenario, we reduced the POM sinking rate by ~22%\nin the PTTM experiments…”\nLines 201-202: “…sensitivity experiments under the “microbial metabolism” scenario reveal…”\nLines 444-445: “…as for “microbial metabolism” scenario, we set the default sinking rate…”.\nLine 200: is “reveal that the reduced POC sinking rate” correct, if so I think this is clearer.\nCorrection made.\nReviewer #4 (Remarks to the Author):\nThe authors used available carbon isotopes values from the literature and numerical simulations to\nreconstructed the response of the biological pump during in 3 past intervals of warming. This\nresponse is presented to be affected by a reduced vertical organic carbon transfer efficiency as a\nresponse to warming, owing to a positive temperature-dependent microbial remineralization in the\nocean interior.\nWe sincerely thank the reviewer for their thorough and constructive evaluation of our manuscript.\nThe comments provided valuable perspectives that helped us to clarify key concepts. The\nsuggestions concerning the conceptual framing and the treatment of proxy data were especially\nhelpful. We have revised the manuscript accordingly, and a detailed, point-by-point response to each\nof your comments is provided below.\nWith regards to the conceptualization:\nThe reduction of vertical fluxes with warming would not in itself result in a biological pump that is\nless capable of storing carbon in the ocean interior, if we focus only in the ocean-atmosphere system.\nThat is because this alone would not affect the ratio of preformed over remineralized nutrients in\nthe ocean interior. What would happen with a weaker vertical flux is that more carbon and nutrients\nare remineralized in shallower depth, but ultimately these nutrients now residing in the upper water\ncolumn would be eventually consumed in an ocean that is mostly oligotrophic (if they are\nremineralized in the uppermost layers), restoring the flux to its original value.\nWe agree that increased remineralization of carbon and nutrients at shallower depths facilitates their\nreturn to the surface. Since phosphate is the only limiting nutrient in our simulation, the recycled\nphosphate can be efficiently utilized to enhance primary productivity. This increase in productivity\nindeed results in enhanced surface export production fluxes, as simulated by cGENIE (Fig. 3a).\nSome relative discussions are available in previous version (Lines 215-223).\nHowever, our results indicate that this compensatory flux is insufficient to fully offset the reduction\nin mid-water POC flux caused by warming-enhanced microbial remineralization (Fig. 3b).\nSpecifically, warming substantially increases microbial respiration rates throughout the upper and\nmid-water columns. The temperature-dependent increase in remineralization rates is non-linear and\nappears to dominate mid-water POC fluxes in our simulations. Accordingly, our findings highlight\nthat temperature-dependent microbial metabolism may significantly influence BCP change and thus\nclimate regulation. We have revised this part to make it clearly.\nLines 236-239: “Despite the increased POC export production in the surface water, our results\nindicate that this compensatory flux is insufficient to fully offset the reduction in mid-water POC\nflux due to warming-induced increases in microbial respiration, leading to a significant decrease in\nPOC transfer efficiency (Fig. 3b, c)”.\nThis can be conceptualized with a 2-box model, where if the upper ocean box is kept oligotrophic\nat steady state (a fair assumption for most of the ocean). Sarmiento & Toggweiler showed in the\n1980s that adding a third box, with partial nutrient consumption in high latitude regions (where\nnitrate/phosphate > 0) does have a real leverage on atmospheric.\nIf productivity/export changes in high-latitude regions (that ventilate the ocean interior, and thus\nchange its preformed nutrient content and thus carbon storage), the mechanism proposed here could\nvery well work. But the idea that warming alone per se causes the biological pump to represent a\npositive feedback on the initial CO2 release is not very compelling to me.\nWe appreciate the reviewer’s insightful comment on the classical view that high productivity\nenhances the oceanic carbon sink. For example, Sarmiento and Toggweiler (1984) showed that a 70\nppmv increase in CO at the end of the last ice age was primarily driven by a decrease in net high-\n2\nlatitude productivity, possibly coupled with change in thermohaline overturning. However, some\nresearchers have argued that quantifying productivity or carbon export alone is insufficient to assess\nthe biological pump’s impact on atmospheric CO , as its efficiency also depends on the residence\n2\ntime of exported carbon in the ocean (e.g. DeVries et al., 2012; Siegel et al., 2023). Carbon in the\nupper ocean typically returns to the atmosphere on sub-centennial timescales, limiting its\ncontribution to long-term CO storage. In contrast, carbon transported to the deep ocean or buried\n2\nin sediments can remain sequestered for significantly longer periods (Siegel et al., 2023). As noted\nin your following comment, temperature-dependent microbial metabolism can influence the rate of\norganic matter burial. Accordingly, in the context of ancient hyperthermals, the transfer efficiency\nof POC to the deep ocean, rather than primary production, is a more appropriate metric for assessing\nthe long-term sequestration capacity of the biological carbon pump.\nTo clarify this issue, we have added discussions highlighting the importance of high productivity in\nenhancing the oceanic carbon sink, as shown in previous studies (e.g., Sarmiento and Toggweiler,\n1984). Furthermore, we clarify that in our study, transfer efficiency, rather than export production,\nis a more appropriate metric for estimating the long-term sequestration capacity of the biological\ncarbon pump on geological timescales.\nLines 226-232: “Productivity or export production is an important metric for evaluating the\nsequestration capacity of the biological carbon pump15, 48. However, quantifying export production\nalone is insufficient to assess the biological pump’s impact on atmospheric CO , as its efficiency\n2\nalso depends on the residence time of exported carbon in the ocean17-19. Specifically, carbon in the\nupper ocean (<1,000 m) typically returns to the atmosphere on sub-centennial timescales, limiting\nits contribution to long-term CO storage, whereas carbon transported to the deep ocean (>1,000\n2\nm) or buried in sediments can remain sequestered for much longer periods.”\nWe acknowledge that transfer efficiency is an indirect proxy, as our simulations cannot directly\ncalculate how a 2–6-fold reduction in transfer efficiency decreases the amount of CO sequestered\n2\nin the ocean by the BCP for a geological timescale. Because CO₂ is used as a forcing to simulate\ntemperature rise in our model, it cannot be modified by feedback processes within the simulation.\nWe have clearly stated this limitation of our study on lines 295-298: “Furthermore, a key limitation\nof our simulations is the uncertainty in how the 2–6-fold reduction in transfer efficiency decreases\nthe amount of CO sequestered in the ocean by the BCP for a geological timescale, which should\n2\nbe addressed in future work.”.\nIf we are talking about the implication of metabolic theory on rates of organic matter burial in the\nseafloor, then what the authors discuss could be relevant for atmospheric CO2. But even then, such\nmechanism could mostly relevant for the Permo-Trias and the Early Toarcian (Million years\ntimescales), I am not sure if it would be relevant on the timescales of the PETM, which onset is\nestimated to have lasted from 2000-20000 years. Carbon burial today, in a much colder ocean\nbaseline, is already pretty small (0.15 Pg C per year, Kandasamy & Nath, 2016, Table 3), in\ncomparison to the 3000 PgC of carbon expected to have been released during the PETM (Zeebe et\nal., 2009).\nYes, we agree that temperature-dependent remineralization can potentially influence the rate of\norganic matter burial, a critical control that has often been overlooked in previous studies which\nprimarily focused on factors such as redox conditions. We also recognize the importance of the\ntimescale over which this feedback develops and manifests. Accordingly, we have clarified in the\nrevised manuscript that temperature-dependent remineralization likely had a limited impact during\nthe PETM due to the relatively short duration of this event.\nLines 285-295: “Our findings highlight that temperature-dependent microbial metabolism may\nsignificantly influence BCP change and thus climate regulation. Specifically, a temperature increase\nof 4–10℃ during three ancient hyperthermals could reduce transfer efficiency of BCP by 2–6\ntimes, thereby lowering the flux of POC into the ocean interior and potentially counteracting OC\nburial (Fig. 3). This feedback likely resulted in reduced CO sequestration in the ocean interior,\n2\namplifying ancient warming (Fig. 4). Conversely, a cooling-induced decline in microbial\nmetabolism enhanced the efficiency of BCP over the past 15 Myr34. It should be noted, however, this\ncritical feedback becomes more effective over longer timescales, as organic carbon burial operates\non a million-year timescale 58. Consequently, the amplifying effect during the PETM may have been\nrelatively limited due to its shorter duration compared to the PTTM and ETTM (Fig. 1), which\nrestricted the establishment of sustained positive feedbacks”.\nExport production proxies are unfortunately ridden with caveats and cannot help us much. For\nexample, there is no commonly accepted evidence that export fluxes declined during the PETM. A\nlot of superficial assessments have been made using Ba fluxes being higher during peak warmth,\nwhich are now been shown to be artifacts of variable sedimentation rates (Torfstein et al., 2010).\nBarium isotopes (flux independent) show the diametric opposite trend to Ba fluxes (Bridgestock et\nal. 2019), which does not support\nThank you for share the conflicting interpretations of export production during the PETM based on\nbarium fluxes and barium isotopes. This discrepancy reflects a broader emerging debate regarding\nchanges in primary productivity during the PETM. Our cGENIE simulations indicate an increase in\nexport production, which is consistent with the trend suggested by barium isotope studies. We have\nnow cited Bridgestock et al. (2019) in the revised manuscript (line 213) to support this interpretation.\nWith regards to the observations:\nThe paper argue for bigger carbon isotope gradients during warming intervals. I did look at\npublished results from the PETM in Sites 527, 401 and 1209 for planktic and benthic and compiled\nthem (see picture attached to this review), The figure seems to show that:\n1) the vertical 13C gradient is smaller, not larger, during the peak-PETM.\n2) Most of the change is driven by planktic (upper ocean) drop in d13C, much less in the deep-\ndwelling benthic organisms (where enhanced remineralization should drive most of the signal).\nWhile this represent a more superficial exploration of the available records in comparison to the\nmore meticulous work done by the authors, I do find these discrepancies very puzzling.\nThank you for pointing out this issue. We also appreciate the well-compiled data and figures that\nyou provided. We have noted that the difference in δ13C between planktonic and deep-ocean benthic\nforaminifera becomes smaller during the PETM. In other words, the magnitude of the carbon isotope\nexcursion (CIE) in planktonic foraminifera is larger than that in deep-ocean benthic foraminifera, a\npattern also observed in previous studies (Zachos et al., 2007; McInerney and Wing, 2011). One\nexplanation is that the smaller CIE magnitude returned by pelagic sites (such as Sites 527, 401 and\n1209) might result from incomplete preservation stemming from carbonate dissolution, rendering\nthem unsuitable for reconstructing the vertical gradient. Furthermore, the strongest vertical δ13C\ngradient usually occurs in the upper ocean (>1500 m), whereas the gradient is less pronounced\nbelow this depth (Figure R1 from Dentith et al., 2020). As a result, deep-ocean benthic\nforaminifera are unsuitable for reconstructing the vertical gradient.\nFigure R1. Depth profiles of δ13C during the 1990s (Dentith et al., 2020)\nDIC\nAccordingly, in addition to the carbonate records, we used paired δ13C and δ18O data from\nplanktonic foraminifera to reconstruct the vertical δ13C gradient in the upper ocean during the PETM\nin previous version (Supplementary Figure 7 and 8), excluding benthic foraminiferal records. This\napproach follows previous study (Boscolo-Galazzo et al., 2021). Planktonic foraminifera records\nindicate an increasing trend in vertical gradients, from 0.65 ± 0.16‰/hm to 1.21 ± 0.41‰/hm, which\nfurther supports our findings. We have incorporated relevant discussion in the Methods section to\nclarify the limitations associated with benthic foraminiferal records. We hope it helps clarify this\nissue.\nLines 394-398: “Deep-ocean benthic foraminifera are excluded from this approach, as previous\nstudies have noted that the vertical δ¹³C gradient between planktonic and deep benthic foraminifera\ndecreasesd during the PETM64, 65. One possible explanation is that the smaller CIE magnitude\nobserved at pelagic sites may result from incomplete preservation due to carbonate dissolution65,\nrendering them unsuitable for reconstructing the vertical gradient”.\nI do appreciate the work done by authors in combining modelling and observations in the context\nof a very important topic such as the biological pump, evaluated across different warming intervals\nof Earth’s history. I personally would suggest that the aforementioned concerns should be addressed\nbefore this manuscript is published.\nWe really appreciate your positive comment on our works. We have tried to address all issues and\nimprove our manuscript. Thanks again.\nReferences\nBoscolo-Galazzo, F., Crichton, K.A., Ridgwell, A., Mawbey, E.M., Wade, B.S. and Pearson, P.N.,\n2021. Temperature controls carbon cycling and biological evolution in the ocean twilight zone.\nScience 371, 1148-1152.\nDentith, J.E., Ivanovic, R.F., Gregoire, L.J., Tindall, J.C. and Robinson, L.F., 2020. Simulating\nstable carbon isotopes in the ocean component of the FAMOUS general circulation model with\nMOSES1 (XOAVI). Geoscientific Model Development 13, 3529-3552.\nDeVries, T., Primeau, F. and Deutsch, C., 2012. The sequestration efficiency of the biological pump.\nGeophysical Research Letters 39.\nMcInerney, F.A. and Wing, S.L., 2011. The Paleocene-Eocene Thermal Maximum: A Perturbation\nof Carbon Cycle, Climate, and Biosphere with Implications for the Future. Annual Review of\nEarth and Planetary Sciences 39, 489-516.\nSarmiento, J.L. and Toggweiler, J.R., 1984. A new model for the role of the oceans"}, "author_response": "Point-by-point response to the reviewers’ comments.\nThose changes appear as highlights in red colour in the manuscript. Below are our responses (blue\ncolour) to reviewers’ comments. Line numbers in our response refer to the revised manuscript.\nREVIEWER COMMENTS\nReviewer #1 (Remarks to the Author):\nEnhanced microbial metabolism during hyperthermal events would have led to diminished Corg\nsequestration in the ocean interior amplifying global warming, but meanwhile during warm periods,\nenhanced stratification caused by slowered ocean circulation and promoted marine primary\nproductivity would have collectively led to expansion of OMZ and favored Corg burial. Thus, the\ncarbon cycling and related feedbacks during hyperthermal events remain unclear. The paper\nattempts to better understand the processes and feedback through data compilation and Earth system\nmodeling for the three representative hyperthermal events during the Phanerozoic. The paper\nfocuses on the role of biological C pump to the global warming, which present two-fold main\nprocesses through Earth system modeling, involving elevated rates of Corg remineralization and\nreduced transfer efficiency of POC from the surface ocean to the deep ocean. The paper also\nexplores relationship between evolution of phytoplankton and BCP and further projects into the\ncurrent and future global warming. This is a well-written, easy-reading manuscript but bears\nsignificant implications for the role of BCP under the global warming conditions, and for readers of\nbroad interests. I therefore recommend publication after minor revisions.\nWe appreciate your positive comments for this manuscript.\nLine 19: What Earth system modeling used in the present study. Would be good to explain it for\nmore information.\nThanks for your suggestion. We have pointed out that the Earth system model cGENIE is used in\nthis study: “Simulations with the Earth system model cGENIE…”.\nLine 52: delete relatively; little is little.\nCorrection made\nIncreased POC transfer efficiency after the rise of modern phytoplankton, but with evolution of\nvarious phytoplankton and various particle size etc would also impact the POC transfer efficiency.\nThus, this process may be illustrated for an individual event, i.e., before and after, regarding the\nmicrobial evolution, but when compared among events, it might be more complicated.\nWe appreciate the reviewer’s insightful comment that different phytoplankton types and particle\nsizes indeed influence POC transfer efficiency, similar to the concerns raised by Reviewer #3.\nFollowing Reviewer #3’s suggestions, we have removed related inferences, as our simulation does\nnot provide robust results to assess the impact of different phytoplankton types on transfer efficiency.\nFigure1. Explain the regression curves and shades in b\nWe have explained the regression curves and shades in the caption of figure. “The blue and red lines\nrepresent the LOESS fit curve for δ13C and SST, respectively, while light blue and red shaded area\nrepresent 68% confidence intervals”.\nFigure 2. Authors may squeeze a-f panels vertically to save some space given the large white space\nin these panels.\nThanks. We have reduced the space between panels in the Figure 2.\nSupplementary Fig. 2. What are the different symbols for? The eight regional areas? If so, please\nindicate it in the Figure caption. Also, how to recognize proxy and modeled data? Need to indicate\nit too. I guess the questions work for Suppl. Fig. 3 and Fig. 4, too.\nThe different symbols represent proxy data and three model scenarios. The legends of symbols have\nbeen added in the Supplementary Figures 2 to 4, which are same as those in the Figure 2.\nSupplementary Fig. 5. h-j, should be g-i?—the caption is also a bit confusing what does peak CIE\nvs pre-CIE mean? Does that mean peak value minus pre-CIE value? Need more explanation.\nCorrection made. “peak CIE vs pre-CIE” has been changed to “peak CIE minus pre-CIE”, and we\nadd more explanation.\nReviewer #2 (Remarks to the Author):\nYuyang Wu et al. with manuscript entitled “Microbial metabolism amplified warming in three\nPhanerozoic hyperthermal events”, address a very important topic of highly interdisciplinary\ninterest for both understanding deep time and future climate change. In particular, the authors\naddress if there were potential change in the mean ocean biological carbon pump (BCP) during three\nmajor global warming events in Earth history (Permian-Triassic, Toarcian and Palaeocene-Eocene\nevents).\nTo test for a change in BCP) the authors focus on carbon isotopes in carbonates from different\ninterpreted ancient water depths. Importantly, the results are then compared to different scenario\ntesting with a global intermediate complexity earth system mode. Based on the large data\ncompilation and model experiments the authors conclude that there likely was a significant\nreduction of BCP during the events.\nMany thanks. We feel that is a fair summary in-line with our intention.\nWhile I am vary sympathetic to the authors approach, method and results, there are three major\nproblems with the implicit assumptions in the interpretation of the observations:\n1) The authors assume they can infer paleowater depth of deposition from sedimentological\nobservations or petrology, while in fact their proxies used are only relative indicators in part related\nto climate oceanographic changes such as storm wave base, and thermal stratification. To put it\nanother way, a marl deposited at one water depth before an event could very well be deposited at a\nshallower or deeper depth during the event. This would compromise the approach.\nWe agree that uncertainty of water depth change should be evaluated, especially the shallower or\ndeeper depth during the event due to the limitations of our proxies. The interval we focused from\npre-CIE to peak CIE is usually short, lasting 0.05~0.2 Myr during three hyperthermal events\n(Burgess et al., 2014; Nordt et al., 2022; Cramwinckel et al., 2018). Thus, the depositional setting\nusually keep constant, but the water depth potentially fluctuates and become shallower or deeper.\nIt is difficult to quantify the change of water depth by proxies, but we could assume multiple possible\nscenarios and evaluate the impact of depth change on the estimated vertical gradient. We consider\ntwo extreme scenarios: (1) the shallow shelf deepening to a deep shelf (“depth increase” scenario)\nand (2) the deep shelf becoming shallower to a shallow shelf (“depth decrease” scenario).\nFor example, in the “depth increase” scenario, we randomly sample several shallow-shelf locations and\nreclassify them as deep-shelf locations at the peak CIE stage. The results of 2000 random sampling\niterations for each hyperthermal show that the vertical gradients still larger obviously at peak CIE\nin most scenarios, except for some shallower scenarios from early Toarcian (Supplementary Figs.\n10 and 11). These analyses suggest that even allowing for extreme scenarios of depositional setting\nchange, increased trend of vertical gradient during the warming, as indicated by the δ13C data, is\ncarb\nstill relatively robust.\nWe have added relevant sensitivity tests (Supplementary Fig. 9) and a description of the methods\n(lines 353–365) to demonstrate that our key conclusion remains robust despite variations in water\ndepth. Some sentences are added on lines 140-149: “Furthermore, palaeodepths of locations may\nhave fluctuated (shallowing or deepening) during the event. Despite these changes could not be\nreliably quantified through existing proxy data, we assume thousand deepening or shallowing\nscenarios during the event to assess vertical gradient responses to depth changes during the event\n(see Methods). Results demonstrate that vertical gradients during the peak CIE stage remained\nsignificantly enhanced in most scenarios (Supplementary Figs. 10 and 11), though gradient\nmagnitudes diminished, as more locations exhibited depth changes. Collectively, the amplified\nvertical δ13C gradients documented across widely distributed sites and multiple geological\nDIC\nintervals strongly suggest that local effects, including palaeodepth variations, mixed DIC sources,\nor diagenetic alterations, do not constitute the dominant control”.\nSupplementary Fig. 10. Effects of “depth increase” scenario on vertical carbon isotope gradients. We\nrandomly select some locations of shallow shelf deepening to a deep shelf at peak CIE. The results of\n1000 random sampling iterations show that the vertical gradients still larger obviously at peak CIE\nin most scenarios. Red lines represent LOESS fit through the data. (a) Permian-Triassic; (b) Early\nToarcian; (c) Paleocene-Eocene.\nSupplementary Fig. 11. Effects of “depth decrease” scenario on vertical carbon isotope gradients. We\nrandomly select some locations of deep shelf becoming shallower to a shallow shelf at peak CIE. The\nresults of 1000 random sampling iterations show that the vertical gradients still larger obviously at\npeak CIE in most scenarios. Red lines represent LOESS fit through the data. (a) Permian-Triassic;\n(b) Early Toarcian; (c) Palaeocene-Eocene.\n2) The authors assume the carbonate sampled at a given interpreted depth where also precipitated at\nthat depth. Carbonate mud or particles may be dispersed from shallow to deep depending on climate.\nAs far as I can tell the authors do not present or discuss indicators that can confirm that the d13C\nsignal at depth is not derived from shallow water and if this redistribution process could be related\nto the climate change during the events. Presumably the carbonate in marls interpreted as their deep\nshelf is coming from micritic carbonate, was this micrite originally produced in shallow water or\nlocally? It is well known that the carbonate sources in platform carbonate factory have quite variable\nd13C_carb composition. Also, could it be that platform carbonate mud from the platform was not\nas efficiently dispersed to greater depth during sea level rise, and therefor indirectly causing an\nincreased gradient in d13C.\nWe acknowledge that the potential mixing of deep and shallow water DIC may influence the\nestimates. Specifically, the observed carbonate δ13C at the deep shelf may reflect a mixing signal of\nboth deep and shallow water DIC (Shen et al., 2011). If this is the case, the actual DIC δ13C of deep\nwater is likely lower than the carbonate δ13C observed at the deep shelf, meaning that our previous\nresults may have underestimated the vertical gradient.\nTo assess the uncertainty in the mixing proportion, we randomly sampled 1,000 different\nmixing proportions of shallow water DIC (pr = 0 to 90%). Using a simple isotopic mixing model of\ntwo end-members, we then calculated the actual δ13C of deep water (δ13C ) and the vertical\nDIC-deep\ngradient.\nδ13C = pr*δ13C + (1−pr)* δ13C\ncarb-deep DIC-shallow DIC-deep\nδ13C = (δ13C −pr*δ13C )/(1−pr)\nDIC-deep carb-deep DIC-shallow\nThe δ13C data from both the shallow and deep shelves are available (e.g. δ13C = 1‰;\ncarb carb-shallow\nδ13C = −1‰), and assuming δ13C could represent the δ13C . According to the\ncarb-deep carb-shallow DIC-shallow\nequation, when pr =0.5, δ13C is −3‰ and DIC gradient is −4‰. First, regardless of the mixing\nDIC-deep\nproportion, the gradient during the peak stage is generally larger than that during the pre-CIE stage,\nand increases as the proportion increases (Supplementary Fig. 9). Even if the proportion differs\nbetween the two stages due to climate change or sea level rise, as suggested by the reviewers, our\nanalysis, which includes 1,000 × 1,000 scenarios (Supplementary Table 3), shows that in most cases,\nthe vertical gradient still increases. For example, 91% of the scenarios at the P-Tr boundary indicate\nan increased gradient (Supplementary Table 3). These analyses suggest that increased gradient\nduring the warming, as indicated by the δ13C data, is largely independent of the mixing proportion\ncarb\nbetween shallow and deep water DIC. Moreover, previous estimates assuming no mixing (pr = 0)\nrepresent a conservative estimate of the gradient increase during warming.\nWe have added some brief description of sensitivity analyses in the Methods (lines 353-357) and\ndetails in supplementary materials. Some sentences discussing the mixing proportion are added in\nmanuscript. Lines 134-140: “The δ13C values recorded in deep shelf and basin environments\ncarb\npotentially represent a mixture of deep and shallow water DIC sources44. To evaluate potential\nmixing effects on gradient amplification, we implemented a simple isotopic mixing model (see\nMethods). Sensitivity tests reveal that increased gradient during the warming, as indicated by the\nδ13C data, show limited dependence on mixing source (Supplementary Fig. 9). The estimates\ncarb\n(shown in Fig. 2) without regard to mixed sources provide conservative values, as the magnitude\nincrease when raising the mixing proportion (Supplementary Fig. 9)”.\nSupplementary Fig. 9. Effects of mixing source of shallow-water DIC on vertical carbon isotope\ngradients. Mixing proportions (0~90%;1000 random samples) of shallow-water DIC to carbonate\nsediments on the deep shelf. The vertical δ13C gradients (mean + s.e.m.) are calculated from the δ13C\nDIC\ndifference between surface water (~40 m) and mid-water (~130 m), expressed as the δ13C change per\nDIC\n100 m of water depth (‰/hm). Blue and red points represent estimates shown in Fig. 2 without\nconsidering a mixing source. Sensitivity tests reveal that the trend of increased gradients at the peak CIE\nstage, as indicated by the δ13C data, is largely independent of the mixing proportion during the\ncarb\nPermian-Triassic (a), Early Toarcian (b) and Palaeocene-Eocene Thermal Maximum (c).\nSupplementary Table 3. Evaluation of potential vertical gradient changes during peak CIE under varying\nmixing source conditions. This comprehensive analysis examines 1,000 × 1,000 possible scenarios\n(Supplementary Fig. 9) comparing vertical gradient increases versus decreases at the peak CIE stage,\nunder varying mixing source conditions.\nAge Number of increased Number of decreased Ratio of increased Ratio of decreased\ngradient scenarios gradient scenarios gradient scenarios gradient scenarios\nPTTM 911,307 88,693 91% 9%\nETTM 722,435 277,565 72% 28%\nPETM 888,295 111,705 89% 11%\n3) The authors assume that any change in the vertical d13C gradient is related mainly to a d13C\nchange at depth. This assumption builds on the zeroing of near surface d13C_carb. I see no\narguments that rule out that the change is actually a change in surface d13C_carb while the deep\nd13C_carb is semi fixed? This would implicate that it it could be a d13C mixed layer&atmosphere\nchange and not a result of change in BCP?\nWe appreciate the reviewer's insightful comment regarding the interpretation of vertical δ13C\ngradient changes. It is necessary to analyze the original carbonate δ13C data (Supplementary Fig. 5\nand Supplementary Table 1). The δ13C at different depth all become negative during the event (i.e.\ncarb\nnegative carbon isotope excursion), because of the massive 13C-depleted carbon injected into\natmosphere from large igneous province (Gutjahr et al., 2017; Wu et al., 2023). While atmospheric\nδ13C changes would indeed produce equivalent isotopic shifts in surface and deep waters, our data\nshow that the CIE magnitude at deep shelf (4.5‰) is larger than those from shallow shelf (3.5‰).\nThis depth-dependent signal cannot be explained solely by atmospheric changes, but need\nadditionally mechanisms involved increased vertical δ13C gradient related to biological pump\nchange (Hülse et al., 2021).\nThe method zeroing surface δ13C is to eliminate the negative CIE in proxy data during peak CIE\ncarb\nphases resulting from carbon emissions not included in the simulations, and also removed any\ndiscrepancy between the modelled δ13C values and proxy δ13C values caused by fractionation\nDIC carb\nbetween oceanic DIC and carbonate (invariant and independent of temperature) and diagenesis.\nTo strengthen the argumentation, we have made the following substantive revisions in the lines 116-\n123:“Compiled global δ¹³Ccarb records exhibit a pronounced negative carbon isotope excursion\n(CIE) at three water depth levels (Supplementary Fig. 6; Supplementary Table 1), which is\ncommonly attributed to the injection of large amounts of ¹³C-depleted carbon into the atmosphere-\nocean system11, 12. Although δ¹³C variations in the ocean-atmosphere carbon inventory would\ntheoretically induce equivalent isotopic shifts in both surface and deep waters, our data indicate\nthat the CIE magnitude in the deep shelf or basin is greater than that in the shallow shelf (e.g.,\nPTTM: 4.5‰ vs. 3.5‰). This discrepancy mainly attributed to an enhanced vertical δ¹³C gradient39,\n40”.\nOn a side note the title is misleading as well as many portions of the text. The authors talk about\nmicrobial metabolism or microbial respiration when in fact they can only talk about remineralization\nrate of sinking organic matter. We can not know if the reminiralizatrion is a result mainly of single\ncell microbes or multicellular eukaryotes that performed the respiration.\nThe remineralization scheme in the cGENIE model follows an Arrhenius-type temperature\ndependency (Crichton et al., 2021). This simplified scheme applies to microbes (e.g., bacteria,\narchaea, and single-celled eukaryotes) rather than complex multicellular eukaryotes, which has been\nused to investigate ocean microbial remineralization in deep time, such as the P-Tr mass extinction\nand Cenozoic (Hülse et al., 2021; Boscolo-Galazzo et al., 2021). Therefore, it is clear that our focus\nis on microbial metabolism rather than multicellular eukaryotes.\nWe have added details of this temperature-dependent scheme in the Methods (lines 386-388):\n“Specifically, the dynamic model incorporates Arrhenius-type temperature-dependent scheme for\nthe microbial remineralization of organic matter to regulates POC flux in mid-waters”.\nIn conclusion, considering point 1-3, the manuscrip is too immature for publication at this stage.\nHowever I do look forward to see further enhancements to the analysis in the future.\nChristian J. Bjerrum\nThank you for your constructive comments. These suggestions really help improve our manuscript.\nReferences\nBoscolo-Galazzo, F., Crichton, K.A., Ridgwell, A., Mawbey, E.M., Wade, B.S. and Pearson, P.N.,\n2021. Temperature controls carbon cycling and biological evolution in the ocean twilight zone.\nScience 371, 1148-1152.\nBurgess, S.D., Bowring, S. and Shen, S.. High-precision timeline for Earth's most severe extinction.\nProceedings of the National Academy of Sciences 111, 3316–3321 (2014).\nCramer, B.S., Toggweiler, J.R., Wright, J.D., Katz, M.E. & Miller, K.G. Ocean overturning since\nthe Late Cretaceous: Inferences from a new benthic foraminiferal isotope compilation.\nPaleoceanography 24, PA4216 (2009).\nGutjahr, M. et al. Very large release of mostly volcanic carbon during the Palaeocene–Eocene\nThermal Maximum. Nature 548, 573-577 (2017).\nHülse, D. et al. End-Permian marine extinction due to temperature-driven nutrient recycling and\neuxinia. Nat. Geosci. 14, 862-867 (2021).\nNordt, L., Breecker, D. & White, J. Jurassic greenhouse ice-sheet fluctuations sensitive to\natmospheric CO dynamics. Nat. Geosci. 15, 54-59 (2022).\noxidation after the Ediacaran Hankalchough glaciation. Geochim. Cosmochim. Acta 75, 1357-\n1373 (2011).\nWu, Y.Y. et al. Volcanic CO degassing postdates thermogenic carbon emission during the end-\nReviewer #3 (Remarks to the Author):\nWu et al., calculate the vertical gradient of ∂13C observations across prominent hyperthermal events\nto show significant increases that may reflect a change in watercolumn cycling of organic carbon\n(the Biological Pump). They further use the cGENIE Earth system model to help interpret these\nsignals in the context of changes in ocean circulation, nutrient inventories, and changes in particulate\norganic carbon flux attenuation (transfer efficiency) with depth due to temperature dependent\nremineralisation. The model results demonstrate that changing transfer efficiency with warming can\nsufficiently explain the ∂13C gradients. The authors conclude that the biological pump may\ntherefore have had more of a role in regulating the carbon cycle during these events than previously\nthought.\nI think the integration of observations and modelling here is robust and informative. In particular,\nthe model is used to help robustly identify TE from a range of potential drivers as an explanation of\nthe ∂13C gradients. I have a few comments on clarifying this approach below but I think this is as\nrobust as a paleo data-model integration can be (I can’t comment on the collation and treatment of\nobservations directly). This outcome is very interesting, potentially significant, and to my\nknowledge has not been shown before.\nHowever, I think the manuscript goes too far in terms of interpretation and implications beyond the\ncore hyperthermal TE and ∂13C result. The work on sinking rates as a function of community\nstructure is too simple and not robust to support the level of interpretation and discussion in the\nmanuscript. Likewise, the discussion of implications for future changes in the biological pump and\ncarbon sequestration is not well supported by this work. The authors need to be much more careful\nabout inferring changes in carbon sequestration from TE alone; not overinterpreting the model w.r.t.\na more complex reality; and consider the differences in future projection timescales (100 years) to\nequilibrium changes presented here. I think the paleo-perspective on future changes is incredibly\nvaluable in general, but it has to be properly and carefully caveated due to the differences in the\nresolution of timescales, spatial scales, and processes. Unfortunately, I don’t think latter part of the\nmanuscript achieves this.\nThanks for your positivity comments. We greatly appreciate your insights and recognize the\nimportance of carefully framing our interpretations and implications. Below, we address your key\nconcerns and outline the revisions we have made (details see the responds to each comment).\n1)The section of “Increased POC transfer efficiency after the rise of modern phytoplankton” and\nrelative Figure 3d have been deleted.\n2) Some strong statements about implications for changes in carbon sequestration and future BCP\nchange have been weakened or removed.\n3) We have performed another three multi-factor experiments varying both [PO ] and temperature\n4\ndependence to evaluate the combined effects by these competitive factors.\n4) We have revised the texts to point out some limitations in cGENIE model.\nGeneral Comments\n∂13C gradient interpretation: There needs to be some clarification on how vertical ∂13C gradients\nare being interpreted. An increase in production, without a change in TE, can lead to an increased\nvertical ∂13C gradient (if measured at a single deep horizon) because both the surface and deep\nvalues are changing. A change in TE will have a similar effect but also changes the distribution over\na range of depth levels. Effectively you show this in Figure 2 g-i as the difference between blue bars!\nIf only one depth level is used then I’m concerned you cannot discern the difference. I was unsure\nif the 0 permil reference standard (lines 105 – 106) are an attempt to control for the production\nscenario? If so, this needs to be much more clearly linked. This is not a major problem because\nFigure 2 shows the production change to empirically to be a smaller effect than with TE. However,\nthis is not clear in the discussion of interpreting vertical gradients. Additionally, there are citations\nthat explore vertical ∂13C gradient interpretations that could be cited here to help such as Morée et\nal., (2018, 2021).\nIn the “eutrophication” scenario, doubling nutrient levels enhance primary production, which can\nresult in a greater vertical δ13C gradient, as shown by the blue bars in Figure 2g–i. We have added\nexpanded discussion of this interpretation. Lines 194-197: “The “eutrophication” scenario\ndemonstrates that doubled nutrient levels stimulate productivity (Supplementary Fig. 13),\nconcurrently increasing the vertical δ13C gradient. However, this effect remains substantially\nDIC\nweaker than the gradient increase induced by the “metabolism” scenario (Fig. 2)”.\nDifferent depth levels had been applied to calculate the gradient to make the conclusion robust,\nincluding δ13C difference between surface water (~40 m) and mid-water (~130 m) showed in Fig.\n2, and δ13C difference between ~40 m and ~230 m (Supplementary Fig. 5). We have added the\ncitation of Supplementary Fig. 5.\nThe method zeroing surface δ13C was consistently applied to both palaeoproxy data and model\noutputs to enable direct comparison. We have clarified it (Lines 110-111): \"Surface water δ13C\nDIC\nvalues were normalized to zero for both proxy records and model scenarios to facilitate consistent\nintercomparisons.\"\nThank you for recommending these relative literatures, which significantly strengthened our\nmechanistic interpretation. Relevant citations have been incorporated throughout the introduction\nand discussion .\n∂13C errors: there are notable uncertainties with the ∂13C values from cGENIE (Figure 2), mainly\nat depth. The authors have subsampled the latitudinal regions most represented by the observations,\nbut does sampling model ∂13C from the nearest location reduce the uncertainty? I this this\nessentially asks if the sampling across shelf vs. open ocean environments in cGENIE adds to the\npresented uncertainty?\nThank you for your suggestions. In the previous manuscript, we had sampled modeled δ¹³C from\nthe nearest locations corresponding to each study region (Supplementary Figs. 2–5). Based on these\nestimates, we recalculated the average gradients and obtained steeper gradients with larger\nuncertainties (Figure R1). Distinguishing between shelf and open-ocean environments did not\nappear to effectively reduce uncertainty, despite slight differences in values. This is likely due to the\nuneven spatial distribution of sampling locations. Furthermore, the relatively low spatial resolution\nof the cGENIE model limits its ability to accurately resolve exact locations, often leading to\nconfusion between shelf and open-ocean environments. Therefore, we have decided to retain our\nprevious approach of sampling model results based on latitudinal regions, as this helps avoid issues\narising from spatial unevenness.\nIn response to the Reviewer, we have revised the Methods section to clarify this approach. Lines\n438-439:“As sampling model results based on latitudinal regions helps avoid issues arising from\nspatial unevenness”.\nFigure R1. Comparison of proxy and simulated vertical DIC δ13C gradients in the upper\nocean during the three hyperthermals. The average of simulated δ13C (mean ± 1 s.d.) is\nDIC\ncalculated based on model results from nearest locations.\nMulti-factor experiments: much of the later text in the manuscript discusses net effects of changing\nexport and TE (see lines 238 – 241). This manuscript only presents single factor experiments which\ndon’t directly quantify this. Multi-factor experiments, e.g., varying both [PO4] and temperature\ndependence, might be illuminating here. I don’t see a reason why it should be assumed that these\ninteractions are linearly additive and can be assumed from the current outputs.\nWe have performed another three multi-factor experiments varying both [PO ] and temperature\n4\ndependence to evaluate the combined effects by these competitive factors. Model results reveal that\nthese competing drivers collectively reduce mid-water POC flux during the PTTM, while enhancing\nfluxes during the the ETTM and PETM, despite maintaining consistently reduced POC transfer\nefficiency across all events (Supplementary Fig. 18).\nWe have expanded the discussion through new multi-factor experiments. Lines 253-260: “We\nfurther perform “multi-factor” experiments varying [PO ] and temperature dependency\n4\nsimultaneously to investigate the interplay between nutrient-dependent productivity and\ntemperature-dependent microbial metabolism on the BCP. Model results reveal that these competing\ndrivers collectively reduce mid-water POC flux during the PTTM, while enhancing fluxes during\nthe the ETTM and PETM, despite maintaining consistently reduced POC transfer efficiency across\nall events (Supplementary Fig. 18). This inverse mid-water POC flux pattern may account for the\ncontrasting OC burial responses observed among different hyperthermals”.\nThe limitations and future directions are clarified. Lines 269-273: “However, our simulated mid-\nwater POC fluxes should not be conflated with OC burial, critical preservation processes of OC are\nnot involved in our experiments. A systematic integration of proxy data (e.g. TOC) with models\nexplicitly incorporating OC burial processes will be essential to elucidate how these warming-\ninduced competitive factors control OC burial across different hyperthermals”.\nSupplementary Fig. 18. Simulated changes in biological carbon pump for metabolism and multi-\nfactor scenarios (varying both [PO ] and temperature dependence) during the three hyperthermal\n4\nevents. (a) Global mean particulate organic carbon (POC) exported at surface water (40 m). (b)\nGlobal mean POC flux at mid-water (1040 m). (c) POC transfer efficiency (i.e., the ratio of POC at\nthe 1040 m to those at the surface water).\nSinking rates: I really struggled with the discussion of the sinking rates to explain the difference in\npre-CIE TEs for PTTM and ETTM in Figure 3. As I understand from the manuscript, a reduced\nsinking rate was applied in the first instance for the PTTM experiment to reflect the absence of\nmodern phytoplankton (lines 152 - 156). The text around lines 250 – 259 seems to then try to\ninterpret the difference in TE despite similar temperatures and concludes it’s a reduced sinking\nrate…which seems circular! Are the authors trying to demonstrate that a reduced sinking rate is\nrequired to help explain the ∂13C gradient? In which case, this is not coming across and would need\nthe equivalent of Figure 2g but from the sensitivity experiments. In general, the “Increased POC\ntransfer efficiency after the rise of modern phytoplankton” section is predicated on the applied\nchange in sinking rate being representative and there is no discussion of this. The temperature-\ndependence findings emerge from the underlying multiple constraints represented by the model.\nThese are therefore much stronger than an applied change in a parameter that directly leads to a\npredictable outcome. See also the specific comment on lines 296 – 299 below.\nAfter thorough evaluation, we acknowledge the limitations in the original analysis of sinking rates\non transfer efficiency. Consequently, we have removed the section titled “Increased POC transfer\nefficiency after the rise of modern phytoplankton” and associated Figure 3d, based on two critical\nconsiderations: 1) single parameter of sinking rate does not consider the various abundances,\ndiameters, and diversities of different plankton community, which could not generate robust\nmodeling results that support the statement in this section. 2) The deleted section constituted\nsupplementary speculation rather than model-supported findings, and its removal does not impact\nour core conclusions regarding temperature-dependent microbial metabolism.\nInstead, we added brief discussion state that sinking rate has slight influence on the gradient based\non sensitivity experiments (Supplementary Fig. 14). Base on this, we revise the statement and\nfigures as follows:\n1)The section of “Increased POC transfer efficiency after the rise of modern phytoplankton” and\nrelative Figure 3d have been deleted.\n2) We added discussions on the effect of the sinking rate on gradients based on sensitivity\nexperiments. Lines 199-202: “Additionally, sensitivity experiments reveal that POC sinking rate in\nthe cGENIE model exhibit limited influence on the gradients (Supplementary Fig. 14), reinforcing\nthe primacy of temperature- dependent microbial metabolism”.\n3) All strong claims about sinking-rate effects in the “Implications for future climate and the BCP”\nsection have been removed or revised (see detailed responds to the specific comments).\nSpecific Comments\nLines 24 – 25: “leading to diminished organic carbon sequestration in the ocean interior”. The\nmanuscript infers this but does not quantify it so you can explicitly say this.\nThere is only quantified transfer efficiency (TE) in our study and no modelled results about how\nreduced TE exacerbated the greenhouse effect actually. We revised text to clarify TE’s potential\nfeedback role rather than asserting direct greenhouse effects, and delete all strong statements.\nLines 24-25: “potentially leading to diminished organic carbon sequestration in the ocean interior”. Also\nsee the consistent revision about comments of Lines 222 – 227.\nLines 42 – 46: Sequestration occurs across the whole water column §and is somewhat evenly\ncontributed from remineralisation above and below 1000m (Siegel et al., 2023). Siegel et al., (2023)\ndemonstrates that you could derive a metric of remineralization weighted mean first-passage time\nto reflect this in one quantity…but…the first-passage time is not easily quantifiable! Instead transfer\nefficiency is an approximate metric of where remineralisation is occurring across the water-column.\nThis is somewhat pedantic but I think an important point in the context of later discussions of TE\nand carbon storage.\nWe thank the reviewer for highlighting this critical distinction. We have revised the text to clarify\nthat transfer efficiency (TE) serves as an approximate proxy for the efficiency of biological carbon\nsequestration. Relative citations are added (e.g. Siegel et al., 2023).\nLines 41-48: “This downward transport process stores carbon within the ocean interior where\ncarbon can be sequestered on centennial and longer timescales17, 18. Transfer efficiency, defined as\nthe fraction of POC transported from the surface ocean (0−100 m) to the ocean interior, serves as\na widely used proxy for estimating the biological carbon pump sequestration. Despite transfer\nefficiency represents an operational approximation rather than a direct sequestration measurement,\nit provides critical insights into natural carbon sequestration process in the ocean.”.\nLines 60 – 63: the citation Kwon et al. (2009) does not comment on vertical gradients of ∂13C\nThanks, correction made. We have cited another relevant paper (Kroopnick et al., 1985) here, which\ndescribes the distributions of δ13C and their vertical gradients.\nKroopnick, P.M., 1985. The distribution of 13C of ΣCO in the world oceans. Deep Sea Research\nSupplementary Figures 2 to 4: It’s not clear which symbols are data and which are model.\nThe symbols have been added to Supplementary Figures 2 to 4, consistent with those in Figure 2.\nLines 137 – 138: see general comment above on productivity impact on gradients – is this what is\nmeant by ocean nutrient conditions?\nIn the “eutrophication” scenario, doubling nutrient levels enhances production, leading to a greater\nvertical δ13C gradient. Relevant revisions can be found in the responses to the general comments.\nLine 142: The “stagnation” terms implies a slow down in overturning– is this the actual case in your\nexperiments? The supplementary figure 10 suggests slightly higher export which to me suggests\nmaybe actually greater overturning. Perhaps the name is not quite appropriate and something like\n“continental and climate” might be more representative.\nWhat you think is correct. The slightly greater overturning causes a slightly larger export in the\n“ocean stagnation” scenario. To avoid this misunderstanding, we have replaced “ocean stagnation”\nwith “circulation change” throughout the manuscript.\nLines 147 – 148: you need to be more explicit about what a “static biological pump” is. This could\nbe construed as fixed export production (not responding to a change in nutrient supply) which I\ndon’t think is the case. I think here this means fixed attenuation of POC with depth only.\nWhat you think is correct. The “static biological pump” refers to a fixed attenuation of POC with\ndepth only. We have clarified it as follows: Lines 167-168: “static biological pump characterized\nby a fixed attenuation of the POC depth profile”.\nLines 193 – 195: “aided by warming-induced increases in photosynthetic rates” – is this in cGENIE\nor a wider statement? You need to be clear whether any temperaturedependence applies to export\nproduction in the model (Crichton et al., 2021 for example introduced this, a temperature term in\nthe production equations exists in earlier productivity schemes). In general, the increased POC\nexport with lower TE is well-known (see Kwon et al., 2009 for example). It might be interesting to\nseparate out the impact of vertical redistribution and export on ∂13C gradients though! (I think this\nis beyond this manuscript, so just an interested comment).\nThank you for raising these important points. In our study, the temperature-dependent\nphotosynthetic rate was incorporated into the “metabolism” scenario using a dynamic model,\nwhereas the static model used in the “circulation change” and “eutrophication” scenarios does not\naccount for temperature dependency. Thus, there is only a negligible change in productivity in the\n“circulation change” scenario (Supplementary Fig. 13). These two approaches align with the STND\nand TDEP schemes proposed by Crichton et al. (2021).\nWe have clarified this in the revised manuscript. Lines 217-221: “This effect is further enhanced by\nwarming-induced increases in photosynthetic rates, leading to higher productivity (Fig. 4). The\n“metabolism” scenario incorporates the temperature dependency of the nutrient uptake rate in the\nmodel45. In contrast, the static model used in the “circulation change” and “eutrophication”\nscenarios does not account for temperature dependency (Supplementary Fig. 13)”.\nYour comment on “separating out the impact of vertical redistribution and export on δ13C gradients”\nis truly inspiring. We will explore this in future work.\nLines 207 – 208: This is fairly unsurprising because you have a fixed sinking rate and a\nremineralisation rate that is varying with local temperature in cGENIE. Since temperature is\nmonotonically decreasing with depth and there are no other processes, such as zooplankton\nmigration, you essentially have a power-law curve (Cael & Bisson 2018). In reality this may not be\nthe case if there are processes that are not monotonic with depth I think.\nThank you for raising these important points. We have added some discussions to clarify the\nlimitations of our model about the non-monotonic processes (e.g., zooplankton migration): “without\nconsidering zooplankton migration” (line 242). More relevant revisions can be found in the next\nresponse.\nLines 210 – 216: There are some details here that feel overstated or not fully stated. 1) The\nindependence of global mean [PO4] and TE is totally unsurprising! At least in cGENIE, there is no\nmechanism linking the distribution of remineralisation rates across depth with production. This is\nnot necessarily always the case, e.g., a sinking rate dependent on ecosystem size structure which\nitself is some partial function of total nutrient supply. 2) The latitudinal distribution of TE in cGENIE\nis solely dependent on temperature remineralisation as it’s the only process resolved for POC that\nleads to spatial variability. TE variability simply falls out of the latitudinal distribution of\ntemperature. 3) DeVries & Weber (2017) is one estimate of TE but it doesn’t agree with all\nobservations, such as deep sediment traps (Henson et al., 2012)!\nWe thank for your insightful comments. We have revised the texts to weaken claims and point out\nsome limitations in cGENIE model on Lines 234-242: “Along with our sensitivity analyses of ocean\ntemperature, the independence of transfer efficiency from [PO ] (Supplementary Fig. 16) and\n4\nspatial pattern of higher transfer efficiency in the high latitudes in our simulations (Supplementary\nFig. 17) further support the idea that transfer efficiency is inversely correlated with ocean\ntemperature. Our spatial pattern is consistent with the results of data-assimilated biological pump\nmodels and modern field observations31, though quantitative discrepancies persist with direct\nsediment trap measurements20. This discrepancy potentially arises from the sole temperature-\ndependent remineralization process incorporated in our model, without considering zooplankton\nmigration46 or size-dependent plankton communities47”.\nLines 222 – 227: The claim that reduced TE exacerbated the greenhouse effect and was an important\ncontributor to persistent warming is dependent on how big this effect actually was but is not\nquantified here! Is there any way to quantify this relative to the other carbon cycle changes? If not,\nthese statements need to be toned down towards TE having a potential role.\nThere is only quantified TE in our study and no modelled results about how much reduced TE\nexacerbated the greenhouse effect actually. We revised this part to clarify the TE only having a\npotential role, and delete some strong statements on lines 248-252: “However, a reduction in POC\ntransfer efficiency due to a warming-induced increase in microbial metabolism—as indicated by\nour model results (Figs. 3 and 4)—potentially lead to less CO sequestration in the ocean interior.\nthree studied hyperthermal events, which persisted for hundreds of thousands to millions of years\n(Fig. 1).”\nLines 250 – 259: SST is not really the best metric because remineralisation is driven by in-situ\ntemperature. The mean temperature across the ocean volume might be more representative.\nWe agree that the mean temperature across the ocean volume might be more representative.\nAnyway, we have removed this section. For more details, please refer to our response to the general\ncomments on sinking rates.\nLines 296 – 299: This conclusion cannot be drawn from this study. You applied a change in the\nsinking rate based on the assumption of plankton community changes. That can reproduced by\nsimply taking an empirical flux curve representing the flux of POC with depth (e.g., Cael & Bisson\n2018) and changing the sinking rate with the same assumption. To support this statement, you’d\nneed to show that the plankton community is driving TE independently of a parameter choice you\nmade and which emerges from your model simulations. In other words, you got out what you put\nin, which doesn’t support this conclusion.\nWe agree what you suggest. We have deleted these strong statements about plankton community\nand focused on impact of the temperature dependent remineralization. For more details, please refer\nto our response to the general comments on sinking rates.\nLines 300 – 302: this is far too strong a statement to be supported by this manuscript. The manuscript\nonly considers temperature dependent remineralisation, whereas there are a wide range of other\nprocesses such as oxygen dependence; ballasting (protection via minerals); organic matter lability.\nSee also above comment on plankton community structure and sinking velocities.\nWe have deleted the overstated statements about phytoplankton community and emphasized a wide\nrange of other processes might influence the transfer efficiency on lines 292-295: “Furthermore,\nthis study does not account for other factors, including trait-based functional types of the plankton\ncommunity45, oxygen-dependent processes51, and mineral ballasting effects52, all of which may\naffect transfer efficiency”.\nLines 302 – 308: Again, these statements are too strong and not consistent with the findings of the\nmanuscript. A ∆Carbon storage (by the biological pump) is proportional to ∆TE only in isolation\nand is also dependent on timescale. Export production and the residence time of regenerated carbon\nare changing simultaneously in projected future changes leading to a net increase in carbon storage\nby 2100 despite the general decrease in TE (Wilson et al., 2022). The increase occurs because of\nthe centennial timescale where the decreases in export and TE haven’t had time to propagate through\nthe ocean circulation. This change in carbon storage will differ over longer timescales (>2100) as\nthose competing physical and biological signals are convolved, but, the magnitude and direction is\nan open question which I don’t think this manuscript is able to provide a clear perspective on.\nThank you for highlighting the importance of timescales in carbon sequestration. We agree that our\noriginal statements were overstated. We have revised the text to clarify the importance of timescales\nin carbon storage and to emphasize some model limitations. These revisions better align our\nconclusions with the scope of the manuscript.\nLines 290-298: “However, this positive feedback on carbon sequestration depends on timescales.\nModel projections suggest a net increase in oceanic carbon storage by 2100 due to short timescales\n(<100 years), despite a decrease in transfer efficiency21. Furthermore, this study does not account\nfor other factors, including trait-based functional types of the plankton community45, oxygen-\ndependent processes51, and mineral ballasting effects52, all of which may affect transfer efficiency.\nNevertheless, our integrated approach combining proxy data compilation with mechanistic\nmodeling highlights that temperature- dependent microbial metabolism may play a significant role\nin modulating the BCP dynamics during both ancient hyperthermal events and future climate\nscenarios”."}, {"round": 2, "reviewer_comments": {"reviewer_1": "The authors have addressed all the issues that I mentioned in the first round of review. I have no further comments on the ms and recommend acceptance, and hope to see the paper in print soon.", "reviewer_3": "Thank you to the authors for addressing my previous comments in detail and engaging in a constructive discussion. I am happy that the authors have addressed the general comments raised on the multi-factor experiments; the rise of modern plankton; and the strength of some statements. Overall, the manuscript is more focussed and stronger in its findings because of this.\nThe final future climate section now feels more appropriate as a discussion-level section. The authors have added some clarifications on the complexities of considering future changes based on my previous comments. However, I am still a bit hesitant on this part of the manuscript because by focusing heavily on implications for the future BCP, the authors are predominantly relying on findings in previous papers and are not informing the discussion with the results presented in this manuscript. This is not totally detrimental to the manuscript but is a lost opportunity to stimulate further discussion in this area!\nIn my view, the strength and novelty of the results presented by the authors is that they strongly confirm that a temperature-dependent biological pump is a notable feedback in the carbon cycle and Earth System. Previous work by John et al., (2013) and John et al., (2014) showed its impact over the long-term geological background but it’s role as a feedback is an assumption based on this – this manuscript confirms we need to include it to help explain features of transient warming events. The discussion-level questions that follow in my mind are: what’s the size or significance of this feedback; what timescale is it relevant over; is it sufficiently resolved in Earth System models? The future is a part of this discussion but not the sole focus.\nLine 45 – 48: I think the preceding statement about TE being a proxy for sequestration is sufficient here and responds well to points raised in the previous review. Therefore, these following lines don’t seem as necessary here.\nLines 142 – 143: “assume thousand deepening or shallowing scenarios” – just needs clarification as this is ambiguous as written. It could be interpreted as 1000m or 1000 different scenarios.\nLine 166: reference to “ocean stagnation” scenario left over from previous manuscript version.\nLines 167 – 168: I think it’s worth being explicit by noting that whilst the remineralisation profile is fixed export production is still dynamic in response to surface ocean nutrient concentrations.\nLines 168 – 169: You should state the metabolic scenario includes both temperature-dependence on production and remineralisation. This is noted later on in the manuscript but it’s of particular relevance here in my opinion.\nLines 170 – 175: Apologies if I missed this previously – is the reduced sinking rate applied in all scenarios or just the metabolic scenario? This may be open to interpretation, but I would see this as a long-term background feature similar to a warmer climate/higher CO2 on top of which the hyperthermal occurs. You seem to have demonstrated this is tested for in lines 200 – 202 but again it’s worth being explicit about the experiment design.\nLine 200: is “reveal that the reduced POC sinking rate” correct, if so I think this is clearer.", "reviewer_4": "The authors used available carbon isotopes values from the literature and numerical simulations to reconstructed the response of the biological pump during in 3 past intervals of warming. This response is presented to be affected by a reduced vertical organic carbon transfer efficiency as a response to warming, owing to a positive temperature-dependent microbial remineralization in the ocean interior.\nWith regards to the conceptualization:\nThe reduction of vertical fluxes with warming would not in itself result in a biological pump that is less capable of storing carbon in the ocean interior, if we focus only in the ocean-atmosphere system. That is because this alone would not affect the ratio of preformed over remineralized nutrients in the ocean interior. What would happen with a weaker vertical flux is that more carbon and nutrients are remineralized in shallower depth, but ultimately these nutrients now residing in the upper water column would be eventually consumed in an ocean that is mostly oligotrophic (if they are remineralized in the uppermost layers), restoring the flux to its original value.\nThis can be conceptualized with a 2-box model, where if the upper ocean box is kept oligotrophic at steady state (a fair assumption for most of the ocean). Sarmiento & Toggweiler showed in the 1980s that adding a third box, with partial nutrient consumption in high latitude regions (where nitrate/phosphate > 0) does have a real leverage on atmospheric.\nIf productivity/export changes in high-latitude regions (that ventilate the ocean interior, and thus change its preformed nutrient content and thus carbon storage), the mechanism proposed here could very well work. But the idea that warming alone per se causes the biological pump to represent a positive feedback on the initial CO2 release is not very compelling to me.\nIf we are talking about the implication of metabolic theory on rates of organic matter burial in the seafloor, then what the authors discuss could be relevant for atmospheric CO2. But even then, such mechanism could mostly relevant for the Permo-Trias and the Early Toarcian (Million years timescales), I am not sure if it would be relevant on the timescales of the PETM, which onset is estimated to have lasted from 2000-20000 years. Carbon burial today, in a much colder ocean baseline, is already pretty small (0.15 Pg C per year, Kandasamy & Nath, 2016, Table 3), in comparison to the 3000 PgC of carbon expected to have been released during the PETM (Zeebe et al., 2009).\nExport production proxies are unfortunately ridden with caveats and cannot help us much. For example, there is no commonly accepted evidence that export fluxes declined during the PETM. A lot of superficial assessments have been made using Ba fluxes being higher during peak warmth, which are now been shown to be artifacts of variable sedimentation rates (Torfstein et al., 2010). Barium isotopes (flux independent) show the diametric opposite trend to Ba fluxes (Bridgestock et al. 2019), which does not support\nWith regards to the observations:\nThe paper argue for bigger carbon isotope gradients during warming intervals. I did look at published results from the PETM in Sites 527, 401 and 1209 for planktic and benthic and compiled them (see picture attached to this review), The figure seems to show that:\n1) the vertical 13C gradient is smaller, not larger, during the peak-PETM.\n2) Most of the change is driven by planktic (upper ocean) drop in d13C, much less in the deep-dwelling benthic organisms (where enhanced remineralization should drive most of the signal).\nWhile this represent a more superficial exploration of the available records in comparison to the more meticulous work done by the authors, I do find these discrepancies very puzzling.\nI do appreciate the work done by authors in combining modelling and observations in the context of a very important topic such as the biological pump, evaluated across different warming intervals of Earth’s history. I personally would suggest that the aforementioned concerns should be addressed before this manuscript is published."}, "author_response": "Point-by-point response to the reviewers’ comments.\nThose changes appear as highlights in red colour in the manuscript. Below are our responses (blue\ncolour) to reviewers’ comments. Line numbers in our response refer to the revised manuscript.\nREVIEWER COMMENTS\nReviewer #1 (Remarks to the Author):\nThe authors have addressed all the issues that I mentioned in the first round of review. I have no\nfurther comments on the ms and recommend acceptance, and hope to see the paper in print soon.\nWe appreciate your positive comments for this manuscript.\nReviewer #3 (Remarks to the Author):\nThank you to the authors for addressing my previous comments in detail and engaging in a\nconstructive discussion. I am happy that the authors have addressed the general comments raised on\nthe multi-factor experiments; the rise of modern plankton; and the strength of some statements.\nOverall, the manuscript is more focused and stronger in its findings because of this.\nWe sincerely thank the reviewer for taking time to evaluate our manuscript and positive feedback\non our previous revision.\nThe final future climate section now feels more appropriate as a discussion-level section. The\nauthors have added some clarifications on the complexities of considering future changes based on\nmy previous comments. However, I am still a bit hesitant on this part of the manuscript because by\nfocusing heavily on implications for the future BCP, the authors are predominantly relying on\nfindings in previous papers and are not informing the discussion with the results presented in this\nmanuscript. This is not totally detrimental to the manuscript but is a lost opportunity to stimulate\nfurther discussion in this area!\nThank you for your positive comments on our revisions. We agree that the final section lacks some\ndiscussions in relation to the results presented in this manuscript. In response, we have made detailed\nrevisions as outlined in the following reply, based on your suggestions.\nIn my view, the strength and novelty of the results presented by the authors is that they strongly\nconfirm that a temperature-dependent biological pump is a notable feedback in the carbon cycle and\nEarth System. Previous work by John et al., (2013) and John et al., (2014) showed its impact over\nthe long-term geological background but it’s role as a feedback is an assumption based on this – this\nmanuscript confirms we need to include it to help explain features of transient warming events. The\ndiscussion-level questions that follow in my mind are: what’s the size or significance of this\nfeedback; what timescale is it relevant over; is it sufficiently resolved in Earth System models? The\nfuture is a part of this discussion but not the sole focus.\nThank you for your constructive suggestions on the final section. They were very helpful in\nimproving the manuscript. In response, we have added further discussion to emphasize the\nimportance of temperature-dependent microbial metabolism. Specifically, we now summarize the\nmagnitude of change in transfer efficiency resulting from a 4–10 °C temperature increase, and we\nalso clarify the relevant timescales associated with this process. In addition, a comparison of transfer\nefficiency across different Earth system models is available in the previous version (lines 301-305).\nThe title of this section has also been revised to “Implications for BCP change and future climate”\nto better reflect the content and its relevance to the results.\nLines 285-298: “Our findings highlight that temperature-dependent microbial metabolism may\nsignificantly influence BCP change and thus climate regulation. Specifically, a temperature increase\nof 4–10℃ during three ancient hyperthermals could reduce transfer efficiency of BCP by 2–6\ntimes, thereby lowering the flux of POC into the ocean interior and potentially counteracting OC\nburial (Fig. 3). This feedback likely resulted in reduced CO sequestration in the ocean interior,\n2\namplifying ancient warming (Fig. 4). Conversely, a cooling-induced decline in microbial\nmetabolism enhanced the efficiency of BCP over the past 15 Myr34. It should be noted, however, this\ncritical feedback becomes more effective over longer timescales, as organic carbon burial operates\non a million-year timescale 58. Consequently, the amplifying effect during the PETM may have been\nrelatively limited due to its shorter duration compared to the PTTM and ETTM (Fig. 1), which\nrestricted the establishment of sustained positive feedbacks. Furthermore, a key limitation of our\nsimulations is the uncertainty in how the 2–6-fold reduction in transfer efficiency decreases the\namount of CO sequestered in the ocean by the BCP for a geological timescale, which should be\n2\naddressed in future work.”.\nLine 45 – 48: I think the preceding statement about TE being a proxy for sequestration is sufficient\nhere and responds well to points raised in the previous review. Therefore, these following lines don’t\nseem as necessary here.\nMany thanks. We have deleted these lines to avoid redundancy.\nLines 142 – 143: “assume thousand deepening or shallowing scenarios” – just needs clarification as\nthis is ambiguous as written. It could be interpreted as 1000m or 1000 different scenarios.\nThank you for pointing out this ambiguity. We agree that the original phrasing was unclear. To clarify,\nthe sentence refers to 1000 different scenarios of deepening or shallowing for each hyperthermal\nevent. We have revised the sentence in the manuscript accordingly. Lines 140-142: “…we assume\n1000 different scenarios of deepening or shallowing for each hyperthermal event to assess the\nvertical gradient responses to depth changes during the event (see Methods)”.\nLine 166: reference to “ocean stagnation” scenario left over from previous manuscript version.\nCorrection made. The“ocean stagnation” is corrected as “circulation change”\nLines 167 – 168: I think it’s worth being explicit by noting that whilst the remineralisation profile\nis fixed export production is still dynamic in response to surface ocean nutrient concentrations.\nYes, export production of static model remains dynamic and responds to surface ocean nutrient\nconcentrations. We have revised relevant sentence to explicitly state that. Lines 167-168:\n“…whereas export production could still vary in response to surface nutrient concentrations”.\nLines 168 – 169: You should state the metabolic scenario includes both temperature-dependence on\nproduction and remineralisation. This is noted later on in the manuscript but it’s of particular\nrelevance here in my opinion.\nWe agree this. We have added relevant description on lines 169-171: “…higher temperatures raise\nthe microbial metabolic rate, thereby altering photosynthesis and POM remineralization depth\nprofile”.\nLines 170 – 175: Apologies if I missed this previously – is the reduced sinking rate applied in all\nscenarios or just the metabolic scenario? This may be open to interpretation, but I would see this as\na long-term background feature similar to a warmer climate/higher CO2 on top of which the\nhyperthermal occurs. You seem to have demonstrated this is tested for in lines 200 – 202 but again\nit’s worth being explicit about the experiment design.\nWe apologize for our mistake. The reduced sinking rate is only applied in metabolism scenario. We\nhave revised relative text in discussion and methods to clarify this.\nLines 170-171: “In the “microbial metabolism” scenario, we reduced the POM sinking rate by ~22%\nin the PTTM experiments…”\nLines 201-202: “…sensitivity experiments under the “microbial metabolism” scenario reveal…”\nLines 444-445: “…as for “microbial metabolism” scenario, we set the default sinking rate…”.\nLine 200: is “reveal that the reduced POC sinking rate” correct, if so I think this is clearer.\nCorrection made.\nReviewer #4 (Remarks to the Author):\nThe authors used available carbon isotopes values from the literature and numerical simulations to\nreconstructed the response of the biological pump during in 3 past intervals of warming. This\nresponse is presented to be affected by a reduced vertical organic carbon transfer efficiency as a\nresponse to warming, owing to a positive temperature-dependent microbial remineralization in the\nocean interior.\nWe sincerely thank the reviewer for their thorough and constructive evaluation of our manuscript.\nThe comments provided valuable perspectives that helped us to clarify key concepts. The\nsuggestions concerning the conceptual framing and the treatment of proxy data were especially\nhelpful. We have revised the manuscript accordingly, and a detailed, point-by-point response to each\nof your comments is provided below.\nWith regards to the conceptualization:\nThe reduction of vertical fluxes with warming would not in itself result in a biological pump that is\nless capable of storing carbon in the ocean interior, if we focus only in the ocean-atmosphere system.\nThat is because this alone would not affect the ratio of preformed over remineralized nutrients in\nthe ocean interior. What would happen with a weaker vertical flux is that more carbon and nutrients\nare remineralized in shallower depth, but ultimately these nutrients now residing in the upper water\ncolumn would be eventually consumed in an ocean that is mostly oligotrophic (if they are\nremineralized in the uppermost layers), restoring the flux to its original value.\nWe agree that increased remineralization of carbon and nutrients at shallower depths facilitates their\nreturn to the surface. Since phosphate is the only limiting nutrient in our simulation, the recycled\nphosphate can be efficiently utilized to enhance primary productivity. This increase in productivity\nindeed results in enhanced surface export production fluxes, as simulated by cGENIE (Fig. 3a).\nSome relative discussions are available in previous version (Lines 215-223).\nHowever, our results indicate that this compensatory flux is insufficient to fully offset the reduction\nin mid-water POC flux caused by warming-enhanced microbial remineralization (Fig. 3b).\nSpecifically, warming substantially increases microbial respiration rates throughout the upper and\nmid-water columns. The temperature-dependent increase in remineralization rates is non-linear and\nappears to dominate mid-water POC fluxes in our simulations. Accordingly, our findings highlight\nthat temperature-dependent microbial metabolism may significantly influence BCP change and thus\nclimate regulation. We have revised this part to make it clearly.\nLines 236-239: “Despite the increased POC export production in the surface water, our results\nindicate that this compensatory flux is insufficient to fully offset the reduction in mid-water POC\nflux due to warming-induced increases in microbial respiration, leading to a significant decrease in\nPOC transfer efficiency (Fig. 3b, c)”.\nThis can be conceptualized with a 2-box model, where if the upper ocean box is kept oligotrophic\nat steady state (a fair assumption for most of the ocean). Sarmiento & Toggweiler showed in the\n1980s that adding a third box, with partial nutrient consumption in high latitude regions (where\nnitrate/phosphate > 0) does have a real leverage on atmospheric.\nIf productivity/export changes in high-latitude regions (that ventilate the ocean interior, and thus\nchange its preformed nutrient content and thus carbon storage), the mechanism proposed here could\nvery well work. But the idea that warming alone per se causes the biological pump to represent a\npositive feedback on the initial CO2 release is not very compelling to me.\nWe appreciate the reviewer’s insightful comment on the classical view that high productivity\nenhances the oceanic carbon sink. For example, Sarmiento and Toggweiler (1984) showed that a 70\nppmv increase in CO at the end of the last ice age was primarily driven by a decrease in net high-\n2\nlatitude productivity, possibly coupled with change in thermohaline overturning. However, some\nresearchers have argued that quantifying productivity or carbon export alone is insufficient to assess\nthe biological pump’s impact on atmospheric CO , as its efficiency also depends on the residence\n2\ntime of exported carbon in the ocean (e.g. DeVries et al., 2012; Siegel et al., 2023). Carbon in the\nupper ocean typically returns to the atmosphere on sub-centennial timescales, limiting its\ncontribution to long-term CO storage. In contrast, carbon transported to the deep ocean or buried\n2\nin sediments can remain sequestered for significantly longer periods (Siegel et al., 2023). As noted\nin your following comment, temperature-dependent microbial metabolism can influence the rate of\norganic matter burial. Accordingly, in the context of ancient hyperthermals, the transfer efficiency\nof POC to the deep ocean, rather than primary production, is a more appropriate metric for assessing\nthe long-term sequestration capacity of the biological carbon pump.\nTo clarify this issue, we have added discussions highlighting the importance of high productivity in\nenhancing the oceanic carbon sink, as shown in previous studies (e.g., Sarmiento and Toggweiler,\n1984). Furthermore, we clarify that in our study, transfer efficiency, rather than export production,\nis a more appropriate metric for estimating the long-term sequestration capacity of the biological\ncarbon pump on geological timescales.\nLines 226-232: “Productivity or export production is an important metric for evaluating the\nsequestration capacity of the biological carbon pump15, 48. However, quantifying export production\nalone is insufficient to assess the biological pump’s impact on atmospheric CO , as its efficiency\n2\nalso depends on the residence time of exported carbon in the ocean17-19. Specifically, carbon in the\nupper ocean (<1,000 m) typically returns to the atmosphere on sub-centennial timescales, limiting\nits contribution to long-term CO storage, whereas carbon transported to the deep ocean (>1,000\n2\nm) or buried in sediments can remain sequestered for much longer periods.”\nWe acknowledge that transfer efficiency is an indirect proxy, as our simulations cannot directly\ncalculate how a 2–6-fold reduction in transfer efficiency decreases the amount of CO sequestered\n2\nin the ocean by the BCP for a geological timescale. Because CO₂ is used as a forcing to simulate\ntemperature rise in our model, it cannot be modified by feedback processes within the simulation.\nWe have clearly stated this limitation of our study on lines 295-298: “Furthermore, a key limitation\nof our simulations is the uncertainty in how the 2–6-fold reduction in transfer efficiency decreases\nthe amount of CO sequestered in the ocean by the BCP for a geological timescale, which should\n2\nbe addressed in future work.”.\nIf we are talking about the implication of metabolic theory on rates of organic matter burial in the\nseafloor, then what the authors discuss could be relevant for atmospheric CO2. But even then, such\nmechanism could mostly relevant for the Permo-Trias and the Early Toarcian (Million years\ntimescales), I am not sure if it would be relevant on the timescales of the PETM, which onset is\nestimated to have lasted from 2000-20000 years. Carbon burial today, in a much colder ocean\nbaseline, is already pretty small (0.15 Pg C per year, Kandasamy & Nath, 2016, Table 3), in\ncomparison to the 3000 PgC of carbon expected to have been released during the PETM (Zeebe et\nal., 2009).\nYes, we agree that temperature-dependent remineralization can potentially influence the rate of\norganic matter burial, a critical control that has often been overlooked in previous studies which\nprimarily focused on factors such as redox conditions. We also recognize the importance of the\ntimescale over which this feedback develops and manifests. Accordingly, we have clarified in the\nrevised manuscript that temperature-dependent remineralization likely had a limited impact during\nthe PETM due to the relatively short duration of this event.\nLines 285-295: “Our findings highlight that temperature-dependent microbial metabolism may\nsignificantly influence BCP change and thus climate regulation. Specifically, a temperature increase\nof 4–10℃ during three ancient hyperthermals could reduce transfer efficiency of BCP by 2–6\ntimes, thereby lowering the flux of POC into the ocean interior and potentially counteracting OC\nburial (Fig. 3). This feedback likely resulted in reduced CO sequestration in the ocean interior,\n2\namplifying ancient warming (Fig. 4). Conversely, a cooling-induced decline in microbial\nmetabolism enhanced the efficiency of BCP over the past 15 Myr34. It should be noted, however, this\ncritical feedback becomes more effective over longer timescales, as organic carbon burial operates\non a million-year timescale 58. Consequently, the amplifying effect during the PETM may have been\nrelatively limited due to its shorter duration compared to the PTTM and ETTM (Fig. 1), which\nrestricted the establishment of sustained positive feedbacks”.\nExport production proxies are unfortunately ridden with caveats and cannot help us much. For\nexample, there is no commonly accepted evidence that export fluxes declined during the PETM.\nA lot of superficial assessments have been made using Ba fluxes being higher during peak warmth,\nwhich are now been shown to be artifacts of variable sedimentation rates (Torfstein et al., 2010).\nBarium isotopes (flux independent) show the diametric opposite trend to Ba fluxes (Bridgestock et\nal. 2019), which does not support\nThank you for share the conflicting interpretations of export production during the PETM based on\nbarium fluxes and barium isotopes. This discrepancy reflects a broader emerging debate regarding\nchanges in primary productivity during the PETM. Our cGENIE simulations indicate an increase in\nexport production, which is consistent with the trend suggested by barium isotope studies. We have\nnow cited Bridgestock et al. (2019) in the revised manuscript (line 213) to support this interpretation.\nWith regards to the observations:\nThe paper argue for bigger carbon isotope gradients during warming intervals. I did look at\npublished results from the PETM in Sites 527, 401 and 1209 for planktic and benthic and compiled\nthem (see picture attached to this review), The figure seems to show that:\n1) the vertical 13C gradient is smaller, not larger, during the peak-PETM.\n2) Most of the change is driven by planktic (upper ocean) drop in d13C, much less in the deep-\ndwelling benthic organisms (where enhanced remineralization should drive most of the signal).\nWhile this represent a more superficial exploration of the available records in comparison to the\nmore meticulous work done by the authors, I do find these discrepancies very puzzling.\nThank you for pointing out this issue. We also appreciate the well-compiled data and figures that\nyou provided. We have noted that the difference in δ13C between planktonic and deep-ocean benthic\nforaminifera becomes smaller during the PETM. In other words, the magnitude of the carbon isotope\nexcursion (CIE) in planktonic foraminifera is larger than that in deep-ocean benthic foraminifera, a\npattern also observed in previous studies (Zachos et al., 2007; McInerney and Wing, 2011). One\nexplanation is that the smaller CIE magnitude returned by pelagic sites (such as Sites 527, 401 and\n1209) might result from incomplete preservation stemming from carbonate dissolution, rendering\nthem unsuitable for reconstructing the vertical gradient. Furthermore, the strongest vertical δ13C\ngradient usually occurs in the upper ocean (>1500 m), whereas the gradient is less pronounced\nbelow this depth (Figure R1 from Dentith et al., 2020). As a result, deep-ocean benthic\nforaminifera are unsuitable for reconstructing the vertical gradient.\nFigure R1. Depth profiles of δ13C during the 1990s (Dentith et al., 2020)\nDIC\nAccordingly, in addition to the carbonate records, we used paired δ13C and δ18O data from\nplanktonic foraminifera to reconstruct the vertical δ13C gradient in the upper ocean during the PETM\nin previous version (Supplementary Figure 7 and 8), excluding benthic foraminiferal records. This\napproach follows previous study (Boscolo-Galazzo et al., 2021). Planktonic foraminifera records\nindicate an increasing trend in vertical gradients, from 0.65 ± 0.16‰/hm to 1.21 ± 0.41‰/hm, which\nfurther supports our findings. We have incorporated relevant discussion in the Methods section to\nclarify the limitations associated with benthic foraminiferal records. We hope it helps clarify this\nissue.\nLines 394-398: “Deep-ocean benthic foraminifera are excluded from this approach, as previous\nstudies have noted that the vertical δ¹³C gradient between planktonic and deep benthic foraminifera\ndecreasesd during the PETM64, 65. One possible explanation is that the smaller CIE magnitude\nobserved at pelagic sites may result from incomplete preservation due to carbonate dissolution65,\nrendering them unsuitable for reconstructing the vertical gradient”.\nI do appreciate the work done by authors in combining modelling and observations in the context\nof a very important topic such as the biological pump, evaluated across different warming intervals\nof Earth’s history. I personally would suggest that the aforementioned concerns should be addressed\nbefore this manuscript is published.\nWe really appreciate your positive comment on our works. We have tried to address all issues and\nimprove our manuscript. Thanks again.\nReferences\nBoscolo-Galazzo, F., Crichton, K.A., Ridgwell, A., Mawbey, E.M., Wade, B.S. and Pearson, P.N.,\n2021. Temperature controls carbon cycling and biological evolution in the ocean twilight zone.\nScience 371, 1148-1152.\nDentith, J.E., Ivanovic, R.F., Gregoire, L.J., Tindall, J.C. and Robinson, L.F., 2020. Simulating\nstable carbon isotopes in the ocean component of the FAMOUS general circulation model with\nMOSES1 (XOAVI). Geoscientific Model Development 13, 3529-3552.\nDeVries, T., Primeau, F. and Deutsch, C., 2012. The sequestration efficiency of the biological pump.\nGeophysical Research Letters 39.\nMcInerney, F.A. and Wing, S.L., 2011. The Paleocene-Eocene Thermal Maximum: A Perturbation\nof Carbon Cycle, Climate, and Biosphere with Implications for the Future. Annual Review of\nEarth and Planetary Sciences 39, 489-516.\nZachos, J.C., Bohaty, S.M., John, C.M., McCarren, H., Kelly, D.C. and Nielsen, T., 2007. The\nPalaeocene–Eocene carbon isotope excursion: constraints from individual shell planktonic\nforaminifer records. Philosophical Transactions of the Royal Society A: Mathematical,\nPhysical and Engineering Sciences 365, 1829-1842."}, {"round": 3, "reviewer_comments": {"reviewer_3": "Thank you to the authors to their constructive responses! The authors have addressed all of my comments really well and I don't have any further comments.", "reviewer_5": "Review Wu et al \"Microbial metabolism amplified warming in three Phanerozoic 2 hyperthermal events\"\nI was invited to assess response to reviewer 4 comments, as this reviewer was no longer available to review. Below I summarize first some of my findings that may help improve the manuscript (there are some inconsistencies), followed by reviewer 4 response checks.\nThis manuscript compares cGENIE model d13C of DIC output with carbon isotope gradients between surface and 'mid-depth' through a data synthesis. While the model results are very interesting, I am not entirely convinced that data supports them.\n-The data show a very high spread; especially during the Permian-Triassic, and the spread of data is bigger than we see in the top 250 meters of the ocean today (-1 to +2) at all three depths for the Permian-Triassic scenaria, and below 100 meters for the early Toarcian scenario. Please discuss why this could be.\n-Are the selected sites representative of continental shelf deposits? Can these reflect the whole ocean (which cGENIE outputs represent)? Please discuss.\nOther comments:\nAbstract\nLine 25: do you mean to say in the deeper ocean interior?\nI am a bit confused about what you are referring to in terms of organic carbon sequestration, and wondering whether you are mixing this up with sequestration of CO2? Organic carbon sequestration refers to the burial of organic material in sediments.\nLines 30-47 are all about the BCP and transfer efficiency of POC, and refers to sequestering of carbon in the ocean interior on centennial and longer timescales.\nThe switch to organic carbon sequestration is abrupt and needs further explanation.\nLines 58-60: this is a generalisation, it is a combination of remineralisation of organic material, ocean circulation and ageing water mass (see a more recent paper I reviewed). See lines 155-159.\nLines 109-111: Planktonic foraminiferal shells, which are more resistant to diagenesis than bulk carbonates, preserve pristine oceanic DIC signals and their δ18O values serve as a paleodepth proxy34. Boscolo-Galazzo et al. (2021) do not suggest that planktic foraminifera d13C pristinely preserve ocean DIC signals, and distinctly mention the influence of the physiology of biomineralization. For this Boscolo-Gallazo et al use constrained shell size fractions to minimize this physiological disequilibrium effect on d13C of planktonic foraminifera. For the time-periods in the current study, they authors likely were not able to use constrained shell size fractions. Please can you include a short section of what the implications of this may be, and how it adds to uncertainty.\nLines 124-125; CIE and ETTM do not appear to be statistically different.\nLines 139-140: Despite – proxy data. Is there a word missing? Sentence currently unclear.\nLine 166 what does a static biological pump mean? Does it mean that productivity cannot vary?\nThe sections comparing data and model output doesn’t seem to take in consideration the range of vast range of d13C values at the different depths. It would be good to see this reflected in the offsets e.g. mean +/- standard deviation or something.\nWhile most of the model results seem to be more or less similar to data during pre-CIE (or at least they fall within the range of d13C values reported for the different depths; although skewed to less depleted values during the Toarcian), it is notable that in the case for Peak CIA at mid-water depth during the Early Toarcian the model seems to be dissimilar. Please explain why this may be\nLine 226-227: Productivity and export production are not the same…\nLines 291-293: rephrase “as organic carbon burial operates on a million-year timescale” to “as the effects of organic carbon burial on carbon cycling typically become apparent on a million-year timescale”\nLines 305-308: where are the cGENIE future projections?\nFigure 3: what is this like for modern; otherwise difficult to understand the implications.\nSupplementary info:\nSupplementary figure 2 what does it show?\nResponse to reviewer 4 comments:\nI believe the authors have addressed these as well as they could."}, "author_response": "Point-by-point response to the reviewers’ comments.\nThose changes appear as highlights in red colour in the manuscript. Below are our responses (blue\ncolour) to reviewers’ comments. Line numbers in our response refer to the revised manuscript.\nREVIEWER COMMENTS\nReviewer #3 (Remarks to the Author):\nThank you to the authors to their constructive responses! The authors have addressed all of my\ncomments really well and I don't have any further comments.\nThank you for reviewing our manuscript.\nReviewer #5 (Remarks to the Author):\nReview Wu et al \"Microbial metabolism amplified warming in three Phanerozoic hyperthermal\nevents\"\nI was invited to assess response to reviewer 4 comments, as this reviewer was no longer available\nto review. Below I summarize first some of my findings that may help improve the manuscript (there\nare some inconsistencies), followed by reviewer 4 response checks.\nWe sincerely thank you for providing constructive comments, which have been helpful in improving\nour manuscript and assisting in evaluating our responses to Reviewer 4’s comments\nThis manuscript compares cGENIE model d13C of DIC output with carbon isotope gradients\nbetween surface and 'mid-depth' through a data synthesis. While the model results are very\ninteresting, I am not entirely convinced that data supports them.\n-The data show a very high spread; especially during the Permian-Triassic, and the spread of data is\nbigger than we see in the top 250 meters of the ocean today (-1 to +2) at all three depths for the\nPermian-Triassic scenario, and below 100 meters for the early Toarcian scenario. Please discuss\nwhy this could be.\nYes, your observation is correct. The δ13C values at the same depth during our study interval\nDIC\nshow large variations of approximately 3–4‰, whereas those in the modern ocean exhibit smaller\nvariations of about 2–3‰ (Supplementary Fig. 10; data from GLODAPv2, Olsen et al., 2016). These\nvariations in both the modern and paleo-ocean primarily reflect differences in ocean circulation,\nnutrient availability, and temperature-dependent metabolic processes among different ocean basins\nor latitude.\nCompared with the modern data, the slightly larger δ13C variability observed in the geological\nDIC\nrecords likely results from additional local effects, including mixed DIC sources, variable water\ndepths, and diagenetic alteration when using the carbonate δ13C as a proxy for δ13C . In the\nDIC\nprevious version, Reviewer 2 raised a similar concern; we therefore performed sensitivity analyses\nto evaluate the potential influence of DIC mixing and fluctuating paleodepths on the estimated δ13C\ngradients; Supplementary Figs. 10 and 11). The results demonstrate that the trend of increased\nvertical δ13C gradients during the peak CIE stage persists despite the uncertainties in δ13C (see\nDIC DIC\ndetailed discussions in Lines 143-159)\n.\nIn this revised version, we further added analyses of modern δ13C data (Supplementary Fig. 2)\nDIC\nfor comparison and included additional discussion.\nLines 138-143: Considerable δ¹³C variability exists even within the same depositional settings (e.g.,\n3.8‰ in shallow shelf during the P–Tr interval), exceeding the 2–3‰ variability observed in\nmodern δ13C (Supplementary Fig. 2). This greater variability may reflect mixing of DIC sources\nDIC\nand fluctuating paleodepths when carbonate δ13C is used as a proxy. Both factors contribute to\nuncertainty in estimates of vertical δ13C gradients, potentially making changes in gradients\nstatistically insignificant.\nSupplementary Fig. 2. Modern DIC carbon isotope depth profiles. Depth profiles of δ13C for the global\nDIC\nocean (a), continental shelf (b), open ocean (c), low latitudes (d), mid-latitudes (e), and high latitudes (f).\nThe continental shelf is defined as areas with a maximum water depth <500 m. The global ocean is\ndivided by latitude into low (0–25°), mid (25–50°), and high (50–75°) latitude zones, as the maximum\nlatitude of available data is ~75°. The vertical δ¹³C gradient between the surface and mid-water (130 m)\nis 0.46‰/hm for the global ocean, 0.49‰/hm for the continental shelf, and 0.46‰/hm for the open ocean.\nThe gradient between the surface and mid-water (∼130 m) differs among low (0.61‰/hm), mid-\n(0.30‰/hm), and high-latitude (0.50‰/hm) regions. Data are from the GLODAPv2 dataset156.\n-Are the selected sites representative of continental shelf deposits? Can these reflect the whole ocean\n(which cGENIE outputs represent)? Please discuss.\nAs for the δ13C data (Fig. 2), nearly all selected sites are from continental shelf settings, as no\ncarb\npelagic sediments are available prior to the Middle Jurassic. Based on δ13C data from the modern\nDIC\nocean (Supplementary Fig. 2a-c), we infer that δ13C values from continental shelves can largely\nDIC\nrepresent those of the global ocean, provided that latitude is consistent.\nFirst, modern δ13C data show that the vertical δ13C gradient between the surface and mid-water\nDIC\n(∼130 m) is nearly identical for the global ocean (0.46‰/hm), continental shelf (0.49‰/hm), and\nopen ocean (0.46‰/hm). This suggests that δ13C data from continental shelf environments in\nDIC\ngeological records can robustly represent the global ocean.\nSecond, modern δ13C data indicate that the vertical δ13C gradient varies with latitude. The\nDIC\ngradient between the surface and mid-water (∼130 m) differs among low (0.61‰/hm), mid-\n(0.30‰/hm), and high-latitude (0.50‰/hm) regions. This highlights the importance of maintaining\nconsistent latitudinal coverage when comparing proxy data and model results. We had implemented\nthis approach in our analysis in the original manuscript, with the model subsampling the latitudinal\nrange (48.6° S–48.6° N) most represented by the observations (see Methods).\nThese clarifications have been incorporated into the main text, Methods, and Supplementary\nInformation.\nLines 96-98: The compiled δ13C data from the modern ocean demonstrate that δ13C values\nDIC DIC\nfrom continental shelf environments provide a robust proxy of the global ocean’s carbon isotopic\ncomposition (Supplementary Fig. 2).\nLines 345-348: Although our compiled δ13C data are derived primarily from continental shelf\ncarb\nsettings and lack pelagic records, modern δ13C data demonstrate that δ13C values from\nDIC\ncontinental shelf environments provide a robust representation of the global ocean (Supplementary\nFig. 2).\nLines 491-492: As vertical δ13C gradients vary across different latitudes according to observations\nfrom the modern ocean (Supplementary Fig. 2d-f).\nOther comments:\nAbstract\nLine 25: do you mean to say in the deeper ocean interior?\nI am a bit confused about what you are referring to in terms of organic carbon sequestration, and\nwondering whether you are mixing this up with sequestration of CO2? Organic carbon sequestration\nrefers to the burial of organic material in sediments. Lines 30-47 are all about the BCP and transfer\nefficiency of POC, and refers to sequestering of carbon in the ocean interior on centennial and longer\ntimescales.\nThe switch to organic carbon sequestration is abrupt and needs further explanation.\nYes, it means carbon sequestration in the deeper ocean interior. In order to avoid confusing, we have\nrevised this part. line 25: “potentially leading to diminished carbon sequestration in the deeper\nocean interior” and also lines 54-55: “…a clear understanding of the feedbacks between oceanic\ncarbon sequestration and climate change across the Phanerozoic”.\nLines 58-60: this is a generalisation, it is a combination of remineralisation of organic material,\nocean circulation and ageing water mass (see a more recent paper I reviewed). See lines 155-159.\nWe appreciate the reviewer’s insightful comment. We agree that the transfer efficiency of POC is\nnot solely governed by temperature-driven microbial metabolism, but rather by the combined effects\nof organic matter remineralization, ocean circulation, and the water mass ventilation. We have\nrevised the sentence accordingly to better reflect this complexity.\nLines 58-62: Such temperature-driven changes in microbial remineralization strongly affect the\ndegradation of sinking organic matter30, 31. Together with ocean circulation and water mass\nventilation, which regulate circulation timescales and pathways of sinking organic carbon, these\nprocesses collectively determine the transfer efficiency of POC21, 32, 33.\nLines 109-111: Planktonic foraminiferal shells, which are more resistant to diagenesis than bulk\ncarbonates, preserve pristine oceanic DIC signals and their δ18O values serve as a paleodepth\nproxy34. Boscolo-Galazzo et al. (2021) do not suggest that planktic foraminifera d13C pristinely\npreserve ocean DIC signals, and distinctly mention the influence of the physiology of\nbiomineralization. For this Boscolo-Gallazo et al use constrained shell size fractions to minimize\nthis physiological disequilibrium effect on d13C of planktonic foraminifera. For the time-periods in\nthe current study, they authors likely were not able to use constrained shell size fractions. Please can\nyou include a short section of what the implications of this may be, and how it adds to uncertainty.\nThank you for highlighting the effect of shell size fractions on planktonic foraminiferal δ13C, which\nwas considered in our previous manuscript. Boscolo-Galazzo et al. (2021) used constrained shell\nsize fractions (250–355 μm) to minimize the influence of physiological disequilibrium on δ13C.\nIn our dataset, 24 of the 32 paired planktonic foraminiferal δ13C and δ18O profiles used a similar\nconstrained size fraction (212–355 μm), while the remaining 8 samples either lack size information\nor used a smaller size fraction (180–212 μm) (see details in Supplementary Data 1). Smaller shells\ntypically show lower δ13C and δ18O values due to vital effects (Birch et al., 2013), which could lead\nto lower estimated δ13C and shallower inferred water depths. Consequently, smaller shell sizes\nDIC\nmay slightly exaggerate vertical δ13C gradients, although this bias is likely limited given the\nDIC\nconstrained size fractions applied to most samples.\nWe have clarified that constrained size fractions were used to minimize uncertainties associated\nwith the vital effect, and we have added corresponding descriptions in the Methods section.\nLines 114-116: Planktonic foraminiferal shells can record oceanic δ13C when using constrained\nDIC\nshell size fractions to minimize physiological effects (see Methods), because smaller-size shells\ntypically show lower δ13C and δ18O values due to vital effects34.\nLines 416-422: In our dataset, the shells of 24 profiles used a constrained size fraction (212–355 μm)\nfollowing previous work34, while the remaining 8 samples either lack size information or used a\nsmaller size fraction (180–212 μm). Smaller shells typically show lower δ13C and δ18O values due\nto vital effects34, which could lead to lower estimated δ13C and shallower inferred water depths.\nDIC\nConsequently, the use of smaller shell sizes may slightly exaggerate vertical δ13C gradients,\nDIC\nalthough this bias is likely limited given the constrained size fractions applied to most samples.\nLines 124-125; CIE and ETTM do not appear to be statistically different.\nWe thank the reviewer for this point. Indeed, given the relatively large uncertainties resulting\nfrom mixing of DIC sources and fluctuating paleodepths, the changes in vertical δ13C gradients\nDIC\nare not statistically significant. However, the observed increased trends are consistent with our\nmodel results (Fig. 2) and previous studies (Meyer et al., 2011; Song et al., 2013). We have\nperformed relative sensitivity analyses of mixing of DIC sources, as well as fluctuating paleodepths\nto discuss the uncertainty in estimates of vertical carbon isotope gradients (see lines 143-159).\nWe have revised the text to clarify this limit caused by large uncertainties.\nLines 141-143: Both factors contribute to uncertainty in estimates of vertical δ13C gradients,\npotentially making changes in gradients statistically insignificant.\nLines 139-140: Despite – proxy data. Is there a word missing? Sentence currently unclear.\nMany thanks. This section addresses how we evaluate the potential uncertainties in carbon isotope\ngradients resulting from fluctuating palaeodepths (i.e., shallowing or deepening) at different\nlocations. We have rephrased this part to make the discussion clearer.\nLines 150-154: Furthermore, in order to evaluate the potential uncertainty resulting from\npalaeodepth fluctuations (shallowing or deepening) at locations, we simulated 2000 different\npalaeodepth scenarios for each hyperthermal event to assess the vertical gradient responses to\ndepth changes during the event (see Methods), as these changes could not be reliably quantified\nthrough existing proxy data.\nLine 166 what does a static biological pump mean? Does it mean that productivity cannot vary?\nProductivity in both static and dynamic biological pumps can vary depending on nutrient availability.\nHowever, the static biological pump assumes a fixed attenuation of the POM depth profile, whereas\nthe dynamic biological pump is characterized by temperature-dependent POM remineralization. We\nhave revised the relevant descriptions to clarify these differences further.\nLines 176-182: “Both the “circulation change” and “eutrophication” scenarios were run using a\nstatic biological pump, in which global POM export production and remineralization are almost\nentirely controlled by nutrient availability (here, [PO ])26. In contrast, the “microbial metabolism”\n4\nscenario used a constant 1× [PO ] and a dynamic biological pump, controlled by both nutrient\n4\navailability and temperature-dependent metabolic rates. In particularly, the dynamic model\nincorporates a temperature-dependent scheme for microbial remineralization of organic matter to\nregulate POC flux in mid-waters”.\nThe sections comparing data and model output doesn’t seem to take in consideration the range of\nvast range of d13C values at the different depths. It would be good to see this reflected in the offsets\ne.g. mean +/- standard deviation or something.\nThank you for pointing this out. We have performed error propagation calculations when comparing\nproxy data and model results. Relative standard deviations have been included in the descriptions,\nand all results are provided in Supplementary Table 4.\nLine 191-192: In this scenario, the differences between proxy and model data at peak CIE phase\nare −0.3 ± 2.7‰/hm for the PTTM, 0.7 ± 2.4‰/hm for the ETTM, and −0.04 ± 2.9‰/hm for the\nPETM, while those for the other two scenarios are notably larger (1.0–2.3 ‰/hm).\nWhile most of the model results seem to be more or less similar to data during pre-CIE (or at least\nthey fall within the range of d13C values reported for the different depths; although skewed to less\ndepleted values during the Toarcian), it is notable that in the case for Peak CIA at mid-water depth\nduring the Early Toarcian the model seems to be dissimilar. Please explain why this may be\nWe appreciate this comment. There is indeed a large discrepancy between the proxy data and the\nmodel results during the peak CIE phase of the ETTM. Most proxy records are derived from the\nwestern Tethys European shelf, a region characterized by restricted ocean circulation. The high\norganic matter burial observed there indicates elevated productivity within this constrained shelf\nsystem (Kemp et al., 2022), which may explain the more depleted mid-water carbon isotope values\nobserved in our proxy data (Fig. 2). In contrast, the relatively low spatial resolution of the model’s\npaleogeographic configuration does not adequately capture the restricted conditions of the European\nshelf, and therefore fails to reproduce the corresponding carbon isotope signatures.\nIn response, we have now added a discussion of the explanation.\nLines 199-204: In particular, there are large discrepancies between the proxy data (most from\nEuropean shelf) and model results during the peak CIE phase of the ETTM (Fig. 2b, e). This\nmismatch likely arises from the relatively low spatial resolution of the cGENIE model, which does\nnot adequately capture the restricted conditions of the European shelf and therefore fails to\nreproduce the 13C-depleted DIC at mid-water depths observed in the proxy data.\nLine 226-227: Productivity and export production are not the same…\nThank you for your comment. We agree that productivity and export production are distinct concepts.\nIn our revised manuscript, we have focused solely on export production, which is the relevant\nquantity for our modeling result. References to productivity have been removed in this sentence to\navoid confusion.\nLine 244: Export production is an important metric for evaluating the sequestration capacity of the\nbiological carbon pump.\nLines 291-293: rephrase “as organic carbon burial operates on a million-year timescale” to “as\nthe effects of organic carbon burial on carbon cycling typically become apparent on a million-year\ntimescale”\nCorrection made.\nLines 305-308: where are the cGENIE future projections?\nThe cGENIE future projections are taken from a previous study (Crichton et al., 2023) rather than\nfrom our work. To avoid ambiguity, we have rephrased this sentence.\nLines 324-327: “Based on this, previous studies using the cGENIE model project that globally\naveraged transfer efficiencies at a depth of 1,040 m will decline from approximately 11% to 4–9%\nunder various carbon emission scenarios22”.\nFigure 3: what is this like for modern; otherwise difficult to understand the implications.\nThank you for the suggestions. We have added the modeled pre-industrial export production, POC\nflux at mid-water, and POC transfer efficiency based on the cGENIE model. These results are taken\nfrom a previous cGENIE simulation (Crichton et al., 2023).\nFig. 3 Simulated changes in biological carbon pump during the three hyperthermal events. (a)\nGlobal mean particulate organic carbon (POC) exported at surface water (40 m) from PTTM (red),\nETTM (blue) and PETM (green). (b) Global mean POC flux at mid-water (1,040 m). (c) POC transfer\nefficiency (i.e., the ratio of POC at the mid-water to those at the surface water). The pre-industrial levels\n(gray) based on the cGENIE model are from previous work22.\nSupplementary info:\nSupplementary figure 2 what does it show?\nSupplementary Figure 2 is similar to Figure 2, as both show the comparison of proxy and modeled\nvertical δ13C gradients. The difference is that Figure 2 presents globally compiled data, whereas\nDIC\nSupplementary Figure 2 groups the data by different ocean regions. For example, during the P–Tr\ninterval, eight ocean regions are considered (e.g., East Paleotethys, Central Paleotethys). Similarly,\nSupplementary Figures 3 and 4 show regional comparisons for the EETM and PTEM, respectively.\nWe have revised the three figure captions to make it clear in the supplementary information.\nResponse to reviewer 4 comments:\nI believe the authors have addressed these as well as they could.\nWe really appreciate reviewer 5 to help assess our responses to reviewer 4 comments.\nReferences\nBirch, H., Coxall, H.K., Pearson, P.N., Kroon, D. and O'Regan, M., 2013. Planktonic foraminifera\nstable isotopes and water column structure: Disentangling ecological signals. Marine\nMicropaleontology 101, 127-145.\nBoscolo-Galazzo, F., Crichton, K.A., Ridgwell, A., Mawbey, E.M., Wade, B.S. and Pearson, P.N.,\n2021. Temperature controls carbon cycling and biological evolution in the ocean twilight zone.\nScience 371, 1148-1152.\nCrichton, K.A., Wilson, J.D., Ridgwell, A. and Pearson, P.N., 2021. Calibration of temperature-\ndependent ocean microbial processes in the cGENIE.muffin (v0.9.13) Earth system model.\nGeoscientific Model Development 14, 125-149.\nKemp, D.B., Suan, G., Fantasia, A., Jin, S. and Chen, W., 2022. Global organic carbon burial during\nthe Toarcian oceanic anoxic event: Patterns and controls. Earth-Science Reviews, 104086.\nMeyer, K.M., Yu, M., Jost, A.B., Kelley, B.M. and Payne, J.L., 2011. δ13C evidence that high\nprimary productivity delayed recovery from end-Permian mass extinction. Earth and Planetary\nScience Letters 302, 378-384.\nOlsen, A., Key, R.M., van Heuven, S., Lauvset, S.K., Velo, A., Lin, X., Schirnick, C., Kozyr, A.,\nTanhua, T., Hoppema, M., Jutterstr Om, S., Steinfeldt, R., Jeansson, E., Ishii, M., P Erez, F.F.\nand Suzuki, T., 2016. The Global Ocean Data Analysis Project version 2 (GLODAPv2) - an\ninternally consistent data product for the world ocean. Earth System Science Data 8, 297-323.\nSong, H.Y., Tong, J.N., Algeo, T.J., Horacek, M., Qiu, H., Song, H.J., Tian, L. and Chen, Z., 2013.\nLarge vertical δ13C DIC gradients in Early Triassic seas of the South China craton: implications\nfor oceanographic changes related to Siberian Traps volcanism. Global and Planetary Change\n105, 7-20."}]

Decision letter