Vinje Eirik

Ulysses Galaasen

Augustin Ninnemann

Nil Kessler

Yair Irvali

Jerry Rosenthal

Nathaëlle Tjiputra

Didier M Bouttes

Helga Roche

David Kleiven

Hodell

Eirik Vinje Galaasen
email: eirik.galaasen@uib.no

Ulysses S Ninnemann

Augustin Kessler

Nil Irvalı

Yair Rosenthal

Jerry Tjiputra

Nathaëlle Bouttes

Didier M Roche

Helga ( Kikki

F Kleiven

David A Hodell

Interglacial instability of North Atlantic Deep Water ventilation

ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

One Sentence Summary: Large and frequent changes in ocean ventilation during recent warm periods question its modern stability.

Main Text: Atlantic meridional overturning circulation (AMOC) and North Atlantic Deep

Water (NADW) ventilation represents a low probability, high impact tipping point (1) in the climate system with implications for the distribution and sequestration of anthropogenic CO 2 and heat, and for Atlantic-wide patterns of climate and sea level [START_REF] Buckley | Observations, inferences, and mechanisms of the Atlantic Meridional Overturning Circulation: A review[END_REF][START_REF] Sabine | The oceanic sink for anthropogenic CO2[END_REF][START_REF] Stocker | Influence of CO2 emission rates on the stability of the thermohaline circulation[END_REF]. While the consequences of any changes are clearly severe, the probability for instabilities in the rate or pathways of NADW ventilation remains highly uncertain. Simple and complex models both suggest large changes are possible, but also that a strong overturning like that found in the modern ocean may be more difficult to disrupt than an overall weaker circulation [START_REF] Stocker | Influence of CO2 emission rates on the stability of the thermohaline circulation[END_REF][START_REF] Stommel | Thermohaline convection with two stable regimes of flow[END_REF][START_REF] Stocker | Rapid transitions of the ocean's deep circulation induced by changes in surface water fluxes[END_REF]. Likewise, most models simulate moderate to no reduction in AMOC in response to future source region buoyancy increases [START_REF] Hodell | Climate Change[END_REF] but may be biased towards stability [START_REF] Hofmann | On the stability of the Atlantic meridional overturning circulation[END_REF] and struggle to reproduce the rich spectrum of variability revealed by a decade of observations [START_REF] Srokosz | Observing the Atlantic Meridional Overturning Circulation yields a decade of inevitable surprises[END_REF][START_REF] Lozier | A sea change in our view of overturning in the subpolar North Atlantic[END_REF]. Testing these physical and conceptual models, and more generally the stability of NADW ventilation in warm climates, requires empirical constraints from beyond the current state of circulation.

Given a background climate similar to today, the modern mode of deep Atlantic ventilation with strong NADW influence (Fig. 1) appears stable on long multi-millennial timescales. Proxy reconstructions indicate that modern NADW ventilation pathways persisted with little multi-millennial variability in recent interglacial periods [START_REF] Galaasen | Rapid Reductions in North Atlantic Deep Water During the Peak of the Last Interglacial Period[END_REF][START_REF] Adkins | Variability of the North Atlantic thermohaline circulation during the last interglacial period[END_REF][START_REF] Kleiven | Reduced North Atlantic Deep Water coeval with the glacial Lake Agassiz freshwater outburst[END_REF][START_REF] Mcmanus | A 0.5-million-year record of millennial-scale climate variability in the North Atlantic[END_REF][START_REF] Oppo | Abrupt climate events 500,000 to 340,000 years ago: Evidence from subpolar North Atlantic sediments[END_REF]. By contrast, pronounced AMOC variability has occurred on timescales of a decade or less in observations [START_REF] Srokosz | Observing the Atlantic Meridional Overturning Circulation yields a decade of inevitable surprises[END_REF][START_REF] Lozier | A sea change in our view of overturning in the subpolar North Atlantic[END_REF] suggesting strong mean overturning is comprised of significant variance. However, little is known about NADW variability on the intermediary timescales, leaving the variability within a long-term 'vigorous' mean ventilation state poorly defined. There are few proxy reconstructions depicting higher-frequency variability and those available are largely confined to the last two interglacials, the Holocene and Marine Isotope Stage (MIS) 5e. During these periods, the largest changes in deep Atlantic ventilation involving reductions of NADW influence occurred on relatively short centennial timescales and were focused around intervals with wasting of continental ice sheet remnants from the preceding glaciation [START_REF] Galaasen | Rapid Reductions in North Atlantic Deep Water During the Peak of the Last Interglacial Period[END_REF][START_REF] Kleiven | Reduced North Atlantic Deep Water coeval with the glacial Lake Agassiz freshwater outburst[END_REF][START_REF] Hodell | Surface and deep-water hydrography on Gardar Drift (Iceland Basin) during the last interglacial period[END_REF]. This includes the century-long NADW reduction at 8.2 thousand years (ky) before present (B.P.) following the freshwater outburst flood from glacial Lake Agassiz [START_REF] Kleiven | Reduced North Atlantic Deep Water coeval with the glacial Lake Agassiz freshwater outburst[END_REF]. Absence of similarly large changes in the last ~eight ky of the Holocene (e.g., 12) has supported the notion of vigorous and stable ventilation as generally representative of interglacial boundary conditions.

Beyond the last two interglacials, little is known about centennial-scale variability in NADW, despite its relevance for delimiting the natural variability of ocean ventilation and the frequency of large NADW reductions under different background climates. The most recent interglacials MIS 5e, 7e, 9e, and 11c are particularly relevant, as these periods had similar climate boundary conditions to the current MIS 1 in addition to episodes of high-latitude warming, Greenland Ice Sheet (GrIS) retreat, and sea level rise relative to the modern [START_REF] Vernal | Natural variability of Greenland climate, vegetation, and ice volume during the past million years[END_REF][START_REF] Berger | Interglacials of the last 800,000 years[END_REF][START_REF] Reyes | South Greenland ice-sheet collapse during Marine Isotope Stage 11[END_REF].

This offers an opportunity to test the robustness of NADW ventilation under source region conditions similar to those projected for the future [START_REF] Hodell | Climate Change[END_REF]. Here we reconstruct northwest Atlantic bottom water δ 13 C to trace NADW influence (Fig. 1) over MIS 7e, 9e, and 11c and provide a detailed perspective on NADW ventilation instability during recent interglacials.

Our epibenthic foraminifera Cibicidoides wuellerstorfi (sensu stricto) δ 13 C record (19) from International Ocean Drilling Program (IODP) Site U1305 (57°29'N, 48°32'W; 3459 m water depth) at the Eirik Drift is situated to monitor lower NADW entering the deep Atlantic (Fig. 1). Due to the potential for uncertainty in δ 13 C reconstructions (e.g., 20), we only consider changes in the running mean of three samples (averaging five data points; see (19)) and signals outside the standard error of data within this window to reflect bottom water δ 13 C variability, which, given negligible influence from organic carbon fluxes [START_REF] Corliss | Influence of seasonal primary productivity on δ13C of North Atlantic deep-sea benthic foraminifera[END_REF], provides a proxy for past changes in the ventilation and distribution of water masses (e.g., [START_REF] Schmittner | Calibration of the carbon isotope composition (δ13C) of benthic foraminifera[END_REF][START_REF] Eide | A global ocean climatology of preindustrial and modern ocean δ13[END_REF]. The Eirik Drift bottom water δ 13 C record indicates large changes in deep Atlantic carbon chemistry during the interglacial δ 18 O plateaus of MIS 7e, 9e, and 11c (Fig. 2). Each interglacial contained abrupt changes in bottom water δ 13 C as large (≤1.0‰) as those of the bordering glacial terminations and inceptions (Fig. 2), and similar to those occurring after freshwater outburst floods such as the ~8.2 ky B.P. event [START_REF] Kleiven | Reduced North Atlantic Deep Water coeval with the glacial Lake Agassiz freshwater outburst[END_REF] and during MIS 5e [START_REF] Galaasen | Rapid Reductions in North Atlantic Deep Water During the Peak of the Last Interglacial Period[END_REF]. Absolute values range from near modern NADW levels (≥0.8‰; Fig. 1) to those typical of the glacial deep Atlantic [START_REF] Mcmanus | A 0.5-million-year record of millennial-scale climate variability in the North Atlantic[END_REF][START_REF] Oppo | Abrupt climate events 500,000 to 340,000 years ago: Evidence from subpolar North Atlantic sediments[END_REF][START_REF] Lisiecki | Atlantic overturning responses to Late Pleistocene climate forcings[END_REF][START_REF] Henry | North Atlantic ocean circulation and abrupt climate change during the last glaciation[END_REF]. We further use a transient interglacial (115-125 ky) simulation ( 19) with the isotopeenabled intermediate complexity iLOVECLIM Earth system model [START_REF] Bouttes | Including an ocean carbon cycle model into iLOVECLIM (v1. 0)[END_REF] to assess potential links between variability in deep Atlantic δ 13 C, NADW distribution, and AMOC. Simulated centennial-scale episodes of NADW shoaling and SSW expansion produce δ 13 C reductions in the deep Atlantic that strongly match the magnitude, rate and duration of the variability observed in our (Fig. 4) and other reconstructions (e.g., 10, 12-14) consistent with the inference that the δ 13 C variability reflect changes in NADW distribution. These large deep Atlantic δ 13 C changes, which are similar in magnitude to glacial millennial-scale changes (e.g., [START_REF] Lisiecki | Atlantic overturning responses to Late Pleistocene climate forcings[END_REF][START_REF] Henry | North Atlantic ocean circulation and abrupt climate change during the last glaciation[END_REF], were achieved without a total collapse of but with a significant (~16 to ~8 Sv) decrease in AMOC strength and accompanied by cooling in the subpolar North Atlantic (Fig. 4).

Our results call for reconsideration of the long-held notion of warm-climate stability in deep Atlantic carbon chemistry and ventilation. This view of stability likely remains true for the (multi-)millennial mean state, as previously depicted by lower-resolution records lacking the fidelity to resolve the shorter timescale that is characteristic of NADW reductions (Fig. 3).

High-resolution records are naturally biased towards the youngest strata and the current interglacial, the Holocene. Yet, when contextualized against the late Pleistocene interglacials, the Holocene stands out as having had the most stable lower NADW ventilation of the last half million years (Fig. 3), which was only strongly curtailed at ~8.2 ky B.P. [START_REF] Kleiven | Reduced North Atlantic Deep Water coeval with the glacial Lake Agassiz freshwater outburst[END_REF]. Bottom water δ 13 C and NADW reductions similar to that at ~8.2 ky B.P. were prevalent features of prior interglacials, occasionally even lasting millennia (Fig. 3). Ventilation patterns changed repeatedly from one similar to the modern (Fig. 1) to one with reduced NADW and incursions of SSW in the deep North Atlantic (~3.4 km), similar to that illustrated by our model simulation (Fig. 4).

The short duration of interglacial NADW reductions may indicate a change in the intrinsic ocean dynamics operating under different background climate states. The interglacial deep Atlantic is clearly better ventilated than the glacial on long equilibrium timescales [START_REF] Adkins | Variability of the North Atlantic thermohaline circulation during the last interglacial period[END_REF][START_REF] Mcmanus | A 0.5-million-year record of millennial-scale climate variability in the North Atlantic[END_REF][START_REF] Oppo | Abrupt climate events 500,000 to 340,000 years ago: Evidence from subpolar North Atlantic sediments[END_REF][START_REF] Lisiecki | Atlantic overturning responses to Late Pleistocene climate forcings[END_REF][START_REF] Henry | North Atlantic ocean circulation and abrupt climate change during the last glaciation[END_REF]. However, the magnitude of ventilation pattern changes that are possible appears similar in (de-)glacial (e.g., [START_REF] Adkins | Variability of the North Atlantic thermohaline circulation during the last interglacial period[END_REF][START_REF] Henry | North Atlantic ocean circulation and abrupt climate change during the last glaciation[END_REF][START_REF] Mcmanus | Collapse and rapid resumption of Atlantic meridional circulation linked to deglacial climate changes[END_REF] and interglacial periods when variability in lower NADW is considered at shorter timescales (Fig. 3). The centennial-scale duration and transient nature of most interglacial NADW reductions (Fig. 3; Fig. 4) suggests the modern ventilation pattern tends to recover quickly when perturbed, and is similar to the AMOC recovery timescale seen in many numerical models forced with buoyancy increases (e.g., 4). With this in mind, the longer-lasting NADW reductions in MIS 7e (~233.5-234.5 ky), 9e (~323-326 ky), and late 11c (~401-408 ky) either required more sustained forcing or suggests that the recovery timescale following perturbations is not fixed. Most interglacial NADW reductions were still short-lived compared to those associated with glacial (Dansgaard-Oeschger) variability (e.g., [START_REF] Henry | North Atlantic ocean circulation and abrupt climate change during the last glaciation[END_REF], suggesting either that NADW ventilation behaved differently or the persistence of any forcing changed, depending on the climate state. One possible explanation for this timescale difference is the extensive glacial expansion of high-latitude sea ice, which could promote a baseline increase in SSW ventilation [START_REF] Ferrari | Antarctic sea ice control on ocean circulation in present and glacial climates[END_REF] and prolong the duration of northern ventilation anomalies [START_REF] Vettoretti | Fast Physics and Slow Physics in the Nonlinear Dansgaard-Oeschger Relaxation Oscillation[END_REF]. A lack of strong sea ice responses could also explain the potentially muted climate variability in interglacial compared to glacial climates [START_REF] Galaasen | Rapid Reductions in North Atlantic Deep Water During the Peak of the Last Interglacial Period[END_REF][START_REF] Mcmanus | A 0.5-million-year record of millennial-scale climate variability in the North Atlantic[END_REF][START_REF] Oppo | Abrupt climate events 500,000 to 340,000 years ago: Evidence from subpolar North Atlantic sediments[END_REF][START_REF] Martrat | Four climate cycles of recurring deep and surface water destabilizations on the Iberian margin[END_REF], despite the presence of NADW variability. More high-resolution climate records spanning past interglacials are however needed to conclusively evaluate the impacts of warm-climate NADW reductions, and delineate its role relative feedbacks such as sea ice responses in driving interglacial climate variability.

Model simulations suggest future warming and freshwater addition from an intensified hydrological cycle and ice sheet melting could all increase source region buoyancy and curtail convective NADW renewal [START_REF] Hodell | Climate Change[END_REF][START_REF] Stocker | Influence of CO2 emission rates on the stability of the thermohaline circulation[END_REF]. The common occurrence of NADW reductions in past interglacials (Fig. 3) clearly demonstrates the potential for large changes in deep Atlantic ventilation and allows us to explore the triggers for perturbations. NADW reductions during the last two interglacial periods were confined mainly to the early warm interglacial phases, concurrent with high northern hemisphere summer insolation and known freshwater outburst floods accompanying the final retreat of residual glacial ice sheets (Fig. 3) [START_REF] Galaasen | Rapid Reductions in North Atlantic Deep Water During the Peak of the Last Interglacial Period[END_REF][START_REF] Kleiven | Reduced North Atlantic Deep Water coeval with the glacial Lake Agassiz freshwater outburst[END_REF][START_REF] Ellison | Surface and deep ocean interactions during the cold climate event 8200 years ago[END_REF][START_REF] Nicholl | A Laurentide outburst flooding event during the last interglacial period[END_REF].

While stratigraphically belonging to interglacial periods, to the extent that these anomalies are related to wasting vestiges of glaciation, they are likely best viewed mechanistically as the final episodes of deglaciation. By contrast, NADW reductions in MIS 7e, 9e, and 11c occurred in the mid and late interglacial phases under low summer insolation (Fig. 3) and after any likely deglacial freshwater influences. This implies that NADW reductions can occur without the excess buoyancy input provided by wasting residual glacial ice sheets, or the influence of large continental ice sheets on atmospheric circulation. Ice sheet activity may still have played a role in regulating the stability of NADW ventilation during some periods. NADW reductions in MIS 7e, 9e, and 11c often coincided with, or were preceded by, input of ice-rafted debris (IRD) at Site U1305 (Fig. 3), indicating supply of icebergs and freshwater proximal to the NADW source region. Furthermore, the prolonged NADW reduction of MIS 9e was associated with elevated southern GrIS sediment discharge and MIS 1 stability with low GrIS activity (32), while particularly strong GrIS retreat in MIS 11c [START_REF] Vernal | Natural variability of Greenland climate, vegetation, and ice volume during the past million years[END_REF][START_REF] Reyes | South Greenland ice-sheet collapse during Marine Isotope Stage 11[END_REF] occurred alongside persistent NADW variability (Fig. 3). These observations are consistent with ice sheet activity and freshwater addition intermittently influencing the formation or downstream density of lower NADW. However, variability in NADW ventilation, IRD input, and GrIS discharge also occurred independently of each other (Fig. 3), implicating additional controls on NADW ventilation and supporting models suggesting convective instability is possible with relatively small buoyancy input if delivered to the convection regions (e.g., 6).

Our results suggest we should consider rapid and large changes in NADW ventilation not only as a possibility [START_REF] Galaasen | Rapid Reductions in North Atlantic Deep Water During the Peak of the Last Interglacial Period[END_REF][START_REF] Kleiven | Reduced North Atlantic Deep Water coeval with the glacial Lake Agassiz freshwater outburst[END_REF][START_REF] Ellison | Surface and deep ocean interactions during the cold climate event 8200 years ago[END_REF] but even as an intrinsic feature of centennial-scale variability in warm climate states. This has implications for constraining the potential for and cause of change in the modern Atlantic. First, it supports concerns that disregarding large variability in simulations may have biased future AMOC projections towards stability [START_REF] Hofmann | On the stability of the Atlantic meridional overturning circulation[END_REF]. The possibility of large natural variability on decadal [START_REF] Srokosz | Observing the Atlantic Meridional Overturning Circulation yields a decade of inevitable surprises[END_REF][START_REF] Lozier | A sea change in our view of overturning in the subpolar North Atlantic[END_REF] to centennial timescales (Fig. 3) also complicates attribution of variability in the deep Atlantic, but the characteristics of this variability may be used to differentiate between natural and anthropogenic change in the coming century. While past changes were predominantly multi-centennial, there are also climate and ocean conditions that can drive longer NADW reductions, as evidenced for example by the ~3000 years-long anomaly in mid-MIS 9e (Fig. 3). Specifically what these conditions were remains unclear, but the triggers for NADW instability have clearly operated across the range of recent interglacial climate conditions. Recognizing this requires moving beyond the notion of vigorous and stable deep Atlantic ventilation as representative of warm climate states (e.g., [START_REF] Hodell | Climate Change[END_REF][START_REF] Stommel | Thermohaline convection with two stable regimes of flow[END_REF][START_REF] Adkins | Variability of the North Atlantic thermohaline circulation during the last interglacial period[END_REF], and towards conceptual and numerical models that can account for pronounced variability across various timescales and climate states.

Materials and Methods

Sample processing and C. wuellerstorfi stable isotopes

The International Ocean Drilling Program (IODP) Site U1305 intervals spanning Marine Isotope Stage (MIS) 7e, 9e, and 11c were identified from the stable isotope stratigraphy of Hillaire-Marcel et al. [START_REF] Hillaire-Marcel | Foraminifer isotope study of the Pleistocene Labrador Sea, northwest North Atlantic (IODP Sites 1302/03 and 1305), with emphasis on paleoceanographical differences between its "inner" and "outer" basins[END_REF] and continuously subsampled at 2-cm spacing at the IODP Bremen Core Repository by the curatorial staff. Bulk sediment samples were kept in deionized water on a shaker for 12 hours for disaggregation before being wet sieved using a 63 µm mesh sieve to separate the fine (<63 µm) and coarse fraction material (>63 µm). Following wet sieving, samples were dried at 50°C. for δ 18 O and δ 13 C, respectively. The standard cleaning step for foraminiferal stable isotope analyses, involving methanol and ultrasonication, was avoided to retain mass and increase the number of measurements possible, as C. wuellerstorfi shells were generally few and small in size, often providing a total weight at the lower limit possible for analysis (~10-15 µg). Instead of performing the standard cleaning step, C. wuellerstorfi tests were visually cleaned using a brush and deionized water ('brush-cleaned') to remove any foreign material on or within them prior to stable isotope analyses. To test the impact of omitting the standard cleaning step, we measured surplus shells from samples (n=22) containing sufficient mass that were cleaned following standard protocols ('standard-cleaned') for removing fine-grained material: adding methanol to reaction vials containing the shells and keeping them in an ultrasonic bath for 10 seconds before extracting the methanol. The δ 18 O and δ 13 C values of the 'standard-cleaned' C.

wuellerstorfi tests are very similar to the sample average of the 'brush-cleaned' tests (Fig. S1).

Further, the intra-sample reproducibility of the 'brush-cleaned' Site U1305 C. wuellerstorfi δ 18 O and δ 13 C data is similar to a large data set of exclusively 'standard-cleaned' C. wuellerstorfi shells from the same region (the MIS 5e section of MD03-2664; location shown in Fig. 1) (Fig. S1). This indicates that the visual 'brush-cleaning' and 'standard-cleaning' steps performed similarly well with no discernible difference in the stable isotope values or in the intra-sample variability. Nonetheless, performing the standard cleaning step for foraminiferal stable isotope analyses is highly recommended when possible. Indeed, in certain regions and time intervals, a cleaning procedure more stringent than the standard one is likely required to accurately determine foraminiferal stable isotope values (40).

Ice-rafted debris

The ice-rafted debris (IRD) counts were performed on the same Site U1305 samples as the benthic foraminiferal stable isotope measurements, but at lower sampling density. IRD counts were performed at 32-cm spacing for MIS 7e and 11c, and 16-cm spacing for MIS 9e (42). Following the sample processing steps outlined above, material in the >150 µm fraction were split, IRD grains visually identified, and IRD calculated as the percent of ≥300 counted entities.

Hole U1305C MIS 11c mcd fine-tuning

Sediment physical properties (e.g., magnetic susceptibility) indicated cm-scale offsets between Hole U1305C and the original (Hole A & B) splice over ~74.5-78.5 mcd (Fig. S3A), corresponding to most of MIS 11c. We fine-tuned the mcd scale for the Hole U1305C MIS 11c interval using magnetic susceptibility, shifting it between -3 cm and -19 cm to align the physical property records (Fig. S3B), and include both the original and corrected core depth scales in the MIS 11c data table. values, we ran to additional 125-115 ky transient simulations using similar initial conditions as above but starting at 125 ky with either: 1) atmospheric δ 13 C fixed at value decreased by ~0.6‰ (at -7.1‰); or 2) with 50% decreased primary productivity in the modeled convection regions off southern Greenland and in the Nordic Seas. We expand on the model results in the discussion section of the supplement below.

The modeled inorganic carbon cycle is represented by dissolved inorganic carbon (DIC) and alkalinity (ALK), while the organic carbon cycle includes phytoplankton, zooplankton, dissolved organic carbon (DOC), slow dissolved organic carbon (DOCs), particulate organic carbon (POC), and CaCO 3 . The phytoplankton is partially remineralized as it sinks through the water column, while all the POC and CaCO 3 is remineralized at depth. The remineralization profile follows an exponential law, adjusted to have less remineralization in the upper layers and more at depth. Carbon fractionation during photosynthesis fixes 12 C in the organic matter, which is added back by the remineralization process at depth. At the air-sea interface, the carbon flux is computed from the CO 2 partial pressure difference between the atmosphere and ocean at a constant gas exchange coefficient of 0.06 mol m -2 yr -1 , where sea surface pCO 2 is a function of temperature, salinity, DIC, and ALK following Millero (44). The modeled δ 13 C distribution is thus affected by air-sea gas exchange and organic matter production/remineralization, and transported by the advection-diffusion scheme of the model.

Unlike 12 C, the atmospheric 13 C is prognostically simulated in response to land and ocean processes.

Supplementary Text

Site U1305 age model

The age models for the MIS 7e, 9e, and 11c intervals of Site U1305 were constructed by correlating our benthic δ 18 O record to and adopting the age model constructed for ODP Site 983 [START_REF] Raymo | Stability of North Atlantic water masses in face of pronounced climate variability during the Pleistocene[END_REF]45) (Fig. 2, Fig. S4). Using Site 983 benthic δ 18 O as a tuning target has advantages over other reference records. First, Site U1305 shares well-defined structures in benthic δ 18 O with Site 983 (Fig. 2; Fig. S4) allowing relatively robust correlation. Further, high-resolution IRD records are available from both sites and allow us to validate the benthic δ 18 O-correlation near deglacial intervals where Site U1305 benthic δ 18 O data were often lacking due to C. wuellerstorfi absence (Fig. 2). Tie points were defined based on major benthic δ 18 O transitions and linearly interpolated between to obtain ages for all core depths (Fig. S4). The Site U1305 and MD03-2664 IRD records were also used to guide the determination of tie points at the start of the interglacial δ 18 O plateaus, as deglacial IRD increases are observed to coincide with transient decreases in benthic δ 18 O values before MIS 5e [START_REF] Galaasen | Rapid Reductions in North Atlantic Deep Water During the Peak of the Last Interglacial Period[END_REF] and MIS 7e, 9e, and 11c in this region (Fig. 2; Fig. S4). Given their transient nature and co-occurrence with deglacial IRD, these deglacial benthic δ 18 O reductions may reflect contamination from low-δ 18 O detrital carbonate commonly deposited during Heinrich-events (40). Consequently, deglacial samples with low C. wuellerstorfi δ 18 O and high IRD were disregarded when we constructed our age model and we defined the first interglacial benthic δ 18 O value, and start of the interglacial plateaus, as the first continuously low δ 18 O value occurring after large deglacial IRD increases.

Corroborating this approach, the Site U1305/MD03-2664 and Site 983 deglacial IRD peaks align using this additional constraint for the benthic δ 18 O tuning (Fig. S4) but would otherwise be offset by a few thousand years.

Age uncertainties involved with benthic δ 18 O-based age models can be considerable.

For example, age differences between major δ 18 O transitions can reach up to a few thousand years between different ocean basins (46). The absolute age uncertainty provided by the Site U1305-Site 983 correlation is likely less than millennial, given i) the relative proximity of these core sites, ii) the similarity of the benthic δ 18 O records-indicating a shared δ 18 O evolution, and iii) the alignment of deglacial IRD peaks (Fig. 2; Fig. S4). Despite a relatively robust correlation, the original age model still carries considerable uncertainty in absolute ages.

