Orbital forcing of the East Antarctic ice sheet during the Pliocene and Early Pleistocene

Integrated Ocean Drilling Program (IODP) core U1361 was recovered from ～3,000 m water depth on the continental rise adjacent to the Wilkes Subglacial Basin of Antarctica, one of the largest marine-based sectors of the East Antarctic ice sheet (EAIS). Landward deepening reverse slope troughs penetrate beneath the EAIS, reaching depths of up to 2 km below sea level in this region¹⁰, heightening the vulnerability of this sector of the EAIS to marine ice sheet instabilities (Fig. 1). U1361 provides a well-dated, continuous geologic archive of Pliocene and Early Pleistocene orbital-scale variability of the EAIS marine margin. Sediment deposition at this site is controlled by the interplay between: downslope marine sediment gravity flows; the rainout of biogenic detritus from surface water plankton; iceberg rafting of terrigenous sediments; and low-energy bottom currents (Supplementary Information).

Also shown is the location of the Miocene–Pleistocene ANDRILL AND-1B core recovered in the northwestern corner of the Ross Ice Shelf, the southern boundary of the Antarctic Circumpolar Current (ACC), the Mertz Glacier tongue and palaeo ice-sheet drainage path (white arrow) extending to the continental shelf edge, and the offshore slope and rise canyon system. Arrows represent downslope currents flowing into the Jussieu channel (black) with non-erosive overbank flow (yellow) towards the U1361 site on the continental rise. Black lines represent seismic reflection profile tracks represented inSupplementary Fig. 1a, b.

The studied interval consists of eighteen sedimentary cycles spanning an age range of 4.3–2.2 Ma, and comprising alternating terrigenous massive to laminated muds and diatom-rich/bearing silty-mud units (cycles 1–18 Fig. 2 and Supplementary Figs 2 and 3). In places, the muds contain packages of well-defined laminae and are consistent with established models of non-erosive overbank hemipelagic deposition onto a channel levee setting via turbidites on the lowermost Antarctic continental rise¹¹. Steeply dipping, seaward prograding wedge sediments are evident in seismic reflection profiles across the continental shelf and extend onto the upper continental slope above seismic unconformity WL-U8 (～4.2 Ma; Supplementary Fig. 1; refs 12, 13, 14). The geometry of these strata is characteristic of grounding zone deposition by repetitive advances of a marine-based ice sheet to the shelf edge during glacial periods¹³. Sediment overloading near the shelf break at submarine canyon heads, in turn triggers turbidity currents down slope channels and leads to overbank deposition at the core site. Low-density turbidity currents in overbank ‘distal’ channel levee environments on the Antarctic continental rise are typically non-erosive¹¹ (Supplementary Information). Thus, turbidite units are associated with periods of glacial advance to the Wilkes Land continental shelf edge, whereas bioturbated, diatom-rich/bearing facies represent warm interglacial periods of relatively ice-free ocean and increased primary productivity when the grounding line had migrated landward away from the shelf edge. Increased productivity during interglacial warm climates may be associated with enhanced upwelling of nutrient-rich Circumpolar Deep Water (CDW; ref. 15) (Supplementary Information), which has been linked to southward expansion of the westerly wind field in response to a reduced pole–equator temperature gradient¹⁶. Presently, this relatively warm nutrient-rich CDW upwells to the surface north of the Southern Boundary Front of the Antarctic Circumpolar Current and is marked by areas of enhanced productivity north of Site U1361 (Fig. 1).

a, XRF-based Ba/Al. b–f, IBRD MAR (b) correlated with January insolation and total integrated summer energy (melt threshold [t] = 400 GJ m^？2) (c), mean annual insolation and total integrated summer energy (melt threshold [t] = 250 GJ m^？2) (d), eccentricity³⁹ (e), and the stacked benthic δ¹⁸O record¹ (f). Also shown are lithofacies, lithological cycles (transitional lithologies are represented by both symbols) and magnetic polarity stratigraphy¹². Grey ellipse denotes alignment between a 1.2-Ma node in obliquity-modulated mean annual insolation and minimum in 400-kyr eccentricity and corresponds with MIS M2, a 1‰ glacial δ¹⁸O excursion culminating with MIS M2 (arrow). g, Atmospheric CO2 reconstructions based on boron isotopes and alkenones^31, 32, 43.

