Research ArticleCLIMATOLOGY

Sea ice variability in the southern Norwegian Sea during glacial Dansgaard-Oeschger climate cycles

See allHide authors and affiliations

Science Advances  06 Mar 2019:
Vol. 5, no. 3, eaau6174
DOI: 10.1126/sciadv.aau6174


The last glacial period was marked by pronounced millennial-scale variability in ocean circulation and global climate. Shifts in sea ice cover within the Nordic Seas are believed to have amplified the glacial climate variability in northern high latitudes and contributed to abrupt, high-amplitude temperature changes over Greenland. We present unprecedented empirical evidence that resolves the nature, timing, and role of sea ice fluctuations for abrupt ocean and climate change 32 to 40 thousand years ago, using biomarker sea ice reconstructions from the southern Norwegian Sea. Our results document that initial sea ice reductions at the core site preceded the major reinvigoration of convective deep-water formation in the Nordic Seas and abrupt Greenland warming; sea ice expansions preceded the buildup of a deep oceanic heat reservoir. Our findings suggest that the sea ice variability shaped regime shifts between surface stratification and deep convection in the Nordic Seas during abrupt climate changes.


Greenland ice core records reveal that the last glacial period was marked by abrupt climate changes, leading to alternating cold Greenland stadials (GS) and warmer Greenland interstadials (GI), known as Dansgaard-Oeschger (D-O) cycles (1, 2). These D-O climate fluctuations involved dramatic changes in northern high-latitude air temperature, hydroclimate, and atmospheric circulation transitioning in a matter of decades (24), and also affected Eurasian (5) and tropical climates (6). There is an emerging consensus that the reconstructed pattern of D-O variability found in various climate proxy records around the world is consistent with the notion of a varying Atlantic meridional overturning circulation (AMOC) and resulting changes in heat transport from the tropics to northern high latitudes (710). Similar to today, a strengthened AMOC with deep-water formation in the Nordic Seas would have enhanced northward heat transport to the Norwegian Sea during warm GI. By contrast, a weakened AMOC during cold GS with deep-water formation south of Iceland would have reduced heat transport to the Nordic Seas.

It is still unclear what triggered the D-O variability, and proposed mechanisms range from coupled ice sheet–ocean circulation variability (710) to volcanic eruptions (11) and solar variability (12). It is, however, clear that the apparent nonlinear behavior of the abrupt D-O climate changes in northern high latitudes must have involved feedbacks in the coupled ocean-climate system (8, 13). Proxy and modeling evidence hints at a potential role of rapid changes in sea ice cover in amplifying and potentially instigating changes in ocean circulation and producing the abrupt D-O climate transitions observed in Greenland ice core records (9, 1420). Moreover, “early-warning” signals reflected by significant variance and autocorrelation increases in the high-frequency band of the δ18O ice core record from Greenland may be linked to sea ice fluctuations before the bifurcation of a D-O warming event (21). However, empirical evidence that reliably resolves the sea ice variability during the millennial-scale climate cycles with a sufficient temporal resolution is lacking until now. Here, we investigate the sea ice changes and its implications across four GS-GI cycles by means of new biomarker sea ice records with an unprecedented multidecadal- to centennial-scale resolution, which allows us to determine the phasing between sea ice and various oceanic and climatic factors.

Sedimentary proxy records were obtained from core MD99-2284 (62°22.48N, 0°58.81W, 1500 m water depth), retrieved from the core of northward inflow of Atlantic waters to the Norwegian Sea (Fig. 1). As suggested by model simulations and available proxy data from this core, this core site is ideally suited to trace glacial changes in sea ice cover and oceanography in the southeastern Nordic Seas during the abrupt GS-GI transitions (9, 14). To investigate the relationships between the ocean and ice core records, we have placed signals in anhysteretic remanent magnetization (ARM) and near-surface temperature in core MD99-2284 and the D-O signals in the δ18O record of the North Greenland Ice Core Project (NGRIP) ice core on a common time scale, and independently verified this by a distinct tephra layer identified in both the sediment and ice core (figs. S1 and S2) (see Materials and Methods). According to this age model, sedimentation rates in MD99-2284 are extremely high [~50 to 400 cm/thousand years (ka)], and enhanced sediment deposition during GI may be linked to a strengthened bottom/deep current activity influencing this sediment drift site (see Materials and Methods). Our new proxy records cover four D-O cycles between ~32 and 40 ka before 2000 CE (ka b2k), with a temporal resolution of 15 to 170 years (table S1). This time interval comprises a suite of well-defined, long- and short-lasting D-O events 5 to 9 that are representative for millennial-scale glacial climate variability and suitable for ocean-ice core comparison.

Fig. 1 Regional context of the study area.