For example, adopting a different age model constructed for Site 983 by Barker et al. ( 45) (e.g., EDC3 and AICC2012), would shift absolute ages up to a few thousand years.

In addition to absolute age uncertainties, relative (sample-to-sample) age uncertainties likely also exist for the Eirik Drift (Site U1305 and MD03-2664) records. Correlation of benthic δ 18 O records is achieved using a limited number of tie points. This is especially true for interglacial δ 18 O plateaus, here defined by two (MIS 5e, 7e, and 9e) or three (MIS 11c) tie points (Fig. 2, Fig. S4). These tie points were linearly interpolated between, assuming constant sedimentation rates. However, interglacial sedimentation rates can vary on a range of timescales in this area [START_REF] Hillaire-Marcel | Foraminifer isotope study of the Pleistocene Labrador Sea, northwest North Atlantic (IODP Sites 1302/03 and 1305), with emphasis on paleoceanographical differences between its "inner" and "outer" basins[END_REF](47)(48)(49). For example, radiocarbon-dated sections indicate that sedimentation rates were higher in the early compared to the late phase of the current interglacial at multiple Eirik Drift locations (e.g., [START_REF] Kleiven | Reduced North Atlantic Deep Water coeval with the glacial Lake Agassiz freshwater outburst[END_REF]48). If this temporal sedimentation pattern persisted in the older interglacial periods, our age models based on the conservative approach of linear interpolation may over-and underestimate the durations of early and late interglacial intervals, respectively.