We have developed a high-resolution record (～3–4-kyr sample spacing) of iceberg-rafted debris (IBRD) mass accumulation rates (MAR) for the U1361 core (Supplementary Table 1 and Fig. 2). Our age model is based on a magnetic polarity stratigraphy constrained by biostratigraphy¹² and is not orbitally tuned (Supplementary Fig. 5). Thus, the age model for our primary spectral analyses (Fig. 3) assumes constant long-term sedimentation rates between polarity reversal tie points (Supplementary Information).

a, Normalized power for the untuned U1361 IBRD-MAR data¹². b, Normalized power for the LR04 δ¹⁸O stack¹. c, Amplitude for the eccentricity-tilt-precession (ETP) solution³⁹. Power has been normalized such that the maximum power in each 500-kyr window is unity. d–g, 2π MTM power spectral estimates for the U1361 IBRD-MAR data using the untuned age model¹² (d,e) and for the δ¹⁸O stack¹ (f,g), with red noise confidence levels using three different approaches (Supplementary Information). In a–c, red indicates large values and blue indicates low values.

In general, the highest intensity of IBRD occurs during transitions from glacial terrigenous mud facies to interglacial diatom-rich/bearing muds up-core until ～47 m below sea floor (mbsf), with most IBRD peaks immediately preceding diatom-rich facies and Ba/Al peaks (Fig. 2). Isotopic Nd and Sr provenance indicators suggest that the terrigenous components in these diatom-rich/bearing muds are associated with periods of glacial retreat of the ice margin back into the Wilkes Land Subglacial Basin during the Early Pliocene¹⁷. The Antarctic ice sheet loses the majority of its mass via iceberg calving and sub-ice-shelf melting¹⁸, which are intimately linked because enhanced ice-shelf melt leads to reduced buttressing, that in turn acts to enhance the flow of glacial tributaries, ultimately enhancing iceberg discharge¹⁹. Ice-sheet reconstructions based on geologic data²⁰ and modelling⁸ suggest that the above processes dominated Antarctic mass loss even during the warmest climates of the past 5 Myr. Thus, we interpret the maxima in IBRD MAR to be the consequence of accelerated calving during glacial retreat from marine terminating outlet glaciers along the Wilkes Land coastline, as well as a contribution from EAIS outlet glaciers entering the western Ross Sea (Supplementary Information). This interpretation is consistent with models and palaeodata that imply the most rapid mass loss of the EAIS margin during the last glacial termination occurred between 17 and 7 kyr ago, and although this retreat was initiated by oceanic warming it was associated with enhanced iceberg discharge^21, 22.

Spectral analysis of the untuned IBRD-MAR time series exhibits a dominant period of ～40 kyr between～4.2 and 3.5 Ma that transitions to strong and significant variance at ～100-kyr and ～20-kyr periods after ～3.5 Ma, with a corresponding decrease in power of the ～40-kyr cycle (Figs 3 and 4). On the basis of the strong orbital relationship exhibited in frequency spectra of the untuned IBRD-MAR time series, the palaeomagnetic age model, and the apparent near-continuous long-term sedimentation (Supplementary Fig. 5), we established a nominal one-to-one correlation between cycles in ice margin variability expressed by our IBRD data and orbitally paced climatic time series (Fig. 2). Between 4.3 and 3.4 Ma, a strong visual correlation can be observed between 40-kyr cycles in IBRD, mean annual insolation and the benthic δ¹⁸O global ice-volume record¹, whereas between 3.3 and 2.2 Ma, IBRD cycles correlate with the ～20-kyr cycles of summer insolation at 65° S (Fig. 2). We acknowledge our visual correlations, although constrained by seven precisely dated palaeomagnetic reversals, do not represent a unique solution but are consistent with the variance in orbital frequencies implied by our spectral estimations.

