Drill core evidence for ice-sheet growth

We reviewed 60 Southern Ocean and Antarctic margin sedimentary successions to reconstruct AIS evolution during the MSC. Most sites reviewed were not included because they either did not recover MSC-aged sediments or because they had insufficient chronostratigraphic control. We identified four drill sites that contain sediments indicative of ice-sheet expansion during the MSC.

Ocean Drilling Program (ODP) sites 1092, 1095 and 1165 (Fig. 1a) were all recovered in relatively deep water (1,000's of metres) and provide a sedimentary archive of Antarctic Circumpolar Current (ACC) evolution in response to AIS growth and decay. Antarctic margin core AND-1B contains a direct record of ice sheet advance and retreat across the continental shelf. Integrated ODP site 1361 recovered a continuous succession spanning the MSC offshore of the Wilkes Land margin in 3,466 m of water. However, no evidence exists at this site for a significant change in sedimentation or ocean circulation, which indicates that this site may not have been sensitive to ice-volume changes during the MSC15.

Figure 1: Maps of core locations. Locations of ODP sites 1092, 1095 and 1165 and the AND-1B succession discussed in this study (a). ODP sites are situated in deep water and, therefore, track the strength of deep-ocean circulation in response to ice-sheet growth. The AND-1B succession is on the continental shelf and was directly influenced by expanding ice sheets that eroded sediments from the shelf and resulted in glacial erosion surfaces in the succession. A continuous succession spanning the Messinian was recovered at site 1361, but there is no evidence of a significant change in sedimentation or ocean circulation at this site. (b) Location map of ODP sites 846 and 926 from which deep-sea δ18O records are derived. Full size image

Sediments at ODP site 1165 were recovered from a water depth of 3,537 m near the East Antarctic coastline at the boundary between the polar gyres and ACC and, critically, near the mouth of Lambert Glacier, which is the largest East AIS outlet glacier, and experienced a ca. 400-km grounding line migration during the latest Miocene16,17. Here, we focus on the interval between 50 and 90 m below seafloor (5–7 Ma), which contains a previously poorly constrained unconformity18 (Fig. 2).

Figure 2: Summary of ODP site 1165 lithology and palaeomagnetic data. ODP site 1165 succession from Prydz Bay, Antarctica, located on the Antarctic Divergence, which comprises a series of gyres at the boundary between the Antarctic Circumpolar Current (ACC) and the Polar Current. The Polar Current and ACC are major surface current systems that extend into Antarctic deep water, and their variable interactions partly control sediment distribution across the continental margin. The ODP site 1165 record between 50 and 90 metres below seafloor (m.b.s.f.). (a) Correlation of the magnetostratigraphy (black (grey) data points are from Hole 1165B (Hole 1165C)) with the geomagnetic polarity timescale (ATNTS2012; ref. 19) as guided by five biostratigraphic constraints (markers adjacent to magnetic susceptibility data). (b) Glacial/interglacial cycles as expressed in sediment density (b, filtered and e, raw), and magnetic susceptibility (c, filtered and f raw) variations are correlated with (d) orbital obliquity66. Spectral analyses of (g) magnetic susceptibility and (h) sediment density indicate statistically significant and (i) coherent cycles in both records. An unconformity at a depth of 67 m is recognized from the abrupt polarity change and from missing obliquity cycles and has an estimated duration of 890 kyr with an upper and lower age of 5.61 and 6.50 Ma, respectively. Solid black lines indicate unambiguous correlations of magnetic polarity intervals with the GPTS and dashed lines indicate correlation of magnetic susceptibility and density data with the obliquity record. Full size image

We developed a precise chronology from a revised magnetobiostratigraphy18, which facilitated correlation with the ATNTS2012 geomagnetic polarity timescale19. Palaeomagnetic analyses reveal a reliable magnetization that results in a well-defined magnetic polarity record with five polarity intervals. Rock magnetic analyses indicate that a mixture of single-domain and pseudo-single-domain magnetite is responsible for the magnetization, and that diagenetic alteration of the magnetite has not occurred. We use biostratigraphic constraints from shipboard observations of first and last appearance datums (FADs and LADs) of diatoms and radiolaria with updated published calibrations to correlate the magnetostratigraphy with ATNTS2012. Marine diatom (MD) datum MD1 is the LAD of Nitzschia donahuensis, which occurs at 56.45 m below seafloor (m.b.s.f.), is calibrated at 5.8 Ma (ref. 20), which results in an unambiguous correlation of the R1-N1 reversal with the C3n.4n-C3r reversal. MD2 is the LAD of Nitzschia miocenica, which occurs at 63.59 m.b.s.f. and is calibrated at between 6.0 and 6.2 Ma (ref. 18). Top Amphymenium challengerae is the LAD of A. challengerae, which occurs at 65.45 m.b.s.f. and is calibrated between 6.1 Ma (ref. 18) and 6.2 Ma (ref. 20). Bottom A. challengerae is the FAD of A. challengerae, which occurs at 72.87 m.b.s.f. and is calibrated at between 6.65 Ma (ref. 18) and 6.8 Ma (ref. 20). MD3 is the FAD of Thalassiosira miocenica, which occurs at 73.00 m.b.s.f. and has multiple calibrations of 6.4 Ma (ref. 18), 5.91 Ma (ref. 20) and 6.25–8.3 Ma (ref. 21). For MD3, we use the original calibration of 6.4 Ma (ref. 18), because it follows the downward progression of FADs and LADs most closely. This results in a correlation of the N3-R2-N2 sequence with C3Bn-C3Ar-C3An.2n.

Below the recorded part of Chron C3r, at 67.02 m.b.s.f., a sharp magnetic polarity transition marks a disconformity (Fig. 2a). We determine precisely the upper and lower ages of the disconformity from time-series analysis and bandpass filtering of orbitally paced sediment density and susceptibility cycles (Fig. 2e,f). Sediment density and susceptibility cycles are driven by alternations of biogenic-rich versus terrigenous dominated sediment that are inferred to correspond to alternating warm and productive periods versus colder periods with greater ice volume22. Spectral analyses (Fig. 2g,h,i) reveal a dominant wavelength of ∼0.77 cycles per metre (above the 95% confidence limit for magnetic susceptibility and above 90% for sediment density). Spectral power is greater in the magnetic susceptibility data probably because it is more sensitive to the terrigenous-to-biogenic ratio in sediments; therefore, magnetic susceptibility data were bandpass filtered to isolate the orbital signal (Fig. 2b,c). In total, 27 obliquity-paced glacial–interglacial cycles were identified and correlated one for one with the orbital timescale. The correlation indicates a break in deposition or removal of 890 kyr of sediment between 5.61 and 6.5 Ma. The disconformity probably resulted from current winnowing and non-deposition by erosive bottom currents during AIS expansion23. The base of the unconformity (6.5 Ma) marks the downward limit of erosion, not the onset of erosion or ice expansion. Our age model indicates that sedimentation resumed at 5.61 Ma, coincident with the warm, interglacial stage TG15 following a reduction in ocean current speed and ice volume.

ODP site 1092 comprises a succession dominated by biogenic carbonate (Fig. 3). We refined the age of a poorly constrained unconformity24,25 using shipboard diatom abundance counts26, improved diatom bioevent calibrations21 and a reassessment of magnetostratigraphic data25. We focus on an interval above the C3An.1n-C3r reversal boundary at 74 metres composite depth (m.c.d.). We used updated LAD and FAD calibrations to assign the normal polarity interval above 70.5 m to chron C3n.4n (Fig. 3d) and to correlate the record with the ATNTS2012 timescale19. This correlation is constrained by the unambiguous FADs of Fragilariopsis lacrima (69.81 m.c.d.) and Thalassiosira inura (68.61 m.c.d.), which have age calibrations of 4.73 and 4.74 Ma, respectively. We use a combination of geomagnetic reversals below the unconformity25 and bioevents higher in the succession26 and obtained a smooth average sedimentation rate that indicates the presence of one or several unconformities between 72.5 and 70.5 m.c.d. (Fig. 3). The palaeomagnetic inclination data in this interval also indicate at least two intervals with inconsistent palaeomagnetic directions and decreased carbonate content that likely indicate the presence of an unconformity caused by strong, corrosive bottom currents. Our best estimate for the interval that contains these unconformities suggests a basal age of between 5.9 and 5.8 Ma and an upper age of 5 Ma. Additional, unrecognized unconformities in the interval between 70 and 73 m could be present, because the site is located on a bathymetric rise and would have been exposed to erosive currents.

Figure 3: Summary of ODP site 1092 lithology with carbonate content and palaeomagnetic data. ODP site 1092 succession from the South Atlantic Ocean. Sediments comprise (a) fine-grained nannofossil ooze and sedimentological analyses indicate an abrupt decrease in (c) carbonate content between ca. 72.5 and 71.5 m. (d) ChRM inclinations indicate that these intervals have unstable magnetization, and correlation with the ATNTS2012 timescale19 indicates that (b) one or more unconformities must occur in the succession. The youngest possible upper age for the uppermost unconformity in our age model is 5.8 Ma. Full size image

Antarctic Peninsula ODP site 1095 contains a continuous succession recovered from a depth of 3,840 m that spans the MSC with a reliable, well-defined chronology27 and no evidence for a significant unconformity (Fig. 4). However, sedimentological analyses (Fig. 4a,c) and anisotropy of magnetic susceptibility (AMS) data shed light on changes in ACC strength28. AMS provides a measure of fabric strength in the sediment (P′, Fig. 4b) where strong fabrics indicate well-aligned grains and an inferred increased current strength. Sedimentological analyses and core logs indicate that the sediments comprise repetitively bedded, weakly laminated silty clays with prominent, graded silt laminae that are interpreted to have been deposited under stronger current regimes. Silty laminae increase in number up-core, which indicates stronger circulation29; sedimentary and AMS analyses indicate increased delivery of terrigenous material and an overall stronger grain alignment at ∼6.3 and ∼5.6 Ma (ref. 28), which indicate greater current speed likely in response to greater ice volume. Circulation was strongest between 6 and 5.6 Ma with weakening currents after 5.6 Ma and minimum circulation at ca. 5.3 Ma coincident with reduced grain size29 and increased biogenic productivity that is interpreted to indicate reduced ice volume28.

Figure 4: Summary of ODP site 1095 lithology with AMS palaeocurrent proxy and palaeomagnetic data. ODP site 1095 succession from offshore of the Antarctic Peninsula. Sediments comprise (a) weakly laminated silty clays with thin-silt laminae throughout the succession. Sedimentological analyses (c) indicate that terrigenous sediments are dominant until ca. 5.5 Ma after which biogenic sediments become more dominant and that coarse-grained laminae reach peak numbers between ca. 6.3 and 5.5 Ma, which indicate a strong erosive current system. (b) AMS (degree of magnetic anisotropy, P′) analyses are indicative of grain alignment in response to ancient current strength with strong grain alignment until ca. 5.5 Ma, which indicate strong, deep currents. AMS data, a decrease in grain size and fewer silt laminae indicate weaker circulation after 5.5 Ma. Full size image

The AND-1B succession (Fig. 5b) is the most ice-sheet-proximal record recovered from the Antarctic margin and contains a record of AIS advance and retreat history that spans the late Miocene to Holocene30. The succession comprises massive to stratified diamictites that represent grounded ice or ice-proximal conditions, muddy units that represent ice distal conditions and diatomite intervals that represent deposition in open-marine conditions, in some cases with minimal sea ice31,32,33. Transitions from warm, low-ice-volume conditions to cool, high-ice-volume conditions are typically separated by surfaces where the ice sheet advanced over the drill site and eroded sediment. The succession contains a glacial erosion surface ‘U8’ at 596.35 m.b.s.f., which has an estimated age range between 5.90 and 5.60 Ma (Fig. 5b)34. Above unconformity U8, a switch from glacially dominated conditions to open-marine conditions is recognized along with the appearance of Shinodiscus tetraoestrupii diatoms that are indicative of warm surface conditions (7–10 °C)32 and likely much lower ice volume.

Figure 5: Regional event summary from drill cores with modelled sea-level and Antarctic ice-volume benthic δ18O. Regional event summary compiled using (a) Southern Ocean ODP successions, (b) the Antarctic Margin AND-1B record and ice-volume/ocean current interpretations derived from sedimentary successions. (c) Interpretation of Antarctic climate events. (d) GIA only modelled r.s.l. at ODP site 1165 (red), at the Strait of Gibraltar (blue) and global eustatic sea level (black). (e,f) The δ18O record38,39, marine isotope stages and isotope-to-ice volume calibration used for sea-level modelling. Full size image

Reconstructing AIS history and sea level

AIS variations exerted a primary control on global sea level on short geological timescales from ∼34 Ma (ref. 35) until expansion of large northern hemisphere ice sheets after ∼2.7 Ma (ref. 36). The presence of large northern hemisphere ice sheets will have likely amplified global eustatic sea-level variations37. AIS history was largely inferred from benthic δ18O records until recovery of well-dated Antarctic margin successions30. These ice-sheet-proximal geological records and subsequent modelling studies reveal that ice sheets grew slowly and retreated rapidly during the late Neogene17,30. To determine realistically how the sea level evolved, we conducted numerical GIA simulations by means of the sea-level equation (SLE)13,14. Solving the SLE requires a solid Earth model for crustal and gravitational response13,14 and an ice-sheet chronology as a forcing function, which we generated by scaling present-day AIS thickness using a δ18O-based ice volume curve38.

AIS volume and thickness reconstructions

For the time interval under consideration (6–5 Ma), Antarctic ice-sheet thickness variations are unavailable from continuous global circulation model studies. Therefore, we reconstructed AIS volume (Fig. 5e) using benthic δ18O records from ODP site 926 (Ceara Rise, 3,598 m water depth, Fig. 1b, ref. 39) between 7 and 6.138 Ma and from ODP site 846 (south of the Galapagos Islands, 3,296 m water depth, Fig. 1b, ref. 38) for the interval between 6.137 and 5 Ma (Fig. 5). We estimate AIS volume using the most conservative approach possible by assuming the modern-day AIS isotopic weight of ∼−53.2‰, total melted water volume of 22.279 × 106 km3 and a 1,335 × 106 km3 global ocean volume40. Therefore, a total loss of the modern AIS would result in a ∼ 0.91‰ inflection of the deep-sea δ18O record. Accordingly, we estimate a ∼58% reduction of Antarctic ice volume during oxygen isotope stage TG5 from a ∼0.53‰ inflection in the δ18O record. However, to test for different glacial to interglacial temperature contributions to the deep-sea δ18O record, we developed two additional ice-volume records that correct for a ca. 20% and 30% temperature contribution (2 and 3 °C), respectively,41 in agreement with Mg/Ca records that indicate a ca. 2 °C temperature variation during the Miocene42. We tested several other methods to convert benthic δ18O to ice volume, including using isotopically heavier ice and the relatively well-understood Pleistocene isotope to sea-level calibration of 0.01‰ m−1. The calculated ice-volume fluctuations are amplified, which results in too many intervals with negative ice volume (15% of the record). The rationale for converting with isotopically heavier ice is plausible because it is likely that Miocene ice sheets were warmer and, therefore, isotopically heavier than the modern AIS (some studies indicate that the Oligocene AIS likely had an isotopic weight of ∼−35‰ (refs 43, 44)). However, we used a more conservative modern ∼−47‰ isotopic weight, which results in more reasonable ice-volume changes. We also tested the Pleistocene to early Pliocene benthic δ18O to sea-level calibration of 0.01‰ m−1 (refs 45, 46), which removes the temperature contribution from the record. We computed ice volumes from the sea-level curve by assuming a modern Antarctic sea-level contribution of 56.6 m (ref. 40), which resulted in unrealistically large ice-volume changes. Including a potential sea-level contribution of 7.3 m from the Greenland Ice Sheet (GIS) produced negligible changes to the ice-volume record. Changes in GIS size and other sources of unidentified northern hemisphere ice would have contributed to variations in the δ18O record; however, the modern GIS comprises only ∼7% of the global ice volume, therefore, the majority of δ18O variations can be attributed to AIS changes. We converted the high-resolution, orbitally tuned δ18O record from ODP site 846 (ref. 38) to ice volume and, accordingly, scaled present-day AIS thickness over time and used this for the GIA simulation.

GIA and sea-level modelling

We performed a numerical GIA simulation to reconstruct the impact of AIS volume changes on r.s.l. at the Antarctic margin and at Gibraltar. Any ice-sheet fluctuation results in local r.s.l. changes (that is, vertical geoid variations with respect to the deforming solid Earth surface) that stem from a complex interplay between gravitational, rotational and solid Earth deformations in response to redistribution of surface ice- and water loads47. Spatial variability of r.s.l. change depends on the distance from the changing ice sheets and on the shape and size of ocean basins48. Second, given the viscous behaviour of the solid Earth on geological timescales, ice-induced r.s.l. change varies in time as a function of mantle viscosity. Local r.s.l. change can, therefore, be significantly different from the globally uniform glacio-eustatic sea-level change. Therefore, correlating r.s.l. change at a given location to a specific ice-sheet volume variation requires precise spatio-temporal discretization of the latter (that is, how much ice thickness changed, and where and when this occurred) and a rheological model for the solid Earth response. These two main factors enter into the gravitationally self-consistent SLE whose solution provides the global r.s.l. change. To determine local r.s.l. at a given point, we solved the SLE using the pseudo-spectral method including consistent time-dependent coastline evolution and rotational feedback to meltwater redistribution49,50,51,52. We employ a radially stratified, spherically symmetric and rotating Earth model characterized by an upper elastic, 100-km-thick lithosphere, a three-layer Maxwell viscoelastic mantle with an inner mantle viscosity of 5 × 1021 Pa s, an outer mantle viscosity of 0.5 × 1021 Pa s and an inner inviscid core. Our initial AIS is smaller than at present, so we also modified the present-day initial global topography model ETOPO1 (ref. 53) by melting the excess mass from the present-day AIS and allowing a 50-kyr viscoelastic relaxation via the SLE. We also decreased the depth of the Gibraltar Strait to 30 m below mean sea level54,55, which is shallower than the modern, deep channel that was eroded at the end of the MSC9.