Site U1304, MD03-2664, MD03-2665, and MD99-2227 age models

To place all proxy records on a common age scale, we revised the age models for all core intervals containing data we compare to the IODP Site U1305 records. We tuned the MIS 5e interval of MD03-2664 and the MIS 5e, 7e, 9e, and 11c intervals of IODP Site U1304 to the same ODP Site 983 reference record (Fig. S4), applying identical tie points and linearly interpolating between as outlined above. The MIS 1 and last deglacial intervals of Site U1304 and MD03-2665 were left on their original age models as presented in Xuan et al. [START_REF] Xuan | Quaternary magnetic and oxygen isotope stratigraphy in diatom-rich sediments of the southern Gardar Drift (IODP Site U1304, North Atlantic)[END_REF] and Kleiven et al. [START_REF] Kleiven | Reduced North Atlantic Deep Water coeval with the glacial Lake Agassiz freshwater outburst[END_REF], respectively.

To compare our records to the MD99-2227 southern Greenland sediment discharge reconstructions [START_REF] Reyes | South Greenland ice-sheet collapse during Marine Isotope Stage 11[END_REF][START_REF] Hatfield | Interglacial responses of the southern Greenland ice sheet over the last 430,000 years determined using particle-size specific magnetic and isotopic tracers[END_REF][START_REF] Colville | Sr-Nd-Pb isotope evidence for ice-sheet presence on southern Greenland during the last interglacial[END_REF], we transferred our benthic δ 18 O-based age models for MIS 5e (MD03-2664), 7e, 9e, and 11c (Site U1305) to MD99-2227 by correlating magnetic susceptibility between this core and MD03-2664 (MIS 5e) and Site U1305 (MIS 7e, 9e, and 11c) (Fig. S5). The strong similarity of the magnetic susceptibility records indicate that these core sites shared sedimentation histories, providing robust relative age control to comparisons 4A), proxy records suggest considerable changes in preformed water mass δ 13 C occurred between and even within past interglacial periods. For example, the preformed δ 13 C of northern source water may have been higher in the Holocene than the late Pleistocene interglacials (see data composites, purple lines, in Fig. 3), while it likely increased across MIS 5e (Fig. 3) consistent with planktic and epibenthic foraminifera δ 13 C records from the Nordic Seas (e.g., 53). This model-data difference should be noted when comparing and contrasting the simulated and reconstructed time series. For the proxy reconstructions, we took this into account by averaging multiple events in order to illustrate common features independent of specific interglacials and preformed δ 13 C values (Fig. 4). Still, the similarity of the modeled and reconstructed bottom water

  While similar in magnitude, the frequency, timing, and duration of these changes differ among individual interglacial periods. Low bottom water δ 13 C values persist for a millennium or more during late MIS 7e (~233.5-243.5 ky; on our age model, 19) and mid-to late MIS 9e (~323.0-326.0 ky), whereas large (~0.5‰) multi-centennial variability punctuated MIS 11c superimposed on multi-millennial trends. Low bottom water δ 13 C values at Site U1305 likely reflect reduced NADW influence and changes in deep Atlantic ventilation patterns. Reduced (high-δ 13 C) NADW influence and incursions of (low-δ 13 C) Southern source water (SSW) explain many features of the observed variability, including the: 1) spatial consistency of intermittently low δ 13 C observed at different deep sites (Site U1304 and U1305; Fig. 3); 2) abruptness of the δ 13 C changes as the NADW-SSW water mass boundary shifted across the core sites; 3) shift of Site U1305 δ 13 C towards the millennially averaged values found near the northern or the southern source regions (Fig. 3); and 4) association of high (low) C. wuellerstorfi δ 13 C with high (low) C. wuellerstorfi B/Ca in selected Eirik Drift samples (Fig. S6) (19).

