
==== Front
Sci Adv
Sci Adv
sciadv
advances
Science Advances
2375-2548
American Association for the Advancement of Science

adp5156
10.1126/sciadv.adp5156
Research Article
Earth, Environmental, Ecological, and Space Sciences
SciAdv r-articles
Climatology
Oceanography
Climatology
Enhanced ocean heat storage efficiency during the last deglaciation
Enhanced deglacial ocean heat storage efficiency
https://orcid.org/0000-0002-9330-4294
Zhu Chenyu Conceptualization Data curation Formal analysis Funding acquisition Investigation Methodology Resources Software Validation Visualization Writing - original draft 1 †
https://orcid.org/0000-0002-7954-6932
Sanchez Saray Conceptualization Data curation Formal analysis Funding acquisition Investigation Methodology Project administration Resources Software Validation Visualization Writing - review & editing 2 †
https://orcid.org/0000-0003-4554-2666
Liu Zhengyu Conceptualization Data curation Funding acquisition Methodology Project administration Supervision Writing - original draft 3 4 *
https://orcid.org/0000-0003-4752-6045
Clark Peter U. Conceptualization Investigation Methodology Supervision Validation Writing - review & editing 2 *
https://orcid.org/0000-0003-1453-9691
He Chengfei Investigation Methodology Writing - review & editing 5
https://orcid.org/0000-0002-7758-8629
Wan Lingfeng Data curation Investigation Resources Software Validation 6
Lu Jiuyou Investigation 7
Zhu Chenguang Visualization 8
https://orcid.org/0000-0002-1818-6080
Li Lingwei Visualization Writing - review & editing 3 9
https://orcid.org/0000-0001-6569-9842
Zhang Shaoqing Writing - review & editing 6
https://orcid.org/0000-0002-9854-0392
Cheng Lijing Conceptualization Funding acquisition Methodology Writing - review & editing 1 10
1 Earth System Numerical Simulation Science Center, Institute of Atmospheric Physics, Chinese Academy of Sciences, Beijing, China.
2 College of Earth, Ocean, and Atmospheric Sciences, Oregon State University, Corvallis, OR 97331, USA.
3 Department of Geography, Ohio State University, Columbus, OH 43210, USA.
4 School of Geography Science, Nanjing Normal University, Nanjing, China.
5 Rosenstiel School of Marine, Atmospheric, and Earth Science, University of Miami, Miami, FL 33149, USA.
6 Frontier Science Center for Deep Ocean Multispheres and Earth System (DOMES), Institute for Advanced Ocean Study (IAOS) and Key Laboratory of Physical Oceanography, MOE (POL), Ocean University of China, Qingdao, China.
7 Laoshan Laboratory, Qingdao, China.
8 School of Ocean Sciences, China University of Geosciences (Beijing), Beijing, China.
9 Institute of Arctic and Alpine Research, University of Colorado Boulder, Boulder, CO 80303, USA.
10 International Center for Climate and Environment Sciences, Institute of Atmospheric Physics, Chinese Academy of Sciences, Beijing, China.
* Corresponding author. Email: liu.7022@osu.edu (Z.L.); peter.clark@oregonstate.edu (P.U.C.)
† These authors contributed equally to this work.

20 9 2024
20 9 2024
10 38 eadp515628 3 2024
13 8 2024
Copyright © 2024 The Authors, some rights reserved; exclusive licensee American Association for the Advancement of Science. No claim to original U.S. Government Works. Distributed under a Creative Commons Attribution NonCommercial License 4.0 (CC BY-NC).
2024
The Authors
https://creativecommons.org/licenses/by-nc/4.0/ This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is not for commercial advantage and provided the original work is properly cited.

Proxy reconstructions suggest that increasing global mean sea surface temperature (GMSST) during the last deglaciation was accompanied by a comparable or greater increase in global mean ocean temperature (GMOT), corresponding to a large heat storage efficiency (HSE; ∆GMOT/∆GMSST). An increased GMOT is commonly attributed to surface warming at sites of deepwater formation, but winter sea ice covered much of these source areas during the last deglaciation, which would imply an HSE much less than 1. Here, we use climate model simulations and proxy-based reconstructions of ocean temperature changes to show that an increased deglacial HSE is achieved by warming of intermediate-depth waters forced by mid-latitude surface warming in response to greenhouse gas and ice sheet forcing as well as by reduced Atlantic meridional overturning circulation associated with meltwater forcing. These results, which highlight the role of surface warming and oceanic circulation changes, have implications for our understanding of long-term ocean heat storage change.

Deglacial ocean heat storage is enhanced by warming of intermediate-depth waters.

http://dx.doi.org/10.13039/100000001 National Science Foundation AGS2202860 http://dx.doi.org/10.13039/100000001 National Science Foundation OPP2103032 http://dx.doi.org/10.13039/100023581 National Science Foundation Graduate Research Fellowship Program 2234662 http://dx.doi.org/10.13039/501100001809 National Natural Science Foundation of China 42106013 http://dx.doi.org/10.13039/501100001809 National Natural Science Foundation of China 42130604 Science and Technology Innovation Project of Laoshan Laboratory LSKJ202203303 Science and Technology Innovation Project of Laoshan Laboratory LSKJ202203303 Science and Technology Innovation Center of Laoshan Laboratory LSKJ202203303
==== Body
pmcINTRODUCTION

As one of the largest heat reservoirs in the climate system, the global ocean absorbs a large fraction of the excess energy generated by the Earth’s energy imbalance (EEI) during transient climate change (1, 2). Now, more than 90% of the excess energy from ongoing anthropogenic warming has accumulated in the ocean (3–5), strongly buffering global warming as well as causing an increase in ocean temperature, ocean heat content, and sea level. Moreover, the lag in the climate response to a forcing is a function of the rate of ocean heat uptake (2), with changes in that rate affecting the EEI (6, 7). In the past century, the greatest warming from the EEI has occurred in the upper 500 m in the ocean (4, 5, 8, 9), with relatively weak warming in the deep ocean because of its long response time (10) and the lingering impact of previous cooling associated with the Little Ice Age (11–13). This response can be expressed as a ratio between changes of global mean ocean temperature (∆GMOT) and global mean sea surface temperature (∆GMSST) or a measure of the ocean heat storage efficiency (HSE), with current HSE (∆GMOT/∆GMSST) << 1 (~0.1).

Paleoceanographic observations, however, suggest that ∆GMOT can be comparable to or even larger than ∆GMSST on millennial and longer (near-equilibration) timescales. This implies that the ocean can serve as a far greater reservoir of energy in the climate system than implied by contemporary observations. Deep ocean temperature (DOT) proxies from marine sediment cores (14–16) and GMOT proxies from Antarctic ice cores (7, 17, 18) suggest that GMOT increased by ~2.8°C from the Last Glacial Maximum [LGM; ~26 to19 thousand years ago (ka)] to pre-Industrial (PI; Fig. 1C), comparable to (slightly smaller than) the concurrent change in GMSST reconstructed from SST proxies of around ~3.0° ± 0.2°C (Fig. 1B) (19–22). Recent air-sea disequilibrium studies suggest that the LGM ocean cooling reconstructed from noble gases in ice cores may be biased by −0.38° ± 0.37°C which is highly sensitive to model high-latitude winds [Fig. 1C, (23)]. Given the good agreement between two estimates for mean LGM DOT cooling relative to PI from marine proxy records [−2.63°C from (15); −2.82° ± 0.17°C from (16)] and the ice core MOT estimate that neglects air-sea disequilibrium [−2.83° ± 0.5°C from the spline reconstruction of (7)], we suggest that the effects of air-sea disequilibrium may be negligible. In any event, the overall GMOT warming magnitude is comparable to the increase in GMSST, leading to a deglacial ocean HSE of ~0.7 to 0.9, which is much closer to 1 (~1) (Fig. 1D) compared to modern HSE. Furthermore, observations show that the HSE exceeds 1 during the millennial-scale deglacial events Heinrich Stadial 1 (HS1; ~18 to 14.7 ka) and the Younger Dryas (YD; 12.9 to 11.5 ka) (Fig. 1D).