Annual mean sea surface temperature anomaly (SSTA) and 10% spring sea ice extent for interstadials (yellow line) compared to stadials (blue line) as well as interstadial sea ice thickness derived from the LOVECLIM MIS3 hindcast simulation (9). Blue shading north of the yellow line reflect interstadial sea ice thickness and red shading south of the yellow line indicate interstadial SSTA. Interstadial and stadial conditions were averaged across 31.05 to 41.05 ka b2k for GI5-9 and GS6-9, respectively (see Materials and Methods). GSR, Greenland-Scotland ridge; NAC, North Atlantic Current (red arrow); ISO, Iceland-Scotland Overflow (blue arrow). Yellow stars, locations of sediment core MD99-2284 (62°22.48N, 0°58.81W, 1500 m water depth) and the NGRIP ice core used in this study. Map at the bottom showing the modern bathymetry was generated using Ocean Data View software ( (59).

Our sea ice reconstruction relies on paired records of open-water phytoplankton biomarkers brassicasterol and dinosterol, as well as the ice algae biomarker IP25 (see Materials and Methods). Brassicasterol and dinosterol are biomarkers produced during spring and/or summer by open-water diatoms and dinoflagellates, respectively (22). The ice proxy IP25 is a highly branched isoprenoid monoene that is produced only by sea ice diatoms during spring (23). Moreover, we combine coeval abundances of brassicasterol/dinosterol and IP25 using the phytoplankton-IP25 index (PBIP25/PDIP25) to derive semiquantitative estimates of past sea ice cover (24) (see Materials and Methods). Previously published PIP25 values of Arctic and subarctic surface sediments reveal a robust correlation with the modern sea ice concentrations during spring/summer, and PIP25 values of <0.1, 0.1 to 0.5, 0.5 to 0.75, and >0.75 are classified to indicate ice-free conditions, a reduced variable sea ice cover, a seasonal sea ice cover, and an extended to perennial sea ice cover, respectively (24, 25).


Biomarker proxy records and sea ice reconstruction

The down-core records of MD99-2284 show distinct variations in open-water phytoplankton biomarker fluxes parallel to Greenland’s D-O cycles ~32 to 40 ka (Fig. 2 and fig. S3). Near-zero dinosterol and brassicasterol fluxes reflect an almost halted open-water phytoplankton production during GS. Highly elevated fluxes of dinosterol and brassicasterol reflect enhanced open-water phytoplankton productivity during GI (least pronounced for GI5), together with increased input of terrigenous material (fig. S3). In contrast to the open-water phytoplankton biomarkers, the IP25 flux reveals a variability slightly offset with respect to the D-O cycles in the NGRIP δ18O record (Fig. 2D). IP25 fluxes increase throughout the latter half of GI and peak at the GI-GS transition, suggesting an increasing sea ice algae production. Subsequently, IP25 fluxes gradually decrease and reach minimum values in mid GS, indicating a reduction of ice algae productivity. IP25 fluxes and thus sea ice algae production remain low across GS-GI transitions, while dinosterol and brassicasterol fluxes reveal rapid increases at the transitions. Distinct IP25 peaks mark the initial increase in ice algae productivity after the early GI peak warm period in Greenland (Fig. 2).

Fig. 2 Biomarker proxy records from core MD99-2284.

(A) NGRIP ice core δ18O plotted as 11-point running average (2, 4). (B) Dinosterol flux. (C) Brassicasterol flux. (D) IP25 flux. (E) PBIP25 and PDIP25. (F) Simulated spring sea ice cover averaged over 60°N-62.5°N, 1.25°W-1.25°E (blue; core site) and 60°N-80°N, 15°W-20°E (gray; Norwegian Sea) (9). All records are plotted on the GICC05 b2k time scale. GI numbered at the top. GS numbered at the bottom (shaded bars). Light blue bars, intra-interstadial sea ice expansion events.

The presence of IP25 in all samples investigated suggest a continuous occurrence of spring sea ice in the southern Norwegian Sea throughout the time period studied (within the sample resolution), with its variations in its abundance pointing to changes in seasonal sea ice cover. Increased IP25 with contemporaneously decreasing brassicasterol and dinosterol indicate most extensive seasonal sea ice conditions at GI-GS transitions and following early parts of GS in our records. These times are also marked by resulting maximum PBIP25 and PDIP25 values of ~0.7 to 0.8 (Fig. 2E). The biomarker signals of enhanced sea ice cover during late GI/early GS are similar to those observed in sediment traps and surface sediments on the modern proximal East Greenland margin (25, 26). Here, the sea ice season extends from late autumn to spring such that ice algae growth is enhanced and open-water phytoplankton growth is reduced due to only short periods of ice-free conditions during summer (2527). In turn, lowered IP25 with contemporaneously increased brassicasterol and dinosterol values, leading to minimum PBIP25 and PDIP25 values of ~0.2, reflects periods of maximum open-ocean conditions in our record (Fig. 2). These maximum open-ocean conditions in the southern Norwegian Sea appear to have persisted only during peak GI, and similar conditions have been reconstructed for the warm early Holocene based on a PIP25 record from a nearby core site (15). These biomarker signals are similar to those observed in the modern Nordic Seas south of the Fram Strait, where sea ice only occurs sporadically and well-mixed surface waters stimulate the phytoplankton production (2426).

The biomarker signals during the younger GS in our records document decreasing IP25 fluxes that may reflect reductions in spring sea ice cover. On the other hand, contemporaneously lowered brassicasterol and dinosterol fluxes indicate reduced open-water phytoplankton production. This could be interpreted as (near-)perennial sea ice cover, which might have been too thick and persistent such that ice algae and phytoplankton growth would cease due to limited light and nutrient availability (25, 26). In that case, one may expect a transition through a seasonal sea ice phase at GS-GI transitions, which the IP25 record does not reveal (Fig. 2D). Therefore, we instead interpret the lowered IP25 fluxes during the latter half of GS as a shortening of the winter sea ice season and thus less spring sea ice occurrence at the core site, probably associated with a northward migration of the sea ice edge into the southern Norwegian Sea. We argue that the open-water phytoplankton productivity was still diminished, probably because sea ice melting during spring and summer freshened the surface ocean, thereby effectively preventing surface-ocean mixing and phytoplankton fertilization, potentially similar to the modern northern Barents Sea. Variable and generally decreasing PBIP25 and PDIP25 values likewise reflect a trend of a reduction in sea ice cover during the latter half of GS (Fig. 2E). This trend is consistent for all late GS in both the PBIP25 and PDIP25 records, except for late GS9 where PDIP25 reaches a value of 1. This is due to the absence of dinosterol in four samples. In the following discussion, we focus on the consistent trends in IP25 and PBIP25. The actual GS-GI transitions are then marked by rapid increases in dinosterol and brassicasterol fluxes and by the final drops in PBIP25 and PDIP25, which indicates a rapid shift to enhanced open-water phytoplankton production (Fig. 2E). These observations imply a further sea ice reduction, that is, a further northward shift of the sea ice edge in the Norwegian Sea.

After the maximum open-ocean conditions during peak GI, interstadials are characterized by increasing IP25 and PIP25, which points to increasing seasonal sea ice conditions before the maximum seasonal sea ice cover is reached at the GI-GS transitions (2427). These conditions are analogous to those prevailing at the distal East Greenland margin today, where seasonal sea ice promotes ice algae growth and IP25 production in spring, while surface ocean mixing in summer and autumn strongly stimulates the open-water phytoplankton production and biomarker fluxes due to enhanced nutrient availability (24, 25). The trend of increasing sea ice cover throughout interstadials is initiated by short-term peaks in IP25 and PIP25, which we interpret as reflecting brief, hitherto unknown intra-interstadial sea ice expansion events following the peak warm periods during GI6 to GI8 (Fig. 2).

Reconstructed and model-simulated sea ice variability during GS-GI cycles

Our biomarker records reflect the sea ice variability in the southern Norwegian Sea during GS-GI cycles in unprecedented detail, implying the presence of the most extensive seasonal sea ice conditions at the onsets and early parts of cold GS, and the most pronounced open-ocean conditions at the onsets of warm GI. Comparable extreme values over all GS-GI cycles in our PIP25 records indicate that the extensive sea ice conditions were relatively similar for all early stadials, and reduced sea ice conditions were similar for all peak interstadials. On the other hand, increased IP25 fluxes point to more severe seasonal sea ice formation during shorter GI6 and GI7 compared to longer-lasting GI8 (Fig. 2). Offset variations in IP25 and open-water phytoplankton biomarkers suggest that the seasonal sea ice diminished through GS and expanded in the mid to late parts of the GI. We interpret this as reflecting shifts of the (spring) sea ice edge, between a position south of core site MD99-2284 during early GS and one north of the core site during peak GI, in a manner that is offset with Greenland warming and cooling events.

This GS-GI sea ice variability is, in part, consistent with evidence from a lower-resolution record by Hoff et al. (15), who inferred enhanced open-water conditions during peak GI and sea ice expansion during later GI for a core site ~150 km west of MD99-2284. Different from our interpretation, however, Hoff et al. (15) argued that there was a near-perennial sea ice cover in the southern Norwegian Sea throughout stadials, which rapidly declined at GS-GI transitions. This might be related to the lower temporal resolution of their record and/or to the more western core location further away from the warm Atlantic inflow to the Norwegian Sea. Trends in spring sea ice occurrence derived from our IP25 record appear comparable with sea ice estimates based on dinoflagellate cysts from a core site also located ~150 km west of MD99-2284 (28). On the basis of this record, Wary et al. (28) proposed intensive winter sea ice formation and summer sea ice melting during GI and reduced sea ice formation during GS. Compared to our reconstruction, however, their record with a markedly lower temporal resolution did not resolve maximum open-water conditions during early peak GI.

Furthermore, our reconstructed trends in sea ice cover agree, in part, with results from a transient D-O numerical experiment performed with the Earth system model LOVECLIM (see Materials and Methods) (9). The model experiment simulates Greenland D-O–like cycles, with a sea ice edge located south of 62°N during GS and largely open-ocean conditions in the Norwegian Sea during GI, which is consistent with our reconstructed maximum and minimum sea ice cover (Fig. 1). While the proxy reconstruction suggests spring sea ice cover during the early part of all GS, the spring sea ice cover does not reach our core site during GS6 and GS7 in the simulation (Fig. 2F). During late GS9 and late GS8, the simulated spring sea ice retreat at the core location precedes the larger-scale sea ice decline in the eastern Nordic Seas at the onsets of GI8 and GI7 by about ~700 and ~300 years, respectively (Fig. 2F). This simulated early sea ice retreat in the southern Norwegian Sea is consistent with our biomarker-based sea ice reconstruction, in terms of the overall trend of reduced spring sea ice cover before the GI and also roughly with respect to the timing of decreasing IP25 fluxes (Fig. 2). In contrast to our IP25 record, however, the simulation does not show any spring sea ice cover at the core site during the later parts of the GI, which indicates that the forcing might be incorrect in detail.

Linkages between glacial ocean circulation, sea ice cover, and climate

A diminished biological export production during GS is also supported by lowered planktic δ13C signals in MD99-2284, which has been interpreted as reflecting poorly ventilated subsurface waters and thus surface stratification in the Norwegian Sea (Fig. 3) (14). Consistent with a reduced ocean-atmosphere heat flux due to surface ocean stratification and the extended, though decreasing, seasonal sea ice cover in the Nordic Seas during GS, reduced benthic δ18O in core MD99-2284 suggests a subsurface warming (Fig. 3), in line with independent low-resolution evidence of intermediate- and deep-water warming from nearby sites (29, 30). Enhanced surface stratification and increased deep-water temperatures in the Nordic Seas suggest that convective deep-water formation was inactive in the Norwegian Sea during all GS investigated here, and probably shifted to an area south of the Greenland-Scotland ridge near the sea ice edge (8, 15, 29, 30). This likely contributed to a general reduction of the meridional ocean circulation across the Greenland-Scotland ridge, as indicated by a lowered ARM (see Materials and Methods) (31) and in agreement with deep North Atlantic evidence of a weakened stadial AMOC (10). We argue that impaired ocean-atmosphere heat exchange under an extended sea ice cover in the Nordic Seas might have played a crucial role for subsurface warming in the Nordic Seas and temperatures in Greenland being equally low during all GS in ice core records.

Fig. 3 Sea ice cover and other paleoceanographic records from core MD99-2284.

(A) NGRIP ice core δ18O plotted as 11-point running average (2, 4). (B) PBIP25. (C) IP25 flux. (D) Near-surface temperature based on transfer function estimates using planktic foraminifera census counts (14). (E) Planktic foraminiferal δ13C (14). (F) Benthic foraminiferal δ18O (14). (G) ARM (14). All records are plotted on the GICC05 b2k time scale. Thick lines in (B), (D), and (E), five-point running averages. Diamonds at the bottom indicate tie points used for the age model. GI numbered at the top. GS numbered at the bottom (shaded bars). Light blue bars, intra-interstadial sea ice expansion events.

We studied the phase relationships between changes in sea ice, near-surface, and deep-water temperatures with respect to Greenland temperature using cross-correlation analyses (see Materials and Methods). According to our age model, which is based on tuning the ARM and a pronounced near-surface temperature overshoot to NGRIP δ18O and verified by the tephra analyses, both ARM and benthic δ18O are positively correlated and in phase with NGRIP δ18O (R = 0.66 and R = 0.71, respectively) (Fig. 4). Using independent temperature estimates from Mg/Ca of benthic foraminifera from a nearby site, Ezat et al. (30) showed that the benthic δ18O chiefly reflects deep-water temperature changes of opposite sign, but in phase with Greenland temperature. Accordingly, increased benthic δ18O during GI indicates a 2° to 5°C cooling of deep waters <1200 m, which was likely linked to active convective deep-water formation in the Nordic Seas, similar to modern conditions (30).

Fig. 4 Relative timing of stadial-interstadial changes in the proxy records.

Cross-correlations of ARM (black), benthic foraminiferal δ18O (dark blue), IP25 (purple), near-surface temperature (red), and PBIP25 from MD99-2284 with NGRIP ice core δ18O (see Materials and Methods). The ARM was tuned to NGRIP δ18O for the age model development, resulting in a positive and in-phase correlation between the two. Benthic δ18O variations, reflecting deep-water temperature changes related to ocean convection, are in phase with NGRIP δ18O. Asymmetric and offset patterns of cross-correlations of IP25, near-surface temperature, and PBIP25 with NGRIP δ18O suggest distinctly different temporal dynamics of changes in seasonal sea ice cover and near-surface temperature as compared with major deep-ocean convection and Greenland climate changes.

The cross-correlations of IP25 and near-surface temperature with NGRIP δ18O (R = 0.56 and R = −0.59, respectively, with a lag of ~200 years) reveal patterns that are distinctly different from the cross-correlation of benthic δ18O with NGRIP δ18O, which seem to mimic each other (Fig. 4). This supports both the asynchronous variability of seasonal sea ice cover with respect to deep-water temperature in the Norwegian Sea and a (causal) link between spring sea ice and near-surface temperature. PBIP25 and NGRIP δ18O are inversely correlated with maximum anticorrelation occurring at a PBIP25 lead time of ~0 to 250 years (R < −0.6), likely reflecting a slight offset between sea ice changes and Greenland climate transitions (Fig. 4). The early sea ice reductions, reflected by reduced IP25 fluxes and decreasing PIP25 values during late GS, may thus have resulted from an increasing influence of warm and saline Atlantic waters flowing into the Norwegian Sea at the surface and/or subsurface. This is supported by salinity estimates derived from the δ18O of planktic foraminifera, reflecting an increasing salinity throughout GS, which point to a gradual weakening of the stadial halocline and surface stratification at the core site (fig. S4) (14).

The benthic δ18O increase is paralleled by an abrupt ~3°C overshoot in near-surface temperature in core MD99-2284 at all GS-GI transitions (Fig. 3) (14). The ~3°C temperature overshoot indicates maximum inflow of warm and saline North Atlantic waters into the Norwegian Sea, consistent with largely open-ocean conditions during early peak GI. It is coeval with reductions in planktic and benthic δ13C signaling admixture of poorly ventilated warmer waters from depth to surface and intermediate waters (Figs. 3 and 5) (14). We interpret these signals as recording an overturning pulse mixing the entire water column after the sea ice retreat at our core site. These overturning pulses most likely released vast amounts of oceanic heat, reflecting a major reinvigoration of convective deep-water formation in the Nordic Seas (14). This convective deep-water formation likely cooled local deep waters, strengthened the Nordic Seas inflow/outflow circulation, and stimulated the phytoplankton productivity, as reflected by the large rises in benthic δ18O, ARM, and brassicasterol (Figs. 3 and 5) (2931).

Fig. 5 High-resolution proxy records covering the GS9-GI8 transition.

(A) NGRIP ice core δ18O (thick line is the 11-point running average) (2, 4). (B) PBIP25 (blue line), IP25 concentration (purple), and brassicasterol concentration (green). (C) Near-surface temperature based on transfer function estimates using planktic foraminifera census counts (14). (D) Planktic foraminiferal δ13C (14). (E) Benthic foraminiferal δ13C (14). (F) Benthic foraminiferal δ18O (14). (A) is plotted on the GICC05 b2k time scale at the top. (B) to (F) are plotted on the MD99-2284 core depth scale at the bottom, which is aligned to the GICC05 time scale using the GS9-GI8 tuning point and is floating in deeper and shallower parts. Dashed lines mark the positions of tephra layers with a statistically similar geochemical composition identified in the NGRIP ice core (green) and core MD99-2284 (blue) (fig. S2). GI8 and GS9 (shaded bar) are indicated for NGRIP at the top and for MD99-2284 at the bottom.

A distinct tephra layer, which geochemically correlates with a tephra layer in the NGRIP ice core, constrains the end of the overturning pulse in MD99-2284 and the early peak warmth in Greenland during peak GI8 (Fig. 5 and fig. S2) (see Materials and Methods). Figure 5 illustrates that the major IP25 drop and parallel onset of PBIP25 decline is situated 36 cm below the near-surface temperature overshoot at the GS9-GI8 transition. This implies that the initial seasonal sea ice reduction in the southern Norwegian Sea preceded the major overturning pulse in the Nordic Seas and the abrupt D-O warming in Greenland. Although particle size–dependent bioturbational mixing could potentially contribute to an offset between the transitions in the biomarker- and foraminifera-based proxy records, the clear separation of the signals and sedimentation rates exceeding 58 cm/ka would strongly minimize such an effect (32). Moreover, brassicasterol increases together with the near-surface temperature overshoot after the IP25 drop, which suggests that the biomarker fluxes chiefly record an autochthonous regional signal of phytoplankton and ice algae production and that a potential additional deposition of laterally transported biomarkers would not affect our main conclusions (Fig. 5). If the IP25 drop at 3120.5 cm was placed to match the GS9-GI8 transition in the NGRIP δ18O, this would lead to an unrealistically high sedimentation rate exceeding 600 cm/ka, given the temporal constraint by the tephra layer. Hence, the time separation of reconstructed parameters within the core and their relation to the warming recorded in the NGRIP δ18O appears robust.

We argue that the rapid reinvigoration of deep convection and concomitant heat loss to the atmosphere were closely tied to an intensification of the large-scale AMOC and the threshold-like abrupt warming of the D-O events recorded in Greenland ice core δ18O, which is dynamically supported by model simulations (8, 9). Furthermore, the sea ice decline starting notably before the GS6-GI5 transition was accompanied by several jumps in near-surface temperature and an early deep-water cooling, which may reflect the flickering of convection close to a threshold point, before the main overturning pulse occurs (Fig. 3). We therefore suggest that both major reinvigoration of convective deep-water formation in the Nordic Seas and abrupt Greenland warming at GS-GI transitions occurred after the sea ice extent had crossed a minimum threshold. This minimum threshold probably implies largely ice-free conditions in the southern Norwegian Sea, consistent with the LOVECLIM model simulation (Fig. 2).

Following the overturning pulse and near-surface temperature overshoot during peak GI, increased planktic δ13C and decreased near-surface and deep-water temperatures document well-ventilated and cooled waters linked to ocean mixing (Fig. 3) (14, 29, 30). As the near-surface temperatures cooled after the overshoot, increased IP25 fluxes suggest enhanced seasonal sea ice formation. The intra-interstadial sea ice expansion events and coeval reductions in Nordic Seas inflow/outflow, as documented by sharp ARM drops, may also have contributed to gradual Greenland cooling after the initial peak GI warmth. This indicates that strong deep-water formation under open-ocean conditions in the Norwegian Sea was only temporarily stable (8). The Scandinavian and British ice sheets surrounding the glacial Nordic Seas formed large reservoirs of freshwater that, when released, potentially affected the stability of the deep convection and thus sea ice formation. Decreasing salinity estimates reconstructed from the planktic δ18O at site MD99-2284 imply a local interstadial surface freshening, supported by increased fluxes of terrigenous material (figs. S3 and S4) (14).

The freshening likely weakened deep convection, reestablished surface stratification, and ultimately fostered sea ice expansion in the Norwegian Sea during the later phases of GI (Fig. 3). This is also in line with proxy and model evidence showing that freshwater release from the Scandinavian Ice Sheet resulted in enhanced sea ice formation in the eastern Nordic Seas just preceding the onset of GS1, the Younger Dryas (33). The increasing seasonal sea ice cover and associated sea ice export during late GI may also have contributed to a freshening of surface waters in the North Atlantic (20). We hypothesize here that once the surface freshening and sea ice cover had reached a threshold, convective deep-water formation ceased in the Nordic Seas and relocated to the south of the Greenland-Scotland ridge. This may have contributed to an AMOC reduction and defined the rapid transition of the North Atlantic climate system to cold GS conditions (8, 9, 34).


Our new biomarker proxy results provide unprecedented, detailed empirical evidence of persistent sea ice variations in the southern Norwegian Sea during D-O climate cycles, slightly offset with Greenland temperature changes. This indicates that sea ice reductions and expansions in the southern Norwegian Sea had distinctly different temporal dynamics as compared to the abrupt Greenland climate transitions during the last glacial. Our relatively consistent proxy-model data comparison suggests that freshwater fluxes to the North Atlantic and/or Nordic Seas are a probable forcing to drive the D-O climate variability. Enhanced northward advection of warm and saline Atlantic surface waters may have acted as trigger, or contributor, of initial sea ice retreat in the Norwegian Sea. Northward propagation of Atlantic surface waters might have been induced by a critical salinity gradient between the saline North Atlantic and the fresher halocline in the Nordic Seas (35), stochastic atmospheric forcing affecting the subpolar gyre circulation (36), wind-driven advection of saline surface waters from the tropical North Atlantic (37), and/or an early AMOC strengthening (10). More detailed proxy and model evidence is needed to better constrain mechanisms that initially triggered the early sea ice reductions in the Norwegian Sea observed in our data. The demise of the sea ice cover at the GS-GI transitions appears connected with the release of oceanic heat to the atmosphere and concomitant convective overturning.

Vettoretti and Peltier (19, 20) used an unforced model simulation performed with the Community Earth System Model version 1 (CESM1), which has revealed D-O–like climate oscillations, to investigate the mechanisms involved in triggering the onset of a D-O event. In this unforced CESM1 model simulation, abrupt warming in Greenland results from oceanic heat release after an initial sea ice disappearance associated with the opening of a super-polynya (19, 20). The rapid disappearance of sea ice in this model simulation occurs first south of Greenland, and subsequent large-scale sea ice retreat in the North Atlantic and Norwegian Sea happens within years, different from our proxy results. The authors of these studies argued that the opening of the super-polynya is triggered by a thermohaline convective instability beneath the sea ice cover (19, 20). Such thermohaline instability can result from an extensive heat reservoir beneath the sea ice and a gradually diminishing vertical salinity gradient enabled by a continuous salinity increase above the pycnocline (19, 20, 38). Although there is proxy evidence of a stadial subsurface heat reservoir in the Norwegian Sea, this mechanism might not solely explain the asynchronous sea ice variability we observe in our records.

Our sea ice record reveals increased seasonal sea ice cover in the southern Norwegian Sea during the latter half of warm GI, probably resulting from a local surface ocean freshening and/or cooling linked to reduced northward advection of Atlantic waters. While freshwater fluxes applied as forcing in the transient LOVECLIM simulation and salinity estimates based on planktic δ18O roughly show similar trends of interstadial freshening (fig. S4), the model does not reveal increasing spring sea ice in the southern Norwegian Sea during GI (Fig. 2) (9). This might suggest that freshwater input from the Scandinavian and British ice sheets into the Nordic Seas is underrepresented in the model experiment, but likely contributed substantially to the buildup of surface ocean stratification and sea ice cover in the Norwegian Sea.

We conclude that the glacial sea ice variability in the Norwegian Sea acted as an important factor that contributed to abrupt glacial changes in ocean circulation and climate in the North Atlantic and over Greenland. We suggest that sea ice reduction and resumption of deep-water formation in the Nordic Seas may have played a persistent role in triggering the abrupt GS-GI climate transitions, terminating both “normal” stadials and “most extreme” stadials that were characterized by severe iceberg surges and collapse of deep-water formation in the North Atlantic and extreme AMOC breakdowns, as, for example, associated with Heinrich Event 4 during GS9 (8, 10). This corroborates previous hypotheses that relate the D-O warming events observed in Greenland ice cores to changes in Nordic Seas sea ice cover. Although the initial sea ice retreat during late GS apparently did not affect Greenland’s temperature, it might explain early-warning signals for D-O events, which have been identified in the δ18O ice core record from Greenland (21). We note that the sea ice retreat during late GS may have been a local feature of the easternmost Atlantic inflow to the southern Norwegian Sea, which likely was not sensitive enough to noticeably affect Greenland’s cold stadial temperatures. Further high-resolution sea ice proxy evidence from the central Norwegian Sea and North Atlantic is desirable to attest to the larger-scale glacial variability in sea ice extent during the millennial-scale D-O climate changes. Last, our findings could also be of relevance concerning the stability of Arctic sea ice and its underlying ocean stratification and warming of the inflow waters in view of the ongoing and future Arctic warming and sea ice demise (39).


Biomarker analyses

Biomarker extraction was carried out on 2.5 to 5.6 g of freeze-dried and homogenized sediment using an accelerated solvent extractor (DIONEX, ASE200; 100°C, 5 min, 1000 psi) with dichloromethane:methanol (2:1, v/v) as solvent. Before extraction, the internal standards 7-hexylnonadecane (7-HND, 0.076 μg, for IP25 quantification), squalane (2.4 μg), and cholesterol-D6 (cholest-5-en-3β-ol-D6, 10 μg, for sterol quantification) were added to each sample. Extracts were separated into hydrocarbon and sterol fractions through open-column chromatography (SiO2) using 5 ml of n-hexane and 6 ml of n-hexane:ethylacetate (2:8, v/v), respectively. The sterol fraction was silylated with 200 μl of BSTFA (bis-trimethylsilyl-trifluoroacetamide) (60°C, 2 hours).

Biomarkers were analyzed by gas chromatography/mass spectrometry using an Agilent 7890B GC coupled to an Agilent 5977A mass selective detector (MSD) in selected ion monitoring (SIM) mode for IP25 and an Agilent 6850 GC coupled with an Agilent 5975C VL MSD for sterols. Component identification was based on comparison of GC retention times with those of reference compounds and published mass spectra for sterols (22, 40) and IP25 (23).

IP25 was quantified using its molecular ion (m/z 350) in relation to the abundant fragment ion m/z 266 of the internal standard 7-HND (26). Dinosterol (4α-23,24-trimethyl-5α-cholest-22E-en-3β-ol), brassicasterol (24-methylcholesta-5,22E-dien-3β-ol), campesterol (24-methylcholest-5-en-3β-ol), and β-sitosterol (24-ethylcholest-5-en-3β-ol) were quantified as trimethylsilyl ethers using molecular ions m/z 500, m/z 470, m/z 472, and m/z 486, respectively, in relation to the molecular ion m/z 464 of the internal standard cholesterol-D6. The detection limit for quantification of IP25 is 0.005 ng μl−1 in SIM mode. Instrument stability and analytical accuracy were controlled by replicate analyses of external standards (7-HND, cholesterol-D6, reference sediment) during each analytical sequence and by replicate analyses for random samples. These replicate analyses indicate an SD (1σ) of 2 ng gSed−1 for a mean sterol concentration of 200 ng gSed−1 and an SD (1σ) of 0.3 ng gSed−1 for a mean IP25 concentration of 5 ng gSed−1, indicating the feasibility of repeatable measurements.

Biomarker concentrations were normalized to the total organic carbon (OC) content that was measured on ~100 mg of sediment from each sample used for biomarker analyses. To obtain biomarker fluxes, component-specific concentrations (per gram sediment) were multiplied with estimates of the sedimentation rate (table S1) and dry bulk density. To generate semiquantitative estimates of the past sea ice cover, we calculated the phytoplankton-IP25 index established by Müller et al. (24)PIP25=[IP25]/([IP25] + ([phytoplankton marker] × c))(1)

The balance factor c corresponds to the ratio of mean IP25 and mean phytoplankton concentration for MD99-2284 records presented in this study (c = 0.033 for dinosterol, c = 0.0056 for brassicasterol).

Foraminifer analyses

Near-surface temperature estimates were based on planktic foraminifera census counts and the maximum likelihood transfer function technique (14). These near-surface temperature estimates were based on a calibration using modern temperatures at 10 m water depth during summer (July-August-September), but trends in its temporal variability might be representative for several tens of meters of the upper ocean (41). Stable oxygen and carbon isotopes (δ18O and δ13C) were measured on planktic, near-surface dwelling Neogloboquadrina pachyderma (sinistral) and epibenthic (or shallow infaunal) Cassidulina neoteretis. More methodological details are given by Dokken et al. (14).

Cryptotephra analyses

Cryptotephra investigation was performed at a sample spacing of 4 cm in the interval between 3020 and 3180 cm in core MD99-2284. This interval comprises the GS9-GI8 transition, for which several, mainly basaltic tephra layers have been identified in Greenland ice cores and geochemically characterized (42, 43). To extract basaltic glass shards, freeze-dried and homogenized sediments were immersed in 10% HCl to remove carbonate material and subsequently wet-sieved in three size fractions (>125, 80 to 125, and 25 to 80 μm) (44). For specific cryptotephra analysis, the 25- to 80-μm fraction was further separated into different density fractions (<2.3, 2.3 to 2.5, and >2.5 g/cm3) to isolate the basaltic glass shards using a heavy liquid technique (45). Furthermore, the >2.5 g/cm3 fraction underwent a magnetic separation procedure to separate paramagnetic basaltic from minerogenic material (46). The samples were then mounted in Canada Balsam on microscope slides, and basaltic glass shards were quantified using optical microscopy.

Geochemical analysis was performed on basaltic shards from 3040 to 3041 cm, where a distinct peak in shard concentration was found (fig. S2). Selected shards were mounted in epoxy resin on a frosted microscope slide, which was then grinded using decreasing grades of silicon carbide paper and polished using diamond suspension. Electron probe microanalysis was performed at the Tephra Analytical Unit, University of Edinburgh. A Cameca SX-100 electron microprobe with five wavelength-dispersive spectrometers was used to obtain the major element composition of 37 individual shards (47). The instrumental drift, precision, and accuracy were determined by running BCR2g basalt secondary standards, and error bars in fig. S2 (B and C) represent two SDs of replicate analyses of the BCR2g reference glass. To compare the tephra layer in MD99-2284 with published tephra layers from Greenland ice cores, the raw geochemical data (provided online) were normalized to an anhydrous basis (i.e., 100% total oxides). The geochemistry of the marine and ice core tephra layers was statistically evaluated using the similarity coefficient (SC) (48) and the statistical distance test (D2) (49).

Transient model simulation

Transient global hindcast simulations were performed with the three-dimensional intermediate complexity Earth system model LOVECLIM (9, 50) to investigate the millennial-scale climate variations of the last glacial [50 to 30 ka ago (ka B.P.)]. Initial conditions for the transient run were based on an equilibrium spin-up simulation with an atmospheric CO2 content of 207.5 ppmv (parts per million by volume) and an orbital forcing and estimated ice sheet orography for 50 ka B.P. (9). In the transient run, these glacial boundary conditions were continuously updated and a time-varying anomalous freshwater flux to the North Atlantic (55°-10°W, 50°-65°N) was applied as forcing (9). The time scale of the model data was originally placed on the GICC05 (B.P.) chronology (9). For the present study, we used simulated stadial-interstadial anomalies of annual mean sea surface temperature, 10% spring sea ice extent, and sea ice thickness for the North Atlantic. To calculate stadial-interstadial anomalies displayed in Fig. 1, we averaged stadial conditions over 39.9 to 38.46, 36.34 to 35.5, 34.5 to 33.7, and 33.1 to 32.45 ka B.P., and interstadial conditions over 40.2 to 40.0, 38.34 to 37.9, 35.35 to 34.9, 33.6 to 33.4, and 32.38 to 32.16 ka B.P. For the model-proxy data comparison, we used simulated time series of spring (April-May-June) sea ice cover averaged over 60°N-62.5°N, 1.25°W-1.25°E and 60°N-80°N, 15°W-20°E. The time scale of simulated time series of sea ice cover shown in Fig. 2F was transferred to the GICC05 (b2k) chronology (51, 52).

Age model

We refined the age model of MD99-2284 previously presented by Dokken et al. (14) for the studied interval. Synchronization with the Greenland ice core is based on tuning of coherent signals in both cores, supported by geomagnetic events and tephrochronology as follows:

1) We tuned the MD99-2284 depth scale to the GICC05 (b2k) ice core chronology (4, 51, 52), using distinct transitions in the sedimentary records of ARM and the near-surface temperature overshoot, as well as pertinent transitions in the NGRIP δ18O record (fig. S1 and table S1) (2). It has been shown that ARM variations in Nordic Seas cores appear to be synchronous with the D-O cycles in Greenland ice core δ18O, resulting from changes in the efficiency of magnetic particle transport linked to deep- and intermediate-ocean currents (31, 53). Increased deposition of magnetic particles at core site MD99-2284 during GI may be linked to enhanced deep convection and bottom currents recirculating within the Nordic Seas (31, 53). We note that the erosional center of the deep overflow streaming southward into the North Atlantic would be located in the narrow and shallow parts of the central Faroe-Shetland Channel, away from site MD99-2284 situated at 1500 m water depth north of the sill. We speculate that increased transport of magnetic particles to core site MD99-2284 in the northeastern exit of the Faroe-Shetland Channel might also be driven by enhanced northward inflow of intermediate Atlantic waters. The inflow currents would need to extend to a water depth of a few hundreds of meters to erode magnetite-containing sediments from the eastern flank of the Faroe-Shetland Channel with its volcanic bedrocks, a situation comparable to today (54). Unfortunately, we cannot disentangle whether increased ARM signals in the glacial section of core MD99-2284 are linked to enhanced deep-water circulation, intermediate inflow, or both. However, since inflow and deep overflow dynamics across the Greenland-Scotland ridge are compensating each other (54), we argue that the sedimentary ARM signals at core site MD99-2284 generally reflect the meridional ocean circulation strength in the Faroe-Shetland Channel.