from

  Northeast Water Polynya plankton tows, sediment traps, and surface sediments. Paleoceanography 11, 679-699 (1996).

Fig. 1 .

 1 Fig. 1. Core locations. IODP Site U1305 (57°29'N, 48°32'W; 3459 m), MD03-2664

Fig. 2 .

 2 Fig. 2. IODP Site U1305 MIS 7e, 9e, and 11c C. wuellerstorfi stable isotope and ice-rafting records. A) Benthic δ 18 O from IODP Site U1305 (thin blue line, sample average of replicate

Fig. 3 .

 3 Fig. 3. Variability in NADW ventilation during interglacials MIS 1-11c. Focused on the interglacial δ 18 O plateaus: A) 65°N insolation at June 21 st (orange) (35); B) core MD99-2227

Fig. 4 .

 4 Fig. 4. Modeled and reconstructed deep Atlantic δ 13 C changes. (Left) The iLOVECLIM simulated δ 13 C distribution and Eirik Drift core location (red circle) along a north-south

  Epibenthic foraminifera Cibicidoides wuellerstorfi (sensu stricto) shells (Fig.S1) were selectively picked from the >150 µm sediment fraction for stable isotope analyses. The stable isotope analyses were performed at the Facility for advanced isotopic research and monitoring of weather, climate and biogeochemical cycling (FARLAB), Department of Earth Science, University of Bergen, Norway, on a Finnigan MAT 253 mass spectrometer coupled to an automated Kiel IV preparation line kept at constant 70°C. Measurements were performed on one to three individual C. wuellerstorfi shells, depending on availability and size, and duplicated per sample when possible (~64% of the samples analyzed). Fig.S2shows the C. wuellerstorfi stable isotope results from Site U1305 including all individual data points. We used Carrera Marble (CM12) as a working standard measured parallel to the foraminifera samples, and all values are reported relative Vienna Pee Dee Belemnite (VPDB) calibrated using National Bureau of Standards (NBS) standard NBS 19 and NBS 18. The long-term reproducibility (1σ) of in-house standards over the analysis period was ≤0.08‰ and ≤0.04‰