Fig. 1. Model-data comparison for deglacial ocean temperature evolution.

(A) Deglacial forcings in iTRACE: June solar insolation (SOLIN) at 45°N (orange), atmospheric CO2 concentration (green) (70), and meltwater fluxes (MWF) in the Northern (black) and Southern (gray) hemispheres. (B) Deglacial GMSST as deviations from PI in iTRACE (red line; centennial mean) and proxy-based reconstructions (blue) (19, 20) and a reanalysis product (black) (31). Reconstructions of the LGM GMSST cooling from Tierney et al. (21) without (orange) and with data assimilation (DA) (green, with 95% confidence interval) and from Clark et al. (22) (purple) are also shown. Shading/error bars: 1σ uncertainty. (C) Deglacial GMOT as deviations from PI in iTRACE (red line; centennial mean) and our proxy-based global ocean temperature stack (orange line, Materials and Methods) as well as by noble gas–based GMOT reconstructions (blue) (7) in spline (line) and samples (dots) version. Also shown are DOT reconstructions for the last deglaciation (green squares) (15) and for the LGM mean (blue and gray triangles) (16, 23). Note that the LGM estimates in (B) and (C) represent averages over the LGM interval. (D) Ocean HSE in proxy reconstructions (blue) and in iTRACE (red). HSE is calculated as ∆GMOT/∆GMSST with a threshold of 0.1°C for both ∆GMSST and ∆GMOT (relative to the LGM state). For proxy-based HSE, GMSST is from Shakun et al. (19) and Marcott et al. (20), and GMOT is from Shackleton et al. (7). ppm, parts per million.

These observations raise several questions: What mechanisms are responsible for ocean heat uptake/storage and how large can the ocean HSE be? In particular, the mechanisms responsible for the large HSE during the last deglaciation (Fig. 1D) are puzzling. Given the large volume of deep and abyssal waters (below 2000 m), the conventional view in paleoceanography is that warming of their source waters in the Southern Ocean and North Atlantic (17, 24–26) would lead to a large increase in GMOT since the LGM. During the last deglaciation, however, proxy records suggest these source waters, especially those for the Antarctic Bottom water, were largely covered by sea ice in winter (27, 28) when deep water forms, such that their temperatures should have remained near the freezing point (~−1.8°C). We would thus expect a much weaker warming of the deep ocean than the surface ocean during the last deglaciation with an HSE much less than 1.

Here, we combine a transient climate simulation from the LGM to mid-Holocene (MH; 6 ka) with new proxy-based DOT reconstructions to show that the deglacial ocean heat storage is enhanced in intermediate-depth waters in response to increasing greenhouse gases (GHGs) and retreating ice sheets as well as reductions in the Atlantic meridional overturning circulation (AMOC) forced by meltwater. The resulting strong warming at these depths then substantially increases ocean HSE throughout the last deglaciation.

RESULTS

Enhanced deglacial ocean heat storage by warming of intermediate waters

To understand the deglacial ocean HSE, we first examine results from our isotope-enabled TRAnsient Climate Evolution (iTRACE) simulation which was performed using the Community Earth System Model version 1.3 (iCESM1) (29). The simulation was forced by reconstructed ice sheets, solar insolation associated with orbital configurations (ORB), GHGs, and meltwater fluxes (MWF) (Materials and Methods) (30). Starting from the LGM, our simulated GMSST and GMOT both increase by ~2.5°C during HS1, plateau during the Bølling-Allerød (14.7 to 12.9 ka), and then increase again by ~1°C during the YD before reaching the relatively stable Holocene state, corresponding to an ocean HSE of ~1 for the LGM-MH transition (Fig. 1, B to D). These changes in the timing and magnitude of GMSST and GMOT largely follow proxy reconstructions (Fig. 1, B and C) (7, 16–18, 19–22) and a reanalysis product (31). The comparable deglacial (MH-LGM) changes of GMOT and GMSST (i.e., HSE of ~1) is robust across climate models as shown in another transient deglacial simulation TraCE-21k (Materials and Methods) (32), a series of time-slice equilibrium simulations for the last deglaciation (33), and equilibrium simulations in most Paleoclimate Modelling Intercomparison Project Phase 3 (PMIP3) models (34) and in iCESM1, regardless of their differing climate sensitivities (fig. S1 and table S1).

Our model results show that the pattern of deglacial warming in the ocean interior is far from uniform (Fig. 2 and fig. S2A). As expected, we find that the global zonal mean MH-LGM temperature change in the abyss is indeed smaller than GMSST change (Fig. 2A) due to source waters being largely covered by sea ice (Fig. 3A). The interior warming is instead dominated by a core region of intermediate waters (300 to 1500 m) extending northward across the equator from the Southern Ocean, where the warming is greater than the GMSST warming by up to two times. This intermediate-depth warming occurs primarily in the Indo-Pacific Oceans which account for ~70% of the global ocean volume (Fig. 2B and fig. S2B). A secondary core region of warming occurs at intermediate to mid-depths (300 to 3000 m) around the northern subpolar region (Fig. 2A) which is caused predominantly by warming in the Atlantic sector (Fig. 2C). While cores of warming in the Indo-Pacific and Atlantic reach different depths, these cores are related in that they both represent the depth of penetration of intermediate/deep water masses sourced at the surface in these regions. The overall warming of the global ocean is stronger than GMSST warming in the upper half of the ocean but similar to or weaker than GMSST warming in the lower half of the ocean, with the resulting GMOT warming being comparable to that of GMSST (HSE of ~1; Fig. 1D). A similar vertical structure of warming with an even greater amplitude can also be seen at some other times, notably at HS1 (relative to LGM) with a greater warming contribution from the Atlantic sector (Fig. 2, D to F), corresponding to an HSE greater than 1 (Fig. 1D). Last, this ocean warming pattern of MH-LGM, accompanied by the large HSE, is also robust across most PMIP3 models (fig. S3 and table S1; Materials and Methods).

Fig. 2. Simulated and reconstructed deglacial subsurface temperature changes.

(A to C) Normalized zonal mean ocean temperature change between the MH (6.5 to 6 ka mean) and LGM (20 to 19 ka mean) states (T*, divided by MH-LGM GMSST change) in the global ocean (A), Indo-Pacific (B), and Atlantic (C). Contours in (A) show overturning stream function in the MH (interval of 4 Sv). (D to F) As (A) to (C) but for HS1 (16 to 15 ka mean) minus LGM. Symbols in (A) to (C) show location of sites used in our regional stack analysis. (G and H) Regional stacks for the Indo-Pacific (G) and Atlantic (H) (shading: 1σ uncertainty; Materials and Methods). Also shown in (G) and (H) is the global stack for comparison. Note that the classification of water depths used in (G) and (H) are approximate since the sites are arranged not only by depth but also by latitude to assess the simulated pattern.

Fig. 3. Simulated deglacial surface temperature changes.

(A) Changes between the MH-LGM in annual mean SAT (green contour, interval of 5°C) and annual mean SST (color shading, °C). Blue and red lines indicate sea ice edge (SIC; defined as 15% annual sea ice coverage) in the LGM (blue) and MH (red), respectively. (B) Zonal mean change of annual SST. (C) As (B) but for SAT. (D to F) As (A) to (C) but for changes between HS1 and the LGM. Blue and red lines in (D) indicate SIC in the LGM (blue) and HS1 (red), respectively. Also shown in (B) and (C) are temperature changes in a reanalysis product (dashed lines) (31).