Cumulative power in the precession (0.04–0.06 cycles kyr^？1), obliquity (0.02–0.035 cycles kyr^？1), and eccentricity (0–0.015 cycles kyr^？1) bands for untuned IBRD-MAR data (a), astronomically tuned IBRD-MAR data (using tie lines in Fig. 2) (b), and LR04 δ¹⁸O stack¹ (c). These results indicate a generally dominant obliquity signal >3,500 kyr in the IBRD-MAR record, with an increase in eccentricity power in all three time series <3,500 kyr attributed to the clipping of precession, which transfers precession power to its eccentricity modulator. Obliquity dominates again after ～2.8 Ma in the LR04 δ¹⁸O stack, but not in the IBRD-MAR data.

Using these correlative relationships (Fig. 2) and our evolutionary spectral analyses on untuned IBRD data (Fig. 3), we then explore the role of longer-period orbital influences on the pattern of iceberg calving. The top of the ～40-kyr-dominated interval is marked by an ～200-kyr-long condensed section between ～3.5 and 3.3 Ma (Supplementary Fig. 5) and corresponds to a long-term + 1‰ glacial δ¹⁸O excursion punctuated by smaller interglacials spanning Marine Isotope Stages (MIS) MG5 and M2 (indicated by arrow in Fig. 2). The stratigraphic condensation, or possible erosion, at Site U1361 is associated with this glacial excursion. Indeed, this glacial excursion has also been associated with southern high-latitude climate cooling and the re-establishment of grounded ice on the middle to outer continental shelf in the Ross Sea following an ～200-kyr period of warm open ocean conditions^7, 9. Previous studies of older Oligocene and Miocene δ¹⁸O glacial excursions have proposed a relationship between intervals of increased glacial amplitude in the δ¹⁸O record with a coincidence of 1.2-Myr nodes in obliquity and 400-kyr minima in long-period eccentricity^23, 24. This orbital configuration, which favours extended periods of cold summers and low seasonality, is considered optimal for Antarctic ice-sheet expansion, and occurs at ～3.3 Ma (Fig. 3c)—the time of the transition from obliquity to precession dominance in the IBRD-MAR time series from U1361 (Figs 2 and 3a).

The observed ～20-kyr-duration IBRD cycles correlate with summer insolation calculated for 65° S for the interval of the core between 3.3 and 2.2 Ma (Figs 2 and 3a) and are embedded within 100-kyr-duration IBRD cycles, with an additional spectral peak at ～400 kyr nearing the 90% confidence level (Fig. 3d). Broad peaks of IBRD maxima are associated with transitions between laminated mudstones to diatom-rich/bearing muds (Fig. 2). A pronounced decrease in the amplitude of ～20-kyr IBRD peaks, and a change to lithofacies associated with non-erosive low-energy bottom currents, is observed at the core site above ～2.5 Ma, broadly coincident with southern high-latitude cooling⁹ and the onset of major Northern Hemisphere glaciations²⁵. We attribute the progressive reduction in calving intensity to cooling and a relative stabilization of the EAIS ice margin (discussed below). Homogenization of the turbidite sediments during glacial maxima by enhanced bioturbation and bottom current activity is observed and probably reflects increased Antarctic sea ice and polynya-style mixing at this time producing oxygenated high-salinity shelf water⁹ that is transferred downslope over Site U1361 to form Antarctic Bottom Water²⁶ (Supplementary Information).

In summary, our correlations of IBRD variations with the benthic δ¹⁸O stack and orbital time series identify up to seventeen ～40-kyr-duration cycles (orange dashed lines in Fig. 2) within six major lithological cycles (cycles 13–18, Fig. 2) during the Early Pliocene (4.3–3.5 Ma). This is followed by forty ～20-kyr-duration cycles (blue dashed lines in Fig. 2), modulated by eccentricity, that occur within twelve lithological cycles (cycles 1–12 Fig. 2).