C

  . wuellerstorfi B/Ca We measured B/Ca ratios in C. wuellerstorfi tests from intervals in MD03-2664 (MIS 5e) and Site U1305 (MIS 7e, 9e, and 11c). To obtain sufficient mass for B/Ca analyses (~250-300 µg), and due to scarcity of C. wuellerstorfi tests in these sediments, it was often necessary to combine tests from up to a maximum of six adjacent samples (or 12 cm of core) that we restricted according to consistent C. wuellerstorfi δ 13 C values. Following the selective picking, the C. wuellerstorfi tests were opened using clean glass slides and transferred into acid-leached vials. The C. wuellerstorfi tests were subsequently cleaned for contaminating phases, including clay removal, reductive and oxidative steps, a weak acid leach, and dissolution in dilute HNO 3 . The B/Ca analyses were performed using the method outlined in Rosenthal et al. (41) on a Finnigan MAT Element XR Sector Field Inductively Coupled Plasma Mass Spectrometer (ICP-MS) at the ICP-MS laboratory at the Institute of Marine and Coastal Sciences, Rutgers, The State University of New Jersey, USA. Site U1305 N. pachyderma (s) δ 13 C records N. pachyderma (s) tests were selectively picked from the 150-212 µm sediment fraction at continuous 4-cm spacing across MIS 7e, 9e, and 11c (with notable gaps only in MIS 9e due to scarcity of N. pachyderma (s) tests). Prior to the stable isotope analyses, the tests were cleaned by adding methanol to the sample kept in reaction vials and ultrasonicating them for ten seconds before removing the supernatant. The stable isotope analyses were performed at FARLAB, University of Bergen, Norway, as outlined for benthic foraminifera C. wuellerstorfi above, and with identical standard reproducibility. Measurements on N. pachyderma (s) were replicated whenever possible (~92% of the samples), and each individual measurement was performed on 6-10 individual N. pachyderma (s) tests.

  iLOVECLIM model simulation We used the iLOVECLIM Earth system model of intermediate complexity to simulate and assess potential relationships between NADW distribution, northwest Atlantic bottom water δ 13 C, and AMOC. The iLOVECLIM model is an isotope-enabled development branch of LOVECLIM version 1.2 (43) and includes an ocean component (CLIO) with 20 vertical layers and 3° by 3° horizontal resolution as well as land and ocean carbon cycle modules (25). We performed a transient simulation for 125-115 ky (corresponding to MIS 5e) using annually interpolated greenhouse gas and orbital forcings following the third phase of the Paleoclimate Modelling Intercomparison Project (PMIP; https://pmip3.lsce.ipsl.fr/), initialized from a quasiequilibrium spin up of 5000 years forced with constant 125 ky boundary conditions. We performed two quasi-equilibrium spin-ups prior to the 125-115 ky simulation, each integrated for 5000 years. The first is based on preindustrial conditions and the second on 125 ky boundary conditions. The preindustrial simulation was used to validate and confirm that iLOVECLIM reproduces the spatial distribution of preindustrial ocean δ 13 C (22), and the 125 ky spin-up to initialize the 125-115 transient simulation. To test the sensitivity of the simulated deep Atlantic δ 13 C variability to changes in surface biological processes and preformed δ 13 C

  of the Site U1305/MD03-2664 and MD99-2227 records. The alignment of glacial and interglacials values in epibenthic foraminifera C. wuellerstorfi δ 18 O from MD03-2664/Site U1305 and planktic foraminifera N. pachyderma (s) δ 18 O from MD99-2227 on the obtained age models supports the magnetic susceptibility correlation (Fig. S5). Carlson et al.'s (50) age model was used for the MIS 1 and last deglacial intervals of MD99-2227. C. wuellerstorfi B/Ca and N. pachyderma (s) δ 13 C To test the fidelity of Eirik Drift C. wuellerstorfi δ 13 C as recorder of bottom water carbon chemistry and water mass tracer, we use epibenthic foraminifera C. wuellerstorfi B/Ca, a proxy for bottom water carbonate ion saturation (Δ[CO 3 2-]) and independent metric of the influence of high-[CO 3 2-] NADW versus low-[CO 3 2-] SSW (e.g., 51, 52). The Eirik Drift C. wuellerstorfi B/Ca data span a range of 185-215 µmol/mol with distinct changes within each of MIS 5e, 7e, 9e, and 11c (Fig. S6). Using Yu & Elderfield's (51) B/Ca to Δ[CO 3 2-] relationship for C. wuellerstorfi, the B/Ca data indicates intra-interglacial changes in bottom water [CO 3 2-] by 25-30 µmol/kg, similar to that expected from shifting between NADW and SSW influence (e.g., 52). The concurrent and coupled changes in C. wuellerstorfi δ 13 C by up to ~0.8‰ within the same sample pool (Fig. S6) is similarly consistent with shifts between NADW and SSW influence. Indeed, the paired change in Eirik Drift C. wuellerstorfi B/Ca and δ 13 C we observe is similar to that recorded in the last glacial to Holocene sections of two North Atlantic cores from similar water depths (Fig. S6) that was previously interpreted to reflect the wellestablished deglacial shift from SSW to NADW influence in the deep North Atlantic (52). We suggest that the Eirik Drift C. wuellerstorfi B/Ca record supports the interpretation of the C. wuellerstorfi δ 13 C variability as reflecting changes in NADW versus SSW influence. An alternative explanation for co-variability in trace element ratios and δ 13 C is contamination by secondary CaCO 3 precipitation and presence of authigenic overgrowths. However, several lines of evidence argue against a role for contamination by secondary CaCO 3 precipitation. The C. wuellerstorfi δ 13 C record suggest no discernible influence by its own merit, showing for example i) absolute values within the range of values observed in the modern ocean or relevant reconstructions (e.g., 22, 23; Fig. 3) and ii) a consistency in the signal and a lack of extremely large fluctuations that would require the mass and isotope value of any contamination to have adjusted itself to balance changes in the mass and isotope value of the foraminifera tests. Further, secondary CaCO 3 precipitation should, if present, influence all foraminifera tests in a given core depth to some degree similarly. That is, if secondary CaCO 3 precipitation drove relatively large changes in C. wuellerstorfi δ 13 C and B/Ca, planktic foraminifera δ 13 C should also show low δ 13 C values in addition to some degree of covariability. The Eirik Drift N. pachyderma (s) δ 13 C records from MIS 5e, 7e, 9e, and 11c do not show values as low as the C. wuellerstorfi δ 13 C records and there is no significant relationship between N. pachyderma (s) and C. wuellerstorfi δ 13 C (Fig. S6). In sum, we suggest secondary CaCO 3 precipitation is either unimportant or entirely absent and consider it as an unlikely explanation for the co-variability in C. wuellerstorfi B/Ca and δ 13 C. Changes in the influence of NADW versus SSW could conversely explain these observations. iLOVECLIM model simulation results We identified persistent centennial-scale variability in δ 13 C and the distribution of NADW in our 125-115 ky transient simulation occurring over a ~6 ky long interval in-between intervals of relative stability during the initial 2-3 and final 1-2 ky. A subsequent study will outline and discuss the results of the model simulation in detail (Kessler et al., in prep.). Here, we use the simulation to help the interpretation of the reconstructed bottom water δ 13 Cvariability by assessing how centennial-scale changes in NADW distribution can impact deep Atlantic δ 13 C. We focus on the simulated episodes of NADW shoaling and recovery to assess the possible rate, duration, and magnitude of water mass distribution and associated bottom water δ 13 C changes and compare these to the reconstructed δ 13 C variability. In the simulation, episodes of NADW shoaling were abruptly initiated, lasted some centuries, and were associated with decreases in AMOC strength and North Atlantic bottom water δ 13 C (Fig.4).The simulated episodes of NADW shoaling/deepening produced distinct North Atlantic bottom water δ 13 C variability as (low-δ 13 C) Southern source water (SSW) expanded/contracted in concert with (high-δ 13 C) NADW contracting/expanding (Fig.4). In the northwest Atlantic region corresponding to the location of Eirik Drift core sites U1305, MD03-2664, and MD03-2665, shoaling of NADW and incursions of SSW manifested as ~0.4‰ decreases in bottom water δ13 C achieved over a few decades, events that ended equally abrupt as NADW deepened some centuries (~100-500 years) later (Fig.4). To illustrate these NADW and δ 13 C distribution changes, Fig. 4A and Fig. 4B displays the North Atlantic mean δ 13 C distribution below 500 m water depth for all simulated years with anomalously strong AMOC (>2σ; n=460 model years; mean AMOC strength: 16.75±0.70 Sv) and anomalously weak AMOC (<2σ; n=63 model years; mean AMOC strength: 8.00±0.42 Sv), respectively. We further selected two simulated intervals of NADW shoaling and recovering to compare to the reconstructed bottom water δ 13 C variability (Fig. 4B), differing from other simulated shoaling episodes only in duration. The magnitude of any bottom water δ 13 C variability driven by such changes in the relative influence of northern versus southern source water could depend strongly on the preformed δ 13 C of, and the gradient between, competing water masses. While iLOVECLIM captures the preindustrial preformed δ 13 C of northern and southern source waters relatively well (e.g., compare Fig. 1 to Fig.

  δ13 C changes could result from the simulation having a specific set of preformed δ 13 C values in NADW and SSW. To assess how different background states and preformed δ 13 C in NADW and SSW could impact the magnitude and character of the simulated bottom water δ 13 C variability, we reran the simulation twice with 1) atmospheric δ 13 C lowered by ~0.6‰ and 2) 50% decreased primary productivity in the simulated deep water formation regions off southern Greenland and in the Nordic Seas. Both experiments shift the absolute values of the simulated Eirik Drift bottom water δ 13 C time series but in both cases the relative magnitudes and character of the variability is only negligibly impacted (Fig.S7). The simulation with perturbed primary productivity, resulting in limited change in deep Atlantic δ 13 C, supports previous studies suggesting that organic carbon fluxes have little influence on the δ 13 C of C. wuellerstorfi (e.g., 21). Thus, the character of the Eirik Drift bottom water δ 13 C variability driven by shifts in the distribution of water masses appears to be relatively stable in face of different preformed δ 13 C values and background biological processes in the model, indicating that the similarity between the simulated and reconstructed variability is not strongly dependent on the specific model configuration. That this character is similar to the reconstructed Eirik Drift bottom water δ 13 C reductions, including in the magnitude, rate, and duration of events (Fig.4; Fig.S7), supports the inference that characteristic deep Atlantic δ 13 C changes can be explained with changes in the distribution and influence of NADW and SSW.