Although existing global DOT (15) and GMOT (7) proxy reconstructions support a large deglacial ocean HSE (Fig. 1D), they are unable to constrain the vertical temperature structure simulated in iTRACE (Fig. 2, A to C). We addressed this issue by reconstructing deglacial ∆DOT from 119 globally distributed, well-dated benthic δ18O (δ18Ob) records archived in the Ocean Circulation and Carbon Cycling (OC3) database (Fig. 2, B and C) (Materials and Methods) (35). The site depths provided with the records cover intermediate to abyssal waters. We derived ∆DOT by accounting for changes in the δ18O of seawater (δ18Osw) associated with the decrease in global ice volume. Independently derived ∆DOT reconstructions that are available at several of our sites show similar changes as in our reconstructions (fig. S4). We also find that our basin-weighted global mean ∆DOT (∆GMDOT) stack is in good agreement with the noble gas–based ∆GMOT reconstructions (Fig. 1C), consistent with a recent assessment of the temperature and land ice seawater contributions to the global δ18Ob record (Materials and Methods) (7). This agreement at local and global scales suggests that any effects of hydrographic changes on δ18Osw were negligible (Materials and Methods).

To assess the vertical temperature structure simulated by iTRACE, we constructed regional ∆DOT stacks from sites occurring in the areas where simulated ocean temperatures show lesser, similar, or greater amounts of change than GMSST (Fig. 2, B and C). The zonal plots in Fig. 2 identify the simulated vertical temperature changes between the MH and LGM (Fig. 2, A to C) and HS1 and LGM (Fig. 2, D to F). Overall, the iTRACE stacks follow the reconstructed deglacial temperature evolution, although temperature variations at intermediate depths caused by AMOC reductions are overestimated in amplitude in the model (fig. S5). In the Atlantic basin, the regional stack showing the greatest warming corresponds to the area of greatest simulated warming, in the intermediate to mid-depths (generally <3000 m and >0°). Here, warming above ∆GMDOT begins at the start of HS1 with the subsequent temperature evolution exhibiting pronounced millennial-scale variability (Fig. 2H). The regional stack representing the remaining Atlantic basin north of 30°S is consistent with iTRACE in showing less warming beyond the core region of greatest warming. Our stacks from the Pacific basin similarly support the vertical structure in the iTRACE simulation, with greater warming than ∆GMDOT at intermediate depths (<1500 m), little difference at mid-depths, and less warming than ∆GMDOT in the abyss (Fig. 2G). Further support for this vertical temperature structure is suggested from regional δ18Ob stacks constructed in (36), with sites from 1000 to 2000 m depth being consistently lighter (i.e., warmer) than deeper sites (fig. S6). We thus conclude that the data support the iTRACE finding that the large HSE is caused by strong warming in the intermediate-depth waters.

Enhanced deglacial HSE by GHGs and ice sheet forcing

Here, we show that substantial warming of intermediate-depth waters can be dynamically related to surface warming in response to GHGs and ice sheet forcing. In iTRACE, annual SSTs warm globally from the LGM to the MH except in the polar regions, with the greatest SST warming (~6°C) centered on ~50°S in the southern oceans and ~45°N in the North Pacific and North Atlantic oceans (Fig. 3, A and B). A similar surface pattern is seen in proxy data (37–39), combined data-model (data assimilation) reanalysis (Fig. 3, B and C) (31, 40, 41), and PMIP3 model simulations (fig. S7).

The deglacial SST warming pattern is forced primarily by increasing GHGs and retreating Northern Hemisphere (NH) ice sheets. This is demonstrated by the strong similarity of the SST warming pattern in iTRACE (Fig. 3, A and B) and in the sensitivity experiment forced just by ice sheets, orbital forcing, and GHGs [continental ice sheets (ICE) + ORB + GHG] (Materials and Methods) (fig. S8). Radiative forcing from GHGs causes global surface warming that is characterized by polar amplification toward both poles in surface air temperature (SAT) (Fig. 3, A and C, and fig. S8, D and F) (42). However, at high latitudes in both hemispheres, the ocean surface is largely covered by sea ice, especially in the cold season (fig. S9) when deep water masses are formed by deep convection, such that SSTs remain around the freezing point with little change during the deglaciation (Fig. 3, A and B, and figs. S8 and S9). This leads to peak SST warming in the mid-to-subpolar latitudes, in contrast to the peak SAT warming toward both poles. In the North Atlantic, warming of newly exposed surface waters in the Greenland-Iceland-Norwegian Seas forms another SST warming peak of 6°C at ~60°N (Fig. 3, A and B, and figs. S7 to S9). At the same time, the lowered elevation and reduced albedo of retreating NH ice sheets induce surface warming over NH continents which is then advected downstream to the North Atlantic and North Pacific, enhancing peak SST warming at 45°N (Fig. 3, A and B, and fig. S10, D to F). This response is confirmed by the sensitivity experiment ICE which is forced only by ice sheet change (fig. S10 and Materials and Methods). Peak SST warming at these latitudes is then ventilated into the ocean interior (fig. S8, A to C) mainly by the mean subduction (43, 44), forming the maximum subsurface warming centered in the intermediate waters, primarily in the Antarctic Intermediate Water. The deeper warming in the North Atlantic associated with the mean North Atlantic Deep Water (NADW) formation does not contribute substantially to the GMOT change due to its small volume (fig. S8). Notably, the strong SST warming remains colocated with strong ventilation regions throughout the deglaciation (fig. S11) and highly correlated with heat flux into the ocean in mid-latitude, wind-driven subduction regions, thus enhancing the intermediate ocean heat storage (43). Even with homogeneous surface warming, as the deep-water formation sites are largely ice covered, the intermediate water masses will enhance the HSE because they represent the main remaining pathway of ventilating the global ocean. Last, our sensitivity experiments show that orbital forcing has only a minor effect on global ocean temperature change during the last deglaciation (Fig. 4B).

Fig. 4. Simulated and reconstructed deglacial evolution of AMOC and ocean temperature changes.

(A) AMOC transport simulated in iTRACE (red) and implied in observed 231Pa/230Th from Bermuda Rise sediments (black) (46). (B) GMSST (dotted line) and GMOT (solid line) change (°C, relative to the LGM state) in iTRACE (ICE + ORB + GHG + MWF) and the other three factorized-forcing sensitivity experiments. (C) As (B) but for ocean HSE. HSE is calculated as ∆GMOT/∆GMSST with a threshold of 0.1°C for both ∆GMSST and ∆GMOT (both relative to the LGM state).

Enhanced deglacial HSE by AMOC weakening

Although we find substantial warming of intermediate waters by surface warming associated with GHGs and ice sheet forcing, our sensitivity experiments demonstrate that reaching an ocean HSE of ≥1 throughout the early deglaciation requires additional subsurface warming that is caused by weakening of the AMOC in response to meltwater forcing (Fig. 4). Paleoceanographic proxies and modeling studies suggest that meltwater forcing from NH ice sheets impeded NADW formation during HS1 and the YD, causing a reduction in the AMOC (Fig. 4A; 32, 45, 46). In response to the weakened AMOC, the northward heat transport was reduced, forming the SST bipolar seesaw pattern of NH cooling and SH warming (Fig. 3, D to F, and figs. S12A and S13A) (19, 32, 47–52). This SST pattern does not change GMSSTs substantially because the NH SST response nearly cancels the opposite SH SST response, as shown by the small increase in GMSST during HS1 and the YD that accompanied an AMOC reduction (~0.5°C; comparing iTRACE with ICE + ORB + GHG) (Fig. 4, A and B). Furthermore, the small GMSST change caused by the AMOC reduction quickly disappears after the end of the meltwater forcing, as seen in GMSSTs being nearly the same after 10 ka in the experiments with (iTRACE) and without (ICE + ORB + GHG) meltwater forcing, such that there is little net effect of meltwater forcing on GMSST in the MH (Fig. 4B).