Although the ANDRILL (AND-1B) marine sediment core from the Ross Sea region provided the first direct evidence that advance and retreat of the WAIS margin across the continental shelf was paced by obliquity during the Pliocene before ～3 Ma (ref. 7), subglacial erosion surfaces in the AND-1B core associated with ice advance have raised the possibility of missing cycles, particularly after 3.1 Ma—for example, the response of WAIS to orbital forcing is more ambiguous after this time. The U1361 and AND-1B records together confirm the dynamic response of both the WAIS and EAIS marine margins to obliquity forcing during the warm Pliocene before the onset of southern high-latitude cooling at 3.3 Ma (ref. 9).

Geologic records^7, 9, 17 and model simulations⁸ of past warm climates highlight the sensitivity of the marine-based portions of the Antarctic ice sheets to ocean warming. However, the mechanism by which the coastal ocean warms and destabilizes marine grounding lines, particularly in response to obliquity forcing, remains elusive. It has been proposed that changes in the intensity and meridional distribution of mean annual insolation controlled by obliquity may have a profound influence on the position and strength of the Southern Hemisphere zonal westerly winds, with implications for marine ice-sheet instability⁷. Indeed, an aerosolic dust record from the Southern Ocean is dominated by ～40-kyr cycles before ～0.8 Ma (ref. 27). Moreover, before ～3.3 Ma the southward expansion of the westerly wind field over the Antarctic circumpolar convergence zone under a reduced meridional temperature gradient has been associated with a reduced sea-ice field⁹. This may have caused the upwelling of warm, CO2-rich CDW onto the continental shelf, with consequences for the migration of marine grounding lines^16, 28, 29.

The reduction in obliquity power revealed by our data after ～3.5 Ma, and comparably stronger influence of precession and eccentricity (Figs 3a and 4), is interpreted to reflect a declining sensitivity of the EAIS to oceanic forcing as the southern high latitudes cooled. Both model and geologic reconstructions imply that past Antarctic ice-sheet expansion is closely linked with development of the sea-ice field³⁰potentially resulting in northward migration of westerly winds and Southern Ocean fronts⁹. Furthermore, sea-ice expansion after 3.5 Ma (ref. 9) potentially restricted upwelling and ventilation of warm CO2-rich CDW at the Antarctic margin, acting to further enhance climate cooling. Under such a scenario, during the warmer climate state of the Early to mid-Pliocene with higher atmospheric CO2 concentration^31, 32, the increased duration and intensity of austral summer surface warming produced a pattern approaching mean annual insolation regulated by obliquity (Fig. 2d), rather than precession that cancels out over the course of a seasonal cycle⁵. Late Pliocene cooling raised the melt threshold such that the duration of the melt season was restricted to times of austral summer insolation maxima controlled by precession (Fig. 2c), and extending the winter sea-ice growth season³³. We propose this decrease in CO2 radiative forcing resulted in extensive sea-ice cover extending into much of the summer season, limiting the influence of upwelling CDW on marine grounding line stability.

The significant variance at ～20-kyr precession and coincident increase in the ～100-kyr eccentricity frequency bands between 3.5 and 2.5 Ma is intriguing (Figs 3a, d and 4), and is similar to the pattern of orbital response expressed in ～100-kyr-duration glacial–interglacial cycles in the benthic δ¹⁸O and greenhouse gas ice core records of the past 600 kyrs (refs 34, 35). Thus the relative dominance of eccentricity in EAIS variability at the Wilkes Land marine margin after 3.5 Ma may reflect a threshold response of the ice sheet to orbital forcing under a colder climate regime, as has been proposed to explain the relative dominance of ～100-kyr cycles in global Late Pleistocene ice volume, albeit with different mass balance forcings^33, 36, 37. Here we propose that the precession response of perennial sea ice after 3.3 Ma was modulated by eccentricity (Fig. 4), such that the areal extent of sea ice was significantly reduced every ～100 kyr, enhancing the influence of oceanic warming at those times. The observation of strong ～100-kyr and 400-kyr power is consistent with a nonlinear clipped climate response³⁴, whereby precession power is transferred to its eccentricity modulator (Fig. 3d).