Fig. S1 .

 S1 Fig. S1. Eirik Drift C. wuellerstorfi (sensu stricto) replicate stable isotope values for different cleaning procedures. Cross-plots comparing sample replicate measurements where sufficient C. wuellerstorfi tests were present to allow it (n=493): A) Site U1305 C. wuellerstorfi δ 13 C values of tests visually cleaned using a brush and distilled water ('brushcleaned'; sample average; x-axis) plotted versus tests cleaned using the standard protocol involving methanol (yellow circles), Site U1305 C. wuellerstorfi δ 13 C values of brush-cleaned versus brush-cleaned tests (red circles), and MD03-2664 C. wuellerstorfi δ 13 C values of standard-cleaned vs. standard-cleaned tests (blue circles; 10). B) Same as in A) but for C. wuellerstorfi δ 18 O. Note the similarity in stable isotope values for the different cleaning protocols. C) Umbilical (lef), apertural (middle), and spiral view (right) of an example specimen of C. wuellerstorfi sensu stricto from Site U1305 Hole C, 8H-3, 100-102 cm.

Fig. S2 .

 S2 Fig. S2. Stable isotope results. The Site U1305 C. wuellerstorfi δ 13 C (red) and δ 18 O (blue) records for A) MIS 7e, B) MIS 9e, and C) MIS 11c plotted versus core depth (meter composite depth; mcd) showing all individual data points (dots; see inset at bottom for color coding), the sample average values (thin lines), the three-point running mean (bold lines), and the standard error of the mean of three-point running window (gray shading; includes on average five individual data points).