In contrast to this GMSST response, a reduction of the AMOC generates a longer-term subsurface warming in both hemispheres. For example, during HS1, the reduced North Atlantic deep convection reduces the vertical mixing between colder surface and warmer subsurface waters, resulting in a cooling in the North Atlantic surface but a warming in the subsurface, which is a robust feature in reconstructions (50–52) and model simulations (32, 53, 54). The warming signal propagates southward along the Deep Western Boundary Current into the South Atlantic in hundreds of years, before leaking northward into other basins via the Southern Ocean (55). This leads to strong subsurface warming that extends to mid-depths over much of the global ocean that persists throughout HS1 (Figs. 2, D to F, and 4B and figs. S12 and S13) (55). The strong subsurface warming may have also warmed abyssal waters through circulation and mixing processes (figs. S12E and S13E). In iTRACE, the weakening AMOC throughout HS1 contributes 2°C warming to GMOT, while the other combined forcings contribute another 0.5°C warming (Fig. 4B). A similar AMOC-related GMOT warming also occurs during the YD but with a smaller magnitude due to its shorter duration (Fig. 4, A and B). As such, during the early deglaciation, nearly 50% (1.8°C) of the total simulated GMOT increase (3.4°C) is contributed by the AMOC reduction (Fig. 4B). By comparison, the HSE in the ICE + ORB + GHG experiment increases monotonically from ~0.2 during HS1 to ~0.9 during the MH (Fig. 4C), reflecting a weaker ocean interior warming from the LGM to MH than in iTRACE (Fig. 2, A to C, versus fig. S8, A to C). Given its small effect on GMSST but strong and long-term effect on GMOT, a decreased AMOC is thus highly effective in enhancing the HSE during the entire deglaciation, with particularly large increases during HS1 and the YD (Figs. 1D and 4C). Moreover, in iTRACE, the heat accumulated in the global ocean by the AMOC reduction lasts for >4000 years after the termination of major meltwater forcing (Fig. 1A), with a notable residual of 0.5°C in ∆GMOT persisting until the MH (Fig. 4B).

DISCUSSION

Our iTRACE simulations and new proxy reconstructions of ∆DOT demonstrate that the deglacial ocean HSE was substantially enhanced by subsurface warming at intermediate depths through the effects of increasing GHGs and retreating ice sheets on surface warming as well as by a weakening of the AMOC in response to meltwater forcing. We emphasize several implications of these results to our understanding of changes in MOT. First, we note that the warming in the Indo-Pacific Ocean around 2000 m is closest to the global mean value (fig. S2B), suggesting that reconstructions of past ∆DOT from similar locations would closely monitor ∆GMOT. Second, our results showing that millennial-scale AMOC variability and its lingering effect on ocean warming substantially enhanced ocean HSE to >1 suggest that persistent AMOC variability throughout a glacial cycle would do the same, whereas reduced variability would cause HSE to be <1. Third, our analyses suggest that different patterns of surface warming associated with different timescales may cause a substantially different ocean HSE. For example, maximum ocean heat uptake today is occurring at mid-latitudes (5), but under the short-term anthropogenic forcing, the dominant large-scale surface warming is found in subtropical latitudes (43, 56), which dominates the GMSST change given the large surface area of the subtropical ocean. Because transmission of this subtropical warming to depth is limited by the shallow subtropical cell and can only warm to depths <500 m, it has damped ocean heat uptake efficiency (43) and thus contributed little to enhancing HSE at this early stage, leading to an HSE of <<1 (fig. S14). On millennial and longer (near-equilibration) timescales such as during the last deglaciation, if we assume that areas of strong surface warming remained colocated with globally dominant ventilation regions in mid- to subpolar latitudes like in iTRACE, then we should expect ventilation of the thicker intermediate-depth waters to enhance HSE. As a result, the ocean will sequester more heat from the atmosphere, potentially slowing the rate of atmospheric warming. Our study therefore suggests that the pattern effect not only influences climate sensitivity through atmospheric feedbacks (57) but also influences ocean heat uptake/storage and its efficiency which can substantially modulate transient climate change (58) and the Earth’s energy budget (59). In the future, if global warming continues as in the case of the longer-term (centennial to millennial timescale) warming investigated by Rugenstein et al. (10), the HSE could be further enhanced (fig. S14) as the strong surface warming extends further poleward to subpolar latitudes where winter sea ice coverage will be substantially reduced.

MATERIALS AND METHODS

Proxy ∆DOT records

We reconstructed deglacial ∆DOT from 119 globally distributed, well-dated benthic δ18O (δ18Ob) records. We selected 108 of these records from the 287 records that are archived in the OC3 database (35) based on their resolution and 14C age control. We also used the following δ18Ob records that are not in the OC3 database: ENAM93-21 (60), EW9302-24GGC (61), EW9302-25GGC (61), EW9302-26GGC (61), EW9302-2JPC (50), MD95-2010 (62), RAPID-10-1P (62), U1308 (62), KNR166-2-26JPC (63), KNR197-3-46CDH (63), and M35003-4 (64).

The average resolution of the 119 records is 520 ± 290 years, with 80 records having a resolution <600 years, and we interpolated all data to a 500-year resolution to reflect the average resolution of the data. Age models of all records used here are based on calibrated 14C ages as opposed to some records in the database with age models derived by alignment to the LR04 δ18Ob stack. For most records, where the age models of records are based on several 14C calibration schemes, we used the one based on the most recent calibration. For the following cores, we used the earlier Jonkers age model in the OC3 database as it was more consistent with other constraints: ODP-162-984, GIK23415_9, and SO213-2-82-1.

Benthic δ18O records (δ18Ob) reflect some combination of the temperature of foraminifera calcification (δ18OT) and the δ18O of seawater (δ18Osw), which in turn reflects global land ice volume and salinity. We derived δ18OT by subtracting changes in δ18Osw associated with the decrease in global ice volume from δ18Ob and assessed possible additional salinity contributions by comparison to independently derived temperature reconstructions (see below). We used the LGM constraints of a 1‰ increase in δ18Osw as established by pore water measurements in deep-sea cores (65) and a sea-level lowering of −134 m (66) to derive the relationship 0.0075‰ m−1 and applied it to a modeled sea-level reconstruction for the last 21 kyr (66). The uncertainties on these constraints are small, and so we do not propagate them into our results as they would not affect our conclusions. We then converted the resulting δ18OT to temperature (∆DOT) using the calibration from Marchitto et al. (67) for Cibicidoides (n = 118) and Uvigerina (n = 1):

T = [0.245–SQRT (0.045461 + 0.0044 δ18OT)]/0.0022.

The ∆DOT time series was then normalized to 20 ka. We calculated the uncertainty on individual site records through a Monte Carlo simulation by using the average uncertainty (300 years) of the calibrated age models at each time interval (0 to 20 ka at 500-year resolution). We generated 1000 age models for each site by randomly adding or subtracting a Gaussian distribution of 300 years to each time interval and then created 1000 temperature records using the 1000 random age models. We then took the average and SD (1σ) of the 1000 temperature records to calculate the Monte Carlo version of our temperature record.

To assess the vertical temperature structure simulated by iTRACE, we constructed regional ∆DOT stacks from sites occurring in the areas where simulated ocean temperatures show lesser, similar, or greater amounts of change than GMSST (Fig. 2, B to F). We used the jackknife procedure to recalculate each regional stack 1000 times after randomly removing 50% of the records in each iteration, with the SD of the resulting stacks being taken as the 1σ age uncertainty. We applied an analytical uncertainty of 0.04‰ to each stack, converted to temperature using 0.25 ‰ °C−1 (68), resulting in an uncertainty of 0.16°C. Similarly, we applied a δ18Osw error of 0.1‰ (65) to each stack corresponding to a 0.4°C uncertainty. We summed these uncertainties in quadrature to derive the total uncertainty on each regional stack.

We also calculated regional ∆DOT stacks and their uncertainties (as described above) for the Atlantic and Indo-Pacific Oceans (there are no data from the Southern Ocean), weighted each by the respective volume of their ocean basins and normalized relative to the total ocean volume (26 and 74%, respectively) (fig. S2B), and then stacked them to generate the ∆GMDOT stack. Similarly, to generate the global stack uncertainty, we weighted the jackknife regional uncertainties for each regional stack and added these in quadrature to the δ18Ob analytical error of 0.16°C and the δ18Osw of porewater error of 0.4°C to derive the uncertainty on the entire stack.