The steep ARM increases parallel near-surface temperature overshoots, which provides further evidence of warm Atlantic waters and oceanic reorganizations that we tied to atmospheric warming of the D-O events (fig. S1). We applied a change point analysis (55) to objectively derive change point probabilities for all records used for the tuning. The final age model was based on 11 tie points and linear interpolation between them (table S1).

2) The tuning-based age model is further supported by the positions of the Laschamp and Mono Lake geomagnetic excursions (fig. S1), which derive from minima in natural remanent magnetization (NRM)/ARM (56) in MD99-2284 and maxima in ice core cosmogenic nuclides (57).

3) Moreover, we identified a pronounced basaltic cryptotephra layer at 3040 to 3041 cm in core MD99-2284 (fig. S2A). Both its stratigraphic position (peak GI8) and its geochemical composition (i.e., Grimsvötn source) suggest that this tephra layer falls within the Faroe Marine Ash Zone III previously found in other marine cores (58) and the high-resolution NGRIP ice core (fig. S2B) (42, 43). If compared with geochemical signatures of several well-investigated tephra layers from the GS9-GI8 transition in the NGRIP ice core (42, 43), the homogeneous geochemical composition of the MD99-2284 tephra layer at 3040 to 3041 cm is statistically most similar to that of the NGRIP tephra layer at 2065.65 m (38.081 ka b2k) (SC = 0.98, D2 = 3.57). The tuning-based age model of MD99-2284 results in an age of 37.955 ka b2k for the tephra layer at 3040 to 3041 cm, which is in good agreement with the age of 38.081 ka b2k given for the pertinent NGRIP tephra layer (42, 43), thus independently validating our age model. In essence, this distinct tephra layer forms the ultimate support that the parallel near-surface temperature overshoot and final ARM rise are most likely closely tied to the abrupt Greenland warming (fig. S2).