Our data are also consistent with a reduction of obliquity variance in the benthic δ¹⁸O (refs 1, 38) at～3.5 Ma (Fig. 3b), coincident with an obliquity node in the astronomical solution³⁹. The gradual re-emergence of a 40-kyr signal at ～3 Ma in benthic δ¹⁸O record³⁸ (Figs 3b and 4c) most probably reflects a similar re-emergence of strong obliquity forcing in the orbital records and possibly also a direct response of the developing Northern Hemisphere continental ice sheets to obliquity forcing. However, it remains possible that precession-driven, anti-phase oscillations in both hemispheric ice-volume histories may have cancelled out in globally integrated proxy records after 3–2.8 Ma (for example, ref. 4).

The strong obliquity signal seen in our IBRD record also suggests that the intensity of summer insolation was a not a direct control on the mass balance of the EAIS before 3.5 Ma. The geometry of strata on the Wilkes Land continental shelf indicate that the EAIS expanded towards the continental shelf edge during glacial maxima in the Pliocene¹², indicating most Antarctic ice-volume variance at this time was growth and retreat of the marine-based ice sheets (Supplementary Fig. 1). During the Early Pliocene, when the sea-ice field was reduced and the Wilkes Land margin of the EAIS was in more direct contact with oceanic influences, iceberg calving occurred more regularly within both facies. Based on the significant decrease in IBRD after 2.5 Ma (Fig. 2) and Southern Ocean records inferring decreased SSTs (refs 9, 40), we also infer the marine margins of the EAIS became less sensitive to ocean-induced melting compared to the WAIS (ref. 29). After this time, East Antarctic ice volume probably fluctuated with a magnitude similar to that of Late Pleistocene cycles, with a minimal contribution to the δ¹⁸O signal^8, 21.

Notwithstanding this relative stability of the marine-based sectors of the EAIS, model results suggest that～20-kyr-duration fluctuations in the total Antarctic ice volume (that is, including the WAIS) may have contributed up to 15 m of the total amplitude during Late Pleistocene glacial–interglacial cycles⁸. Given the δ¹⁸O composition of Antarctic ice, this magnitude of variability could have offset a larger out-of-phase precessional change in Northern Hemisphere ice volume (for example, 20–40 m), resulting in an enhanced obliquity signal in globally integrated sea-level and ice-volume proxy records after 3.0–2.8 Ma (for example, ref. 4, Figs 3b and 4). Alternatively, direct obliquity forcing of the Northern Hemisphere ice sheet is supported by proxy evidence including ice-rafted debris records²⁵ and a recent dust flux record⁴¹, suggesting that Northern Hemisphere ice-sheet variability (marine-based margins) and climate primarily responded to obliquity under a relatively warm Northern Hemisphere climate state. Northern Hemisphere cooling and ice-sheet growth across the mid-Pleistocene transition ～0.8 Ma has been implicated in a similar switch to precession- and eccentricity-paced glaciations, albeit by different mass balance forcing mechanisms⁴². In contrast, however, our results imply that Southern Ocean sea-ice feedbacks caused a fundamentally different response of the marine-based sectors of the EAIS under a cooler Late Pliocene/Early Pleistocene climate state, characterized by a dominance of precession-paced variability.

We conclude that, before 3.5 Ma, under a warm climate state, EAIS demonstrates high sensitivity on orbital timescales to a relatively small increase in atmospheric CO2 concentration and mean global surface temperature. With atmospheric CO2 concentrations and global surface temperatures projected to remain above 400 ppm and >+ 2 °C beyond 2100⁴³, our results have implications for the equilibrium response of the Antarctic ice sheets, and suggest that the marine margins of the EAIS, as well as the marine-based WAIS, may become increasingly susceptible to ocean warming, with the potential for widespread mass loss raising sea level by metres over the coming centuries to millennia.