We find that independently derived ∆DOT reconstructions that are available at several of our sites show similar changes as in our reconstructions (fig. S4). Given our inferred contribution of δ18Osw from global ice sheets, this agreement suggests a negligible contribution from local salinity changes at these sites. We also find that our basin-weighted ∆GMDOT stack is in good agreement with the noble gas–based ∆GMOT reconstructions (Fig. 1C).

Shackleton et al. (7) converted the noble gas–based ∆GMOT reconstructions to δ18OT and estimated the δ18Osw component from land ice as Δδ18Oice volume = δ18Oice/(D/ΔSL − 1), where δ18Oice is their assumed mean ice sheet δ18O (−30 ± 2‰), D is the modern mean ocean depth, and ΔSL is the global sea-level change from Lambeck et al. (66). When applied to the last deglaciation, Δδ18Oice volume is 1.02‰ at 20 ka or the same as the pore water reconstruction (65) used here. Summing their δ18OT and Δδ18O reconstructions gave a δ18O record (predicted) that is the same, within uncertainty, as the LR04 δ18Ob record (observed) from 0 to 20 ka. The agreement of our independently derived global δ18OT and δ18Osw reconstructions with those in Shackleton et al. (7) thus indicates that any additional contribution from salinity changes at the global scale is negligible.

Model and simulations

iTRACE is a transient simulation of global climate and water isotopes during the last deglaciation (30) conducted with iCESM1 (29). The resolution of atmosphere and land surface is a nominal 2° (1.9° in latitude and 2.5° in longitude) with 30 vertical levels in the atmosphere; the resolution of the ocean and sea ice is a nominal 1° with 60 vertical levels in the ocean. Starting from the LGM (20 ka), four forcing factors, namely, the ICE, solar insolation associated with ORB, GHG, and MWF, were applied additively to simulate the deglacial climate change. Continental ice sheet configuration following the ICE-6G reconstruction (69) was modified every 1000 year during the simulation. GHG concentrations (CO2, CH4, and N2O) were prescribed based on ice core reconstructions (70). MWF was applied on the basis of sea-level reconstructions, largely following the scheme used in TraCE-21k (32). Assuming linearity of the responses, the difference between different pairs of the factorized-forcing simulations can be used to estimate the impact of each individual forcing approximately. For example, (ICE + ORB + GHG + MWT) − (ICE + ORB + GHG) isolates the effects of meltwater forcing (fig. S13). The set of iTRACE simulations have been shown to be able to reproduce many major features of the observed change of global climate and water masses during the early deglaciation (from 20 to 11 ka) (30, 55, 71). Here, we used the model output of 20 to 6 ka, which is now available.

Besides iTRACE, another transient simulation TraCE-21k, the predecessor of iTRACE (32, 72), and multiple equilibrium model simulations are also investigated. TraCE-21k was performed with a lower-resolution version of the Community Climate System Model version 3 (CCSM3). Similar to iTRACE, TraCE-21k is also forced by deglacial changes in ice sheets, solar insolation, GHGs, and MWF and captures many important features of deglacial climate in its all forcing experiment (19, 32). More details of TraCE-21k simulation can be found in (72). For equilibrium simulations, we use “PI,” “MH,” and “LGM” simulations performed with iCESM1 and seven PMIP3 models (34). LongRunMIP simulations (10) are also investigated to show the long-term ocean HSE under future anthropogenic warming.

Calculation of GMOT, GMSST, and AMOC intensity

For model simulations, the three-dimensional ocean temperature fields are interpolated from their corresponding horizontal grids to 1° rectilinear latitude-longitude grids. The GMOT is calculated as the volume-weighted average over the entire ocean. GMSST is calculated as the area-weighted average over the ocean surface. Both GMOT and GMSST calculations use all ocean cells available including open sea cells and those covered by sea ice. AMOC intensity is defined as the maximum overturning stream function below 300 m over the section of 20° to 50°N in the Atlantic.

Further model-data comparison

A further temperature-salinity diagram comparison of sites from Adkins et al. (73) indicates that our iTRACE model faithfully reproduces the homogeneous glacial deep ocean temperature that around freezing point and larger deglacial warming (>4°C) in the deep Atlantic than in the deep Pacific and Southern Ocean (<2.5°C), similar to iTRACE (Fig. 2, B and C) and TraCE-21k (74).

There are some model-data discrepancies. Compared with observations, the evolution of the simulated GMSST experiences sharper peaks and troughs that are related to strong AMOC variability (Figs. 1B and 4A). This high-frequency surface temperature variability, however, has little imprint in the broad global ocean temperature change, as indicated by the consistent GMOTs in iTRACE and observations. For GMSST records, the absence of these sharp peaks and troughs can be caused, partly, by the low temporal resolution and chronology uncertainty.

For PMIP3 models, we note that there are differences in the detailed ocean warming pattern, especially in the North Atlantic (fig. S3). This can be related to model differences in simulating the deep convection and the AMOC transport. Specifically, AMOC strength (relative to present day) at the LGM is highly uncertain among PMIP models (75, 76) which will influence the temperature response at the LGM. However, as the North Atlantic only occupies a small fraction in volume, such difference makes only little differences to GMOT change.

Acknowledgments

We acknowledge the high-performance computing support from NCAR’s Computational and Information Systems Laboratory (CISL), Laoshan Laboratory, and the National Key Scientific and Technological Infrastructure project “Earth System Numerical Simulation Facility” (EarthLab).

Funding: This work was supported by National Natural Science Foundation of China (42106013 to Chenyu Zhu and 42130604 to L.W.), Science and Technology Innovation Project of Laoshan Laboratory (LSKJ202203303 to Chenyu Zhu, L.W., and J.L.), US National Science Foundation (AGS2202860 to Z.L. and OPP2103032 to P.U.C.), and US National Science Foundation Graduate Research Fellowship (2234662 to S.S.).

Author contributions: Z.L. and P.U.C. conceived the study. Chenyu Zhu and S.S. performed the analysis. Chenyu Zhu drafted the manuscript with Z.L. Chenyu Zhu, C.H., L.W., and J.L. performed the iTRACE simulations. S.S. and P.U.C. performed the proxy data analysis. Chenguang Zhu and L.L. assisted in the model analysis. All authors discussed and revised the manuscript.

Competing interests: The authors declare that they have no competing interests.

Data and materials availability: All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. iTRACE (20 to 11 ka) and TraCE-21k data used in the present study are publicly available at https://earthsystemgrid.org/dataset/ucar.cgd.ccsm4.iTRACE.html and https://earthsystemgrid.org/project/trace.html, respectively. PMIP3 model outputs for LGM, MH, and PI simulations can be downloaded at https://aims2.llnl.gov/search/cmip5/. LongRunMIP simulation dataset is available at http://longrunmip.org/. The new regional and global ocean temperature stacks are provided in the Supplementary Materials. iTRACE simulation data for Figs. 1 to 4 are available at https://figshare.com/account/home#/projects/195203.

Supplementary Materials

The PDF file includes:

Figs. S1 to S14

Table S1

Legend for data S1

References

Other Supplementary Material for this manuscript includes the following:

Data S1
==== Refs
REFERENCES AND NOTES