Statistical analyses

For the statistical analyses, all records were linearly interpolated to produce evenly spaced time series at intervals of 30 years and subsequently detrended and high-pass–filtered to remove frequencies higher than 1/100 years and thus reduce the noise (fig. S4). Cross-correlations were calculated using the R-Studio software and the package “stats”. The stadial-interstadial tuning points used for the age model were based on change points determined for the original (unsmoothed) records of δ18O in NGRIP and ARM and near-surface temperature in MD99-2284, using the R-Studio software and the “bcp” package (55).


Supplementary material for this article is available at

Table S1. Tie points for the age model of core MD99-2284.

Fig. S1. Age model of core MD99-2284 32 to 40 ka ago.

Fig. S2. Tephra data supporting the age model of core MD99-2284.

Fig. S3. Supplementary biomarker proxy records from core MD99-2284.

Fig. S4. Proxy-model data comparison.

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.


Acknowledgments: We thank the IMAGES program and the R/V Marion Dufresne crew on leg MD114. Moreover, we are grateful to W. Luttmer for assistance with biomarker analysis andto C. Hayward for assistance with geochemical tephra analysis using the electron microprobe at the Tephrochronology Analytical Unit, University of Edinburgh. We thank D. Naafs and two anonymous reviewers for their constructive comments leading to improvement of the manuscript. Funding: The research leading to these results has received funding from the European Research Council under European Union’s Seventh Framework Programme (FP7/2007-2013)/ERC grant agreement number 610055. L.M. is supported by Australian Research Council DECRA project DE150100107. A.T. received funding from the Institute for Basic Science, South Korea. Author contributions: T.M.D., E.J., and H.S. designed the study. H.S. conducted biomarker measurements with guidance from K.F. and R.S., analyzed the data, produced the figures, and drafted the original manuscript. S.M.P.B. counted and analyzed the tephra shards. F.M. contributed to the statistical analyses. L.M. and A.T. provided and interpreted the model output data. All authors contributed to the final 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. Additional data related to this paper may be requested from the authors. Data reported in this study can be found on the PANGAEA database (

Stay Connected to Science Advances

Navigate This Article