Fig. S3 .

 S3 Fig. S3. Hole U1305C MIS 11c meter depth scale tuning. A) Magnetic susceptibility (54) of the original splice (gray) and our sampled intervals of U1305A and U1305B (blue) and U1305C (red) on the original mcd scale. Note the offset between U1305C and the original splice over ~74.5-78.5 mcd. B) Magnetic susceptibility (54) of the original splice (gray) and our sampled intervals of U1305A and U1305B (blue) on mcd, and our sampled interval of U1305C (red) on the tuned meter depth scale.

Fig. S4 .

 S4 Fig. S4. Age model correlation. Plotted on LR04 age (ky): A) Benthic foraminifera δ 18 O from MD03-2664 (dark blue; MIS 5e) (10), Site U1305 (dark blue; MIS 7e, 9e, and 11c), Site U1304 (light blue) (15, 38), and age model correlation target Site 983 (gray) (33) with tie points denoted by triangles; B) Site 983 ice-rafted debris (IRD) as grains gram -1 (45) compared to NW Atlantic IRD (black with gray shading) from MD03-2664 (MIS 5e) (55) and Site U1305 as percent IRD grains in total entities >150 µm shown only for the deglacial intervals where it was used to guide the age model construction. Dashed line in A) marks a prolonged interval lacking Site 983 δ 18 O data across the MIS 10 to 9e transition. Yellow shading denotes the interglacial δ 18 O plateaus as in the main text.