1 J. M. Gregory, R. J. Stouffer, S. C. B. Raper, P. A. Stott, N. A. Rayner, An observationally based estimate of the climate sensitivity. J. Climate 15 , 3117–3121 (2002).
2 J. Hansen, L. Nazarenko, R. Ruedy, M. Sato, J. Willis, A. Del Genio, D. Koch, A. Lacis, K. Lo, S. Menon, T. Novakov, J. Perlwitz, G. Russell, G. A. Schmidt, N. Tausnev, Earth’s energy imbalance: Confirmation and implications. Science 308 , 1431–1435 (2005).15860591
3 S. Levitus, J. I. Antonov, T. P. Boyer, O. K. Baranova, H. E. Garcia, R. A. Locarnini, A. V. Mishonov, J. R. Reagan, D. Seidov, E. S. Yarosh, M. M. Zweng, World ocean heat content and thermosteric sea level change (0–2000 m), 1955–2010. Geophys. Res. Lett. 39 , 10.1029/2012GL051106 (2012).
4 M. Rhein, S. R. Rintoul, S. Aoki, E. Campos, D. Chambers, R. A. Feely, S. Gulev, G. C. Johnson, S. A. Josey, A. Kostianoy, C. Mauritzen, D. Roemmich, L. D. Talley, F. Wang, “Observations: Ocean” in Climate Change 2013: The Physical Science Basis. Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change, T. F. Stocker, D. Qin, G.-K. Plattner, M. Tignor, S. K. Allen, J. Boschung, A. Nauels, Y. Xia, V. Bex, P. M. Midgley, Eds. (Cambridge Univ. Press, 2013), pp. 255–315.
5 L. Cheng, J. Abraham, K. E. Trenberth, T. Boyer, M. E. Mann, J. Zhu, F. Wang, F. Yu, R. Locarnini, J. Fasullo, F. Zheng, Y. Li, B. Zhang, L. Wan, X. Chen, D. Wang, L. Feng, X. Song, Y. Liu, F. Reseghetti, S. Simoncelli, V. Gouretski, G. Chen, A. Mishonov, J. Reagan, K. von Schuckmann, Y. Pan, Z. Tan, Y. Zhu, W. Wei, G. Li, Q. Ren, L. Cao, Y. Lu, New record ocean temperatures and related climate indicators in 2023. Adv. Atmos. Sci. 41 , 1068–1082 (2024).
6 E. D. Galbraith, T. M. Merlis, J. B. Palter, Destabilization of glacial climate by the radiative impact of Atlantic meridional overturning circulation disruptions. Geophys. Res. Lett. 43 , 8214–8221 (2016).
7 S. Shackleton, A. Seltzer, D. Baggenstos, L. Lisiecki, Benthic δ18O records Earth’s energy imbalance. Nat. Geosci. 16 , 797–802 (2023).
8 K. von Schuckmann, A. Minière, F. Gues, F. J. Cuesta-Valero, G. Kirchengast, S. Adusumilli, F. Straneo, M. Ablain, R. P. Allan, P. M. Barker, H. Beltrami, A. Blazquez, T. Boyer, L. Cheng, J. Church, D. Desbruyeres, H. Dolman, C. M. Domingues, A. García-García, D. Giglio, J. E. Gilson, M. Gorfer, L. Haimberger, M. Z. Hakuba, S. Hendricks, S. Hosoda, G. C. Johnson, R. Killick, B. King, N. Kolodziejczyk, A. Korosov, G. Krinner, M. Kuusela, F. W. Landerer, M. Langer, T. Lavergne, I. Lawrence, Y. Li, J. Lyman, F. Marti, B. Marzeion, M. Mayer, A. H. MacDougall, T. McDougall, D. P. Monselesan, J. Nitzbon, I. Otosaka, J. Peng, S. Purkey, D. Roemmich, K. Sato, K. Sato, A. Savita, A. Schweiger, A. Shepherd, S. I. Seneviratne, L. Simons, D. A. Slater, T. Slater, A. K. Steiner, T. Suga, T. Szekely, W. Thiery, M.-L. Timmermans, I. Vanderkelen, S. E. Wjiffels, T. Wu, M. Zemp, Heat stored in the Earth system 1960–2020: Where does the energy go? Earth Syst. Sci. Data 15 , 1675–1709 (2023).
9 S. K. Gulev, P. W. Thorne, J. Ahn, F. J. Dentener, C. M. Domingues, S. Gerland, D. Gong, D. S. Kaufman, H. C. Nnamchi, J. Quaas, J. A. Rivera, S. Sathyendranath, S. L. Smith, B. Trewin, K. von Schuckmann, R. S. Vose, “Changing state of the climate system” in Climate Change 2021: The Physical Science Basis. Contribution of Working Group I to the Sixth Assessment Report of the Intergovernmental Panel on Climate Change, V. Masson-Delmotte, P. Zhai, A. Pirani, S. L. Connors, C. Péan, S. Berger, N. Caud, Y. Chen, L. Goldfarb, M. I. Gomis, M. Huang, K. Leitzell, E. Lonnoy, J. B. R. Matthews, T. K. Maycock, T. Waterfield, O. Yelekçi, R. Yu, B. Zhou, Eds. (Cambridge Univ. Press, 2021), pp. 287−422.
10 M. Rugenstein, J. Bloch-Johnson, A. Abe-Ouchi, T. Andrews, U. Beyerle, L. Cao, T. Chadha, G. Danabasoglu, J.-L. Dufresne, L. Duan, M.-A. Foujols, T. Frölicher, O. Geoffroy, J. Gregory, R. Knutti, C. Li, A. Marzocchi, T. Mauritsen, M. Menary, E. Moyer, L. Nazarenko, D. Paynter, D. Saint-Martin, G. A. Schmidt, A. Yamamoto, S. Yang, LongRunMIP—Motivation and design for a large collection of millennial-length GCM simulations. Bull. Am. Meteorol. Soc. 100 , 2551–2570 (2019).
11 Y. Rosenthal, B. K. Linsley, D. W. Oppo, Pacific Ocean heat content during the past 10,000 years. Science 342 , 617–620 (2013).24179224
12 G. Gebbie, P. Huybers, The Little Ice Age and 20th-century deep Pacific cooling. Science 363 , 70–74 (2019).30606843
13 A. Bagnell, T. DeVries, 20th century cooling of the deep ocean contributed to delayed acceleration of Earth’s energy imbalance. Nat. Commun. 12 , 4604 (2021).34326319
14 H. Elderfield, P. Ferretti, M. Greaves, S. Crowhurst, I. N. McCave, D. Hodell, A. M. Piotrowski, Evolution of ocean temperature and ice volume through the mid-Pleistocene climate transition. Science 337 , 704–709 (2012).22879512
15 J. D. Shakun, D. W. Lea, L. E. Lisiecki, M. E. Raymo, An 800-kyr record of global surface ocean δ18O and implications for ice volume-temperature coupling. Earth Planet. Sci. Lett. 426 , 58–68 (2015).
16 E. J. Rohling, G. L. Foster, T. M. Gernon, K. M. Grant, D. Heslop, F. D. Hibbert, A. P. Roberts, J. Yu, Comparison and synthesis of sea-level and deep-seatemperature variations over the past 40 million years. Rev. Geophys. 60 , e2022RG000775 (2022).
17 B. Bereiter, S. Shackleton, D. Baggenstos, K. Kawamura, J. Severinghaus, Mean global ocean temperatures during the last glacial transition. Nature 553 , 39–44 (2018).29300008
18 M. Häberli, D. Baggenstos, J. Schmitt, M. Grimmer, A. Michel, T. Kellerhals, H. Fischer, Snapshots of mean ocean temperature over the last 700,000 years using noble gases in the EPICA Dome C ice core. Clim. Past 17 , 843–867 (2021).
19 J. D. Shakun, P. U. Clark, F. He, S. A. Marcott, A. C. Mix, Z. Liu, B. Otto-Bliesner, A. Schmittner, E. Bard, Global warming preceded by increasing carbon dioxide concentrations during the last deglaciation. Nature 484 , 49–54 (2012).22481357
20 S. A. Marcott, J. D. Shakun, P. U. Clark, A. C. Mix, A reconstruction of regional and global temperature for the past 11,300 years. Science 339 , 1198–1201 (2013).23471405
21 J. E. Tierney, J. Zhu, J. King, S. B. Malevich, G. J. Hakim, C. J. Poulsen, Glacial cooling and climate sensitivity revisited. Nature 584 , 569–573 (2020).32848226
22 P. U. Clark, J. D. Shakun, Y. Rosenthal, P. Köhler, P. J. Bartlein, Global and regional temperature change over the past 4.5 million years. Science 383 , 884–890 (2024).38386742
23 A. M. Seltzer, P. W. Davidson, S. A. Shackleton, D. P. Nicholson, S. Khatiwala, Global ocean cooling of 2.3°C during the Last Glacial Maximum. Geophys. Res. Lett. 51 , e2024GL108866 (2024).
24 C. Emiliani, Temperatures of Pacific bottom waters and polar superficial waters during the Tertiary. Science 119 , 853–855 (1954).17751537
25 J. Hansen, M. Sato, G. Russell, P. Kharecha, Climate sensitivity, sea level and atmospheric carbon dioxide. Phil. Trans. R. Soc. A. 371 , 20120294 (2013).24043864
26 J. E. Hansen, M. Sato, L. Simons, L. S. Nazarenko, I. Sangha, P. Kharecha, J. C. Zachos, K. von Schuckmann, N. G. Loeb, M. B. Osman, Q. Jin, G. Tselioudis, E. Jeong, A. Lacis, R. Ruedy, G. Russell, J. Cao, J. Li, Global warming in the pipeline. Oxford Open Clim. Change 3 , kgad008 (2023).
27 N. Maffezzoli, P. Vallelonga, R. Edwards, A. Saiz-Lopez, C. Turetta, H. A. Kjær, C. Barbante, B. Vinther, A. Spolaor, A 120 000-year record of sea ice in the North Atlantic? Clim. Past 15 , 2031–2051 (2019).
28 X. Crosta, K. E. Kohfeld, H. C. Bostock, M. Chadwick, A. Du Vivier, O. Esper, J. Etourneau, J. Jones, A. Leventer, J. Müller, R. H. Rhodes, C. S. Allen, P. Ghadi, N. Lamping, C. B. Lange, K.-A. Lawler, D. Lund, A. Marzocchi, K. J. Meissner, L. Menviel, A. Nair, M. Patterson, J. Pike, J. G. Prebble, C. Riesselman, H. Sadatzki, L. C. Sime, S. K. Shukla, L. Thöle, M.-E. Vorrath, W. Xiao, J. Yang, Antarctic sea ice over the past 130 000 years—Part 1: A review of what proxy records tell us. Clim. Past 18 , 1729–1756 (2022).
29 E. C. Brady, S. Stevenson, D. Bailey, Z. Liu, D. Noone, J. Nusbaumer, B. L. Otto-Bliesner, C. Tabor, R. Tomas, T. Wong, J. Zhang, J. Zhu, The connected isotopic water cycle in the Community Earth System Model version 1. J. Adv. Model. Earth Syst. 11 , 2547–2566 (2019).
30 C. He, Z. Liu, B. L. Otto-Bliesner, E. C. Brady, C. Zhu, R. Tomas, P. U. Clark, J. Zhu, A. Jahn, S. Gu, J. Zhang, J. Nusbaumer, D. Noone, H. Cheng, Y. Wang, M. Yan, Y. Bao, Hydroclimate footprint of pan-Asian monsoon water isotope during the last deglaciation. Sci. Adv. 7 , eabe2611 (2021).33523950
31 M. B. Osman, J. E. Tierney, J. Zhu, R. Tardif, G. J. Hakim, J. King, C. J. Poulsen, Globally resolved surface temperatures since the Last Glacial Maximum. Nature 599 , 239–244 (2021).34759364
32 Z. Liu, B. L. Otto-Bliesner, F. He, E. C. Brady, R. Tomas, P. U. Clark, A. E. Carlson, J. Lynch-Stieglitz, W. Curry, E. Brook, D. Erickson, R. Jacob, J. Kutzbach, J. Cheng, Transient simulation of last deglaciation with a new mechanism for Bølling-Allerød warming. Science 325 , 310–314 (2009).19608916
33 J. R. Alder, S. W. Hostetler, Global climate simulations at 3000-year intervals for the last 21 000 years with the GENMOM coupled atmosphere–ocean model. Clim. Past 11 , 449–471 (2015).
34 P. Braconnot, S. P. Harrison, M. Kageyama, P. J. Bartlein, V. Masson-Delmotte, A. Abe-Ouchi, B. Otto-Bliesner, Y. Zhao, Evaluation of climate models using palaeoclimatic data. Nat. Clim. Chang. 2 , 417–424 (2012).
35 J. Muglia, S. Mulitza, J. Repschläger, A. Schmittner, L. Lembke-Jene, L, Lisiecki, A. Mix, R. Saraswat, E. Sikes, C. Waelbroeck, J. Gottschalk, J. Lippold, D. Lund, G. Martinez-Mendez, E. Michel, F. Muschitiello, S. Naik, Y. Okazaki, L. Stott, A. Voelker, N. Zhao, A global synthesis of high-resolution stable isotope data from benthic foraminifera of the last deglaciation. Sci. Data. 10 , 131 (2023).36899009
36 L. E. Lisiecki, J. V. Stern, Regional and global benthic δ18O stacks for the last glacial cycle. Paleoceanography 31 , 1368–1394 (2016).
37 CLIMAP Project Members, The surface of the ice-age Earth: Quantitative geologic evidence is used to reconstruct boundary conditions for the climate 18,000 years ago. Science 191 , 1131–1137 (1976).17781630
38 M. Kucera, M. Weinelt, T. Kiefer, U. Pflaumann, A. Hayes, M. Weinelt, M. Chen, A. C. Mix, T. T. Barrows, E. Cortijo, J. Duprat, S. Juggins, C. Waelbroeck, Reconstruction of the glacial Atlantic and Pacific sea-surface temperatures from assemblages of planktonic foraminifera: Multi-technique approach based on geographically constrained calibration datasets. Quat. Sci. Rev. 24 , 951–998 (2005).
39 MARGO Project Members, Constraints on the magnitude and patterns of ocean cooling at the Last Glacial Maximum. Nat. Geosci. 2 , 127–132 (2009).
40 J. D. Annan, J. C. Hargreaves, A new global reconstruction of temperature changes at the Last Glacial Maximum. Clim. Past 9 , 367–376 (2013).
41 J. D. Annan, J. C. Hargreaves, T. Mauritsen, A new global surface temperature reconstruction for the Last Glacial Maximum. Clim. Past 18 , 1883–1896 (2022).
42 L. Back, K. Russ, Z. Liu, K. Inoue, J. Zhang, B. L. Otto-Bliesner, Global hydrological cycle response to rapid and slow global warming. J. Climate 22 , 8781–8786 (2013).
43 E. Newsom, L. Zanna, S. Khatiwala, J. M. Gregory, The influence of warming patterns on passive ocean heat uptake. Geophys. Res. Lett. 47 , e2020GL088429 (2020).
44 L. Zanna, S. Khatiwala, J. M. Gregory, J. Ison, P. Heimbach, Global reconstruction of historical ocean heat storage and transport. Proc. Natl. Acad. Sci. U.S.A. 116 , 1126–1131 (2019).30617081
45 J. F. McManus, R. Francois, J. Gherardi, L. D. Keigwin, S. Brown-Leger, Collapse and rapid resumption of Atlantic meridional circulation linked to deglacial climate changes. Nature 428 , 834–837 (2004).15103371
46 J. Lippold, F. Pöppelmeier, F. Süfke, M. Gutjahr, T. J. Goepfert, P. Blaser, O. Friedrich, J. M. Link, L. Wacker, S. Rheinberger, S. L. Jaccard, Constraining the variability of the Atlantic meridional overturning circulation during the Holocene. Geophys. Res. Lett. 46 , 11338–11346 (2019).
47 W. S. Broecker, Paleocean circulation during the last deglaciation: A bipolar seesaw? Paleoceanography 13 , 119–121 (1998).
48 P. U. Clark, N. G. Pisias, T. F. Stocker, A. J. Weaver, The role of the thermohaline circulation in abrupt climate change. Nature 415 , 863–869 (2002).11859359
49 T. F. Stocker, S. J. Johnsen, A minimum thermodynamic model for the bipolar seesaw. Paleoceanography 18 , 1087 (2003).
50 S. A. Marcott, P. U. Clark, L. Padman, G. P. Klinkhammer, S. R. Springer, Z. Liu, B. L. Otto-Bliesner, A. E. Carlson, A. Ungerer, J. Padman, F. He, J. Cheng, A. Schmittner, Ice-shelf collapse from subsurface warming as a trigger for Heinrich events. Proc. Natl. Acad. Sci. U.S.A. 108 , 13415–13419 (2011).21808034
51 S. Weldeab, T. Friedrich, A. Timmermann, R. R. Schneider, Strong middepth warming and weak radiocarbon imprints in the equatorial Atlantic during Heinrich 1 and Younger Dryas. Paleoceanography 31 , 1070–1082 (2016).
52 S. Barker, P. Diz, Timing of the descent into the last Ice Age determined by the bipolar seesaw. Paleoceanography 29 , 489–507 (2014).
53 J. Mignot, A. Ganopolski, A. Levermann, Atlantic subsurface temperatures: Response to a shutdown of the overturning circulation and its recovery. J. Climate 20 , 4884–4898 (2007).
54 C. He, Z. Liu, J. Zhu, J. Zhang, S. Gu, B. L. Otto-Bliesner, B. Esther, C. Zhu, Y. Jin, J. Sun, North Atlantic subsurface temperature response controlled by effective freshwater input in “Heinrich” events. Earth Planet. Sci. Lett. 539 , 116247 (2020).
55 C. Zhu, J. Zhang, Z. Liu, B. L. Otto-Bliesner, C. He, E. Brady, R. Tomas, Q. Wen, Q. Li, C. Zhu, S. Zhang, L. Wu, Antarctic warming during Heinrich Stadial 1 in a transient isotope-enabled deglacial simulation. J. Climate 35 , 3753–3765 (2022).
56 H. Yang, G. Lohmann, C. Stepanek, Q. Wang, R. X. Huang, X. Shi, J. Liu, D. Chen, X. Wang, Y. Zhong, Q. Yang, Y. Bao, J. Müller, Satellite-observed strong subtropical ocean warming as an early signature of global warming. Commun. Earth Environ. 4 , 178 (2023).
57 T. Andrews, J. M. Gregory, D. Paynter, L. G. Silvers, C. Zhou, T. Mauritsen, M. J. Webb, K. C. Armour, P. M. Forster, H. Titchner, Accounting for changing temperature patterns increases historical estimates of climate sensitivity. Geophys. Res. Lett. 45 , 8490–8499 (2018).
58 M. Winton, K. Takahashi, I. M. Held, Importance of ocean heat uptake efficacy to transient climate change. J. Climate 23 , 2333–2344 (2010).
59 D. Baggenstos, M. Häberli, J. Schmitt, S. A. Shackleton, B. Birner, J. P. Severinghaus, T. Kellerhals, H. Fischer, Earth’s radiative imbalance from the Last Glacial Maximum to the present. Proc. Natl. Acad. Sci. U.S.A. 116 , 14881–14886 (2019).31285336
60 T. L. Rasmussen, E. Thomsen, T. C. E. van Weering, L. Labeyrie, Rapid changes in surface and deep water conditions at the Faeroe Margin during the last 58,000 years. Paleoceanography 11 , 757–772 (1996).
61 D. W. Oppo, W. B. Curry, J. F. McManus, What do benthic δ13C and δ18O data tell us about Atlantic circulation during Heinrich Stadial 1? Paleoceanography 30 , 353–368 (2015).
62 T. Dokken, E. Jansen, Rapid changes in the mechanism of ocean convection during the last glacial period. Nature 401 , 458–461 (1999).
63 D. W. Oppo, W. Lu, K.-F. Huang, N. E. Umling, W. Guo, J. Yu, W. B. Curry, T. M. Marchitto, S. Wang, Deglacial temperature and carbonate saturation state variability in the tropical Atlantic at Antarctic Intermediate Water depths. Paleoceanography 38 , e2023PA004674 (2023).
64 C. Rühlemann, S. Mulitza, G. Lohmann, A. Paul, M. Prange, G. Wefer, Intermediate depth warming in the tropical Atlantic related to weakened thermohaline circulation: Combining paleoclimate data and modeling results for the last deglaciation. Paleoceanography 19 , PA1025 (2004).
65 D. P. Schrag, J. F. Adkins, K. McIntyre, J. L. Alexander, D. A. Hodell, C. D. Charles, J. F. McManus, The oxygen isotopic composition of seawater during the Last Glacial Maximum. Quat. Sci. Rev. 21 , 331 (2002).
66 K. Lambeck, H. Rouby, A. Purcell, Y. Sun, M. Sambridge, Sea level and global ice volumes from the Last Glacial Maximum to the Holocene. Proc. Natl. Acad. Sci. U.S.A. 111 , 15296–15303 (2014).25313072
67 T. M. Marchitto, W. B. Curry, J. Lynch-Stieglitz, S. P. Bryan, K. M. Cobb, D. C. Lund, Improved oxygen isotope temperature calibrations for cosmopolitan benthic foraminifera. Geochim. Cosmochim. Acta 130 , 1–11 (2014).
68 N. J. Shackleton, Attainment of isotopic equilibrium between ocean water and the benthonic foraminifera genus Uvigerina: Isotopic changes in the ocean during the last glacial. Colloques Internationaux du C.N.R.S 219 , 203–209 (1974).
69 W. R. Peltier, D. F. Argus, R. Drummond, Space geodesy constrains ice age terminal deglaciation: The global ICE-6G_C (VM5a) model. J. Geophys. Res. Solid Earth 120 , 450–487 (2015).
70 D. Lüthi, M. L. Floch, B. Bereiter, T. Blunier, J. Barnola, U. Siegenthaler, D. Raynaud, J. Jouzel, H. Fischer, K. Kawamura, T. F. Stocker, High-resolution carbon dioxide concentration record 650,000–800,000 years before present. Nature 453 , 379–382 (2008).18480821
71 S. Gu, Z. Liu, D. W. Oppo, J. Lynch-Stieglitz, A. Jahn, J. Zhang, L. Wu, Assessing the potential capability of reconstructing glacial Atlantic water masses and AMOC using multiple proxies in CESM. Earth Planet. Sci. Lett. 541 , 116294 (2020).
72 F. He, “Simulating transient climate evolution of the last deglaciation with CCSM3,” thesis, University of Wisconsin-Madison (2011), pp. 171.
73 J. F. Adkins, K. McIntyre, D. P. Schrag, The salinity, temperature and δ18O of the glacial deep ocean. Science 298 , 1769–1773 (2002).12459585
74 Z. Liu, S. Shin, R. Webb, E. Lewis, B. L. Otto-Bliesner, Atmospheric CO2 forcing on glacial thermohaline and climate. Geophys. Res. Lett. 32 , L02706 (2005).
75 B. L. Otto-Bliesner, C. D. Hewitt, T. M. Marchitto, E. Brady, A. Abe-Ouchi, M. Crucifix, S. Murakami, S. L. Weber, Last Glacial Maximum ocean thermohaline circulation: PMIP2 model intercomparisons and data constraints. Geophys. Res. Lett. 34 , L12706 (2007).
76 C. Zhu, Z. Liu, S. Zhang, L. Wu, Global oceanic overturning circulation forced by the competition between greenhouse gases and continental ice sheets during the last deglaciation. J. Climate 34 , 7555–7570 (2021).
77 S. G. Valley, J. Lynch-Stieglitz, T. M. Marchitto, Intermediate water circulation changes in the Florida Straits from a 35 ka record of Mg/Li-derived temperature and Cd/Ca-derived seawater cadmium. Earth Planet. Sci. Lett. 523 , 115692 (2019).
78 N. E. Umling, D. W. Oppo, P. Chen, J. Yu, Z. Liu, M. Yan, G. Gebbie, D. C. Lund, K. R. Pietro, Z. D. Jin, K.-F. Huang, K. B. Costa, F. A. L. Toledo, Atlantic circulation and ice sheet influences on upper South Atlantic temperatures during the Last Deglaciation. Paleoceanogr. Paleoclimatol. 34 , 990–1005 (2019).