Fig. S5 .

 S5 Fig. S5. Magnetic susceptibility correlation of MD03-2664/Site U1305 and MD99-2227. Plotted on the LR04 age scale obtained from transferring our δ 18 O-based age models of MD03-2664 (MIS 5e) and Site U1305 (MIS 7e, 9e, and 11c) to MD99-2227 by correlating magnetic susceptibility: A) epibenthic foraminifera C. wuellerstorfi δ 18 O from core MD03-2664 (dark blue; MIS 5e) (10) and Site U1305 (dark blue; MIS 7e, 9e, and 11c) and planktic N. pachyderma (s) from core MD99-2227 (light blue) (32); B) magnetic susceptibility for MD03-2664 (red; MIS 5e) (56), Site U1305 (red; MIS 7e, 9e, and 11c) (54), and MD99-2227 (black) (47).

Fig. S6 .

 S6 Fig. S6. Eirik Drift C. wuellerstorfi δ 13 C, C. wuellerstorfi B/Ca, and N. pachyderma (s) δ 13 C. A) The Eirik Drift C. wuellerstorfi B/Ca values (blue diamonds) of tests combined from selected intervals in MIS 5e (core MD03-2664), 7e, 9e, and 11c (Site U1305), the average δ 13 C values of C. wuellerstorfi in the same intervals (red diamonds), and the full C. wuellerstorfi δ 13 C time series shaded in the background for comparison (thin gray line, sample average; bold red line: 3-point running mean) over MIS 5e (MD03-2664; 10), 7e, 9e, and 11c (Site U1305) plotted versus age (ky; LR04). B) Cross-plot of the Eirik Drift C. wuellerstorfi δ 13 C and B/Ca data (yellow diamonds; comparing the red and blue diamonds from A)) compared to the equivalent data of the last glacial-to-Holocene intervals of cores BOFS 5K (50°42'N, 21°54'W; 3547 m w.d.; B/Ca from C. wuellerstorfi and C. mundulus; δ 13 C from Cibicidoides spp.) and BOFS 8K (52°30'N, 22°06'W; 4045 m w.d.; B/Ca from C. wuellerstorfi and C. mundulus; δ 13 C from C. wuellerstorfi) (gray diamonds; B/Ca data from Yu et al. (52), adjusted by +5% to correct for inter-laboratory differences; δ 13 C data from Yu et al. (52) and references therein). C) The Eirik Drift C. wuellerstorfi δ 13 C (coloring and references as in A)) and N. pachyderma (s) δ 13 C from MIS 5e (MD03-2664; 55) and MIS 7e, 9e, and 11c (Site U1305; this study). The N. pachyderma (s) data were corrected by +1.0‰ to account for disequilibrium (57). D) Cross-plot of Eirik Drift N. pachyderma (s) δ 13 C (also corrected by +1.0‰) versus C. wuellerstorfi δ 13 C (sample average values; data and references as in C)) where the solid line shows the linear fit (r 2 =0.045).

Fig. S7 .

 S7 Fig. S7. Reconstructed and simulated Eirik Drift bottom water δ 13 C variability. Plotted versus years relative the final high values preceding distinct bottom water δ 13 C reductions at the Eirik Drift: A) reconstructed bottom water (C. wuellerstorfi) δ 13 C changes shown as averages (bold lines) of individual events from MIS 1, 5e, 9e, and 11c (thin lines) at 30-year steps (obtained by linear interpolation) and binned according to durations of <100 (red; n=5), 101-200 (blue; n=4), and 201-300 years (green; n=3); for the same two simulated NADW shoaling events (showing ten-year running means): B) the original simulation as shown in Fig. 4; C) the original simulation but with atmospheric δ 13 C decreased by ~0.6‰ and fixed (at -7.1‰); and D) the original simulation but with 50% decreased primary productivity in the northern deep water formation regions off southern Greenland and in the Nordic Seas.

Acknowledgments:

We thank the scientific party and crew of R/V JOIDES Resolution IODP Expedition 303, the curatorial staff at the IODP Bremen core repository for core sampling assistance, and Douglas Richmond for assistance with processing MIS 9e samples. Funding: Funded by grants from the Research Council of Norway (RCN) through the project THRESHOLDS (RCN grant 254964) and data generation at FARLAB (RCN grant 245907).

Author contributions: U.S.N. and E.V.G. designed the project and all co-authors helped conceptualize the study; E.V.G. processed the sediment samples and conducted the benthic and planktic foraminifera stable isotope analyses; N.I. conducted the IRD analyses; A.K., J.T.,

Supplementary Materials:

Materials and Methods

Supplementary text