Dated radar-stratigraphy between Dome A and South Pole, East Antarctica: old ice potential and ice sheet history

Abstract

An array of information about the Antarctic ice sheet can be extracted from ice-sheet internal architecture imaged by airborne ice-penetrating radar surveys. We identify, trace and date three key internal reflection horizons (IRHs) across multiple radar surveys from South Pole to Dome A, East Antarctica. Ages of ~38 ± 2.2, ~90 ± 3.6 and ~162 ± 6.7 ka are assigned to the three IRHs, with verification of the upper IRH age from the South Pole ice core. The resultant englacial stratigraphy is used to identify the locations of the oldest ice, specifically in the upper Byrd Glacier catchment and the Gamburtsev Subglacial Mountains. The distinct glaciological conditions of the Gamburtsev Mountains, including slower ice flow, low geothermal heat flux and frozen base, make it the more likely to host the oldest ice. We also observe a distinct drawdown of IRH geometry around South Pole, indicative of melting from enhanced geothermal heat flux or the removal of deeper, older ice under a previous faster ice flow regime. Our traced IRHs underpin the wider objective to develop a continental-scale database of IRHs which will constrain and validate future ice-sheet modelling and the history of the Antarctic ice sheet.

1. Introduction

Widespread radio-echo sounding (RES) data have revealed an extensive archive of the dynamic ice flow and past climate of Antarctica, providing information complementary to ice-core analyses (Siegert and others, 1998; Cavitte and others, 2016; Winter and others, 2019; Ashmore and others, 2020; Bodart and others, 2021). Ice-sheet englacial stratigraphy obtained from internal reflection horizons (IRHs) in RES data (Bingham and Siegert, 2007) enables us to expand our understanding of past accumulation rates (e.g. Leysinger Vieli and others, 2011; Bodart and others, 2022) and ice-flow changes over thousands of years (e.g. Bingham and others, 2007; Winter and others, 2015; Siegert and others, 2019), thousands of kilometres away from rare, point-location ice cores. Connecting ice cores via continuous IRHs can: (i) synchronise ice-core age scales; (ii) reduce uncertainties associated with current ice-core age-depth sequences through comparisons; and (iii) inform the selection of future ice-core sites (MacGregor and others, 2015).

The presence of IRHs are often a result of conductivity variations as a result of differing ice chemistry (most of the ice column); density changes in the ice, typically associated with impurities within firn layers (typically in the upper tens of m of the ice column); or changes in the ice fabric (most commonly in the deepest ice) (Bingham and Siegert, 2007; Cavitte and others, 2016; Holschuh and others, 2018). Continuous IRHs are generally considered isochronal (Whillans, 1976; Siegert and others, 1998; Siegert, 1999), and therefore primarily reflect the burial and advection of palaeo-ice-sheet surfaces (Ashmore and others, 2020). Imaged IRH architecture therefore provides a record of surface mass balance, ice flow and basal melt, while setting critical age tracers in the ice sheet (Sutter and others, 2021).

A continental-scale database of traced and dated IRHs in Antarctica, similar to that produced for the Greenland Ice Sheet (MacGregor and others, 2015), is required to constrain and validate ice-sheet models. Information obtained through englacial architecture can be used to reconstruct the evolution of the East Antarctic Ice Sheet (EAIS). However, to date, englacial stratigraphy with any degree of dating control, has only been obtained over finite areas of the ice sheet (Leysinger Vieli and others, 2011; Cavitte and others, 2016; Winter and others, 2019; Wang and others, 2023). IRHs have previously been traced throughout a 200 km radius around Dome C, East Antarctica; and at this site intersections of the englacial stratigraphy with the EPICA Dome C Ice Core (hereafter EDC) age-depth profile have allowed the construction of a 3-D age-depth profile spanning the last two glacial cycles (Leysinger Vieli and others, 2011; Cavitte and others, 2016; Winter and others, 2019). Some of these IRHs were subsequently traced along flightlines connecting Dome C, Vostok and Dome A (Winter and others, 2019). Englacial architecture between South Pole and the southern flanks of Dome A (Fig. 1) has, however, received little attention. Consequently, we know little about the age-depth relationship or ice-sheet history of this slow flowing region of the EAIS. RES surveys in this region have been collected and interrogated to determine past and future ice-core targets (Brook and others, 2006; Parrenin and others, 2017; Van Liefferinge and others, 2018); geothermal conditions beneath the ice (Jordan and others, 2018); ice-stream onset-zone boundary conditions (Winter and others, 2018); and the occurrence and form of basal ice units (Bell and others, 2011; Wrona and others, 2018). While all of these studies have provided important insights and advances, there remains a clear need for a more holistic and widely spatially-extensive approach to characterising the englacial stratigraphy of the region between South Pole and Dome A, taking advantage of airborne RES surveys acquired over the last two decades that collectively link the two regions.

Figure 1. (a) Regional map of data analysed in this study. Points M1 and M2 denote sites we used for Dansgaard/Jensen age-depth modelling. A and A’ define the extent of the radargram in b. Background is ice thickness from BedMachine v2 (Morlighem and others, 2020). Regional ice divides (Zwally and others, 2012) are noted in light grey and highlighted ice features include: Recovery Glacier (RG); Academy Glacier (AG); and Lambert Glacier (LG). (b) An example of radar data (in chirp mode) from the AGAP survey (flightline A10b A to A’ in a) in which the three distinct IRHs (H1-H3) mapped in this study are marked with coloured symbols. (c) flightline A10b including H1-H3 picks.

Here, we undertake extensive tracing of IRHs between South Pole and Dome A. We focus on three distinct IRHs which traverse central East Antarctica, linking previously traced IRHs at Dome A (Cavitte and others, 2016; Winter and others, 2019) and Titan Dome (Beem and others, 2021), with South Pole Ice Core (SPICEcore) data (Winski and others, 2019). We provide age constraints for each of our traced IRHs through: (1) their intersections with dated englacial stratigraphy around Dome A linked to Vostok and EDC; (2) a direct intersection with SPICEcore; and (3) independent verification with a 1-D model, which provides insight into the appropriate parameterisation of such models elsewhere in East Antarctica where there may be no direct links to any ice-core age-depth profiles. We use our new regional englacial stratigraphy to provide insights towards locating sites of oldest ice on the flanks of Dome A, and on past ice dynamics and elevated geothermal heat flux near South Pole.

2. Data and methods

2.1 Radar datasets

The RES data utilised for this study were collected during two research campaigns: the Antarctic GAmburtsev Province (AGAP) survey conducted in the austral seasons 2007/08 and 2008/09 (Bell and others, 2011; Ferraccioli and others, 2011; Rose and others, 2013; Sanderson and others, 2023) and the PolarGap survey acquired in 2015/16 (Jordan and others, 2018; Winter and others, 2018; Paxman and others, 2019). AGAP as an entire project incorporated comprehensive surveying of the Gamburtsev Subglacial Mountains and Lambert Glacier/Rift – with, crucially for this paper, several connecting flight tracks to South Pole – and for logistical reasons was divided, largely regionally, into surveys conducted respectively from the British Antarctic Survey and Lamont-Doherty Earth Observatory aerogeophysical platforms. For this paper, we use AGAP data that were acquired by the British Antarctic Survey, focusing on the flightlines that explicitly connect Dome A to the South Pole (Fig. 1): and for the rest of this paper we use the term AGAP to refer only to this subset of the whole project's dataset. The RES data were acquired with the British Antarctic Survey's Polarimetric Radar INstrument (PASIN), which operated with a 150 MHz centre frequency, and used a 10 MHz chirp to sound deep into the ice, producing vertical sampling resolution of ~8.4 m. Chirp compression and incoherent 2-D Synthetic Aperture Radar (SAR) processing were applied to the data to enhance the along-track resolution and echo signal to noise. An additional incoherent averaging filter was applied to the chirp data (Corr and others, 2021). The PolarGAP survey used an updated radar, PASIN2, which while still operating with a 150 MHz centre frequency, used a wider-frequency (13 MHz) linear chirp, producing an improved vertical sampling resolution of 6.5 m. As for the AGAP data, chirp compression was applied, but the PolarGAP data were differently processed using a coherent averaging filter (unfocused SAR processing) with Doppler beam sharpening to enhance the signal to clutter ratio of the bed echo, improving visualisation (Ferraccioli and others, 2021). Together, these data sets provide multiple opportunities to link IRHs across intersecting RES lines, and to date IRHs using ages acquired from the South Pole Ice Core and other RES surveys.

2.2 Tracing internal reflecting horizons

To optimise the display and traceability of IRHs, we applied a natural-log filter and a 10-trace horizontal average to both the AGAP and PolarGAP radar data. Data were loaded into the freely available Opendtect Seismic Interpretation Software (https://www.dgbes.com/software/opendtect) for 3-D analysis and maximum-amplitude layer picking. IRH tracing was initiated on AGAP flightline A10b (Fig. 1b) because it contains clearly visible and continuous IRHs imaged over an area of almost stagnant ice flow (<5 ma⁻¹). We traced three particularly clear IRHs, which we name H1-H3, along this A10b control line. In addition to being distinct along flightline A10b, our initial reconnaissance of the dataset suggested that these IRHs occurred at similar depths in the ice column to IRHs traced by Winter and others (2019) between Dome A and Dome C. We then progressively extended the tracing of H1-H3 along intersecting flightlines from the AGAP and PolarGap surveys. IRHs H1 and H2 are easily identified and traced throughout A10b and almost all other intersecting flight lines. H3, being further down the ice column, was harder to trace in places, although we still traced it along a large number of flight lines, down to depths of ~2950 m. We note that many other IRHs are visible within the radar data, but most could not be traced across the study area (e.g. due to reflection bifurcation or convergence) and are therefore not included in this study.

In places where tracing was not possible, due to discontinuity or an absence of reflectors, H1-H3 tracing was extended across the IRHs that were distinctly brighter than the other IRHs in the column and their diagnostic stratigraphic signature (i.e. layer width and amplitude). All IRH two-way-travel (TWT) picks were converted to depth using an electromagnetic wave speed of 168.5 m μs⁻¹ and a spatially constant firn correction of + 12 m (supplementary material, section S1). Following Ashmore and others (2020) and Bodart and others (2021), a conservative uncertainty in IRH tracing depth records, arising from variations in electromagnetic wave speed, firn correction and radar range accuracy results in a conservative vertical uncertainty of ± 17 m for H1, ± 21 m to H2 and ± 27 m to H3 (supplementary material, section S2).

2.3. Applying age constraints to internal reflection horizons

2.3.1. Intersections with ice-core chronologies

We first applied age constraints to H1-H3 by analysing where they intersected flightlines with previously traced and dated englacial stratigraphy connecting to Vostok and EDC. This existing age-depth information was provided by Winter and others (2019) who traced IRHs in radar data collected between Dome A, Vostok and Dome C using a 150 MHz centre-frequency RES system (based on the MCoRDS system (Bell and others, 2011)), and based their chronology on Bazin and others (2013) ice-core chronology for Vostok and EDC. To determine the errors associated with dating IRHs using these previously traced IRHs we used a root-mean-squared analysis of the differences in depth at the crossover points (supplementary material, section S3) and refer to the closest ice core chronology to the crossover intersections.

Secondly, because our dataset traverses South Pole, we were able to provide an independent verification of the age of any of the IRHs occurring down to 1751 m depth, the deepest ice dated by SPICEcore (Winski and others, 2019). Our IRHs traced from the PolarGAP radar survey pass within 86 m of the SPICEcore. Following MacGregor and others (2015), we took the unweighted mean traced pick depth ± 250 m from the closest trace approaching the drill site and used this to assign an age from the chronology that closely matches with our IRH depth. Uncertainties associated with this method are discussed in supplementary material, section S4. An additional independent constraint on H1 was provided by a further intersection with an IRH traced by Beem and others (2021) using various surveys around South Pole and Titan Dome and dated to 37.6 ka at its own intersection with SPICEcore (see supplementary material, section S3) using a coherent 60 MHz centre-frequency radar ice sounder (Peters and others, 2005). We combined the age association from the intersection with Winter and others (2019) and the date obtained from the SPICEcore to achieve a final age association for H1 (supplementary material, section S5).

2.3.2. Age-depth modelling

We calculated an independent validation of IRH ages (Bodart and others, 2021) by applying the Dansgaard and Johnsen (1969) one dimensional vertical strain rate model to our RES data (supplementary material, section S6). We chose to apply the Dansgaard-Johnsen model in this case for its simplicity, allowing us to evaluate published accumulation rates and the impact of the basal shear level thickness on deep, older IRH ages. The model has previously been applied to calculate accumulation rates near ice divides as the model assumes negligible horizontal velocity (Siegert and Payne, 2004; Jacobel and Welch, 2005). We identified two suitable locations (Site M1 and M2 in Fig. 1) on ice divides where the model is most likely to be valid (under the assumption that the ice is at steady state in these locations). We note that other alternatives such as the Nye (1957) (supplementary material, section S7) and Parrenin and others (2006) models, or the more developed quasi-Nye model (MacGregor and others, 2015) exist. However, we do not use these to avoid the additional complexities and potential errors that would be involved. The Dansgaard-Johnsen model is used as secondary analysis and not as a method of IRH age association. We applied accumulation rates from multiple direct ice-core measurements collected at Dome A, ranging from 0.019–0.023 m water equivalent yr⁻¹, (site M2 in Fig. 1) (Minghu and others, 2011). We also applied (i) modelled estimated ages using 1-D ice compression (Dansgaard and Johnsen, 1969); (ii) joint inversion models (Wolovick and others, 2021); and (iii) isochronal IRHs calculated by Siegert (2003) and Wolovick and others (2021) of 0.016 m water equivalent yr⁻¹ and 0.014 m water equivalent yr⁻¹ respectively.

The Dansgaard and Johnsen (1969) model relies on a constant basal shear level thickness, and we applied this constant based on previous studies. Karlsson and others (2014), Ashmore and others (2020) and Bodart and others (2021) all used appropriate ranges for West Antarctica (200–1100 m) but, given our East Antarctic focus, the range we applied considered thicker estimations of the basal shear layer thickness. These estimations are based on work by Schwander and others (2001) who applied a basal shear layer thickness of 0.373H (H = ice thickness) at Dome C. We therefore applied a range of scenario estimates of basal shear layer thickness appropriate for East Antarctica, from 600–1200 m thick (Table S4, S5 and S6. and supplementary material, section S6).

2.4. Locating old ice

We used the variations in the depth of ice below the deepest IRH, H3, to identify potential regions of oldest ice across our extensive survey region. Using the fractional depth (i.e. depth of the IRH in respect to ice thickness), a spatially constant constraint was applied to the deepest IRH traced across the region. This constraint emulates an exercise undertaken by Winter and others (2019) in the Dome C region based on the understanding that the maximum age of undisturbed ice a few metres above the bed is (~800 ka). This method only recognises areas where there is a large proportion of ice deeper (and therefore older) than our deepest IRH; we therefore then excluded areas where old ice might be unexpected based on factors such as faster ice flow and high geothermal heat (Van Liefferinge and others, 2018).

3. Results

3.1. Extent and geometry of internal-reflecting horizons

H1, H2 and H3 were traced along 13 000 km of RES lines across a wide central area of East Antarctica (Fig. 2). This area includes a range of ice drainage basins that flow into Ross Ice Shelf, Filchner-Ronne Ice Shelf, the coastline of Dronning Maud Land, and Amery Ice Shelf. Across the region there is large variability in the depth of H1-H3 below the ice surface and its position as a fraction of ice thickness (Fig. 2), broadly accordant with the large spatial coverage (~582 000 km²) of the radar surveys and the large variability in ice thickness (~1300–4000 m). H1-H3 are all found deeper in the ice column near to South Pole when compared to the other parts of the analysed RES data (Fig. 2). H1, H2 and H3 are all relatively conformable to major undulations in bed topography (Fig. 3) and largely have unbroken continuous profiles.

Figure 2. IRHs traced between South Pole and Dome A displayed in terms of: a-c, depth of H1, H2 and H3 below the ice surface, plotted on bed elevation (BedMachine v2; Morlighem and others, 2020); d-f, fractional depth of H1-H3, plotted on ice thickness (BedMachine v2; Morlighem and others, 2020). Key locations are labelled on panel a, where acronyms are used for Recovery Glacier (RG); Academy Glacier (AG); Lambert Glacier (LG).

Figure 3. (a) Traced IRHs H1 (blue), H2 (red) and H3 (yellow) across intersecting radar profiles (viewed in 2D and 3D) from AGAP and PolarGap RES surveys (location map in panel b). The zoomed inset highlights challenges to IRH tracing across the Gamburtsev Subglacial Mountains. (b) Site map, showing the extent of horizon tracking in AGAP survey (yellow) and PolarGap (red) data. In the background, black lines mark the AGAP/PolarGap RES flightlines not used in this study, and the greyscale map is ice thickness from BedMachine v2 (Morlighem and others, 2020). Panels c and d zoom in on the intersection of IRHs across two radar flight lines (where the location is marked by an X arrow on panels a and b). (c) processed radargrams showing the depth of IRHs H1 (blue triangle), H2 (red circle) and H3 (yellow square). (d) traced IRHs, visible in panel c.

H1, traceable along 90% of the flightlines, demonstrates the least variability in depth, ranging from 313 to 1681 m below the ice surface. The fractional depth of H1 ranges from 0.12 to 0.64 but the mean is 0.29, highlighting that H1 is typically found in the upper half of the ice column. Lower in the ice column, H2, traced along 95% of the flightlines, ranges from depths of 645 to 2266 m and is the most extensively traceable of all three IRHs. The fractional depth of H2 is also the most variable, ranging from 0.25 to 0.94, with a mean of 0.45. H3, traced along 62% of the flightlines, reaches a depth of 2956 m below the surface. The fractional depth ranges from 0.45 to 0.95 with a mean fractional depth 0.65 (i.e. generally there is more than 35% of ice thickness below this deepest traced IRH). H3 is predominantly found in the survey grid-east of South Pole where ice is ~1600–4000 m thick (Fig. 1).

An exception to the general widespread traceability and bed-conformability of the IRHs occurs where they cross the Gamburtsev Subglacial Mountains (Fig. 3). There, exceptionally rough bed topography results in an undulating ice surface and spatially variable accumulation (Wolovick and others, 2021), causing significant dipping of IRHs, back-scatter from the mountains and, in places, loss of IRH visibility, reducing our ability to trace IRHs in this location. Elsewhere, the Recovery Subglacial Highlands also present challenges to IRH tracing. There, it was possible to trace H2 completely across the region, but surface clutter and a weaker or lost signal due to rough bed topography respectively precluded most tracing of H1 and H3 (Fig. 2). H1-H3 were also undetectable east of the onset zone of Lambert Glacier due to ice flow increasing and converging (Sanderson and others, 2023).

Previous studies have noted inconsistencies in the depth of prominent manually traced IRH when using different radar systems, largely as a result of variable central frequency and bandwidth (Ashmore and others, 2020; Bodart and others, 2021). Despite the AGAP and PolarGap data having been collected by different versions of the PASIN radar system, there is very little mismatch in IRHs at crossover locations (e.g. Figures 3c, d). An empirical error analysis of the crossovers of the traced IRHs was performed at ten intersections for AGAP-PASIN crossovers only, and for a further ten intersections between AGAP-PASIN and PolarGap-PASIN2. Root-mean-square error of the differences in H1, H2 and H3 depths generates RMS errors of 17.3 m (AGAP only) and 20.2 m (AGAP/ PolarGap). Crossover analysis for the three traced IRH is therefore low, and falls within the uncertainty range for the surveys (see section 2.2).

3.2. IRH ages

At all eight intersections between our three newly-traced IRHs and the three IRHs traced by Winter and others (2019) from Dome A through Vostok to EDC (their ‘layers H1, H5 and H8’), we found the respective three IRHs to occur at similar depths: the root mean square differences were 19 m between our H1 and their ‘H1’, 42 m between our H2 and their ‘H5’, and 51 m between our H3 and their ‘H8’ (supplementary material, section S3). Considering that the radar range resolutions of the two data sets are similar (~7 m for AGAP-South MCoRDS data (Winter and others, 2019) and ~8 m for AGAP-North PASIN data) these depth uncertainties are comparable, and we therefore conclude that, collectively, both studies have traced the same three palaeo-ice surfaces across central East Antarctica. Winter and others (2019) tied their IRHs to the Vostok and EDC ice-core chronologies, and accordingly used Bazin and other's (2013) ice-core chronology for these sites to preliminarily assign ages of 38.2 ± 0.6, 90.2 ± 1.6 and 161.1 ± 3.5 ka to the three IRHs. Through integration of the root mean square differences at the intersections with the IRH from Winter and others (2019), we therefore assign ages of 38.2 ± 2.0, 90.4 ± 3.6 and 161.9 ± 6.8 ka for H1 to H3 respectively (supplementary material, section S3 and S5).

At South Pole, the mean depths of H1 and H2 respectively are 1476.1 m and 1939.1 m, while H3 could not be traced within 3 km of the SPICEcore site because of a loss of layer visibility at depth. At the crossover with SPICEcore, we determined an age for H1 of 39.8 ka ± 0.89 ka. As a test of the IRH intersections with Winter and others (2019), this age was established independently based solely on the PolarGAP RES data and the ice core. With SPICEcore only extending to a depth of 1751 m, dated to 54 ka (Winski and others, 2019), we could not use the ice core to independently test the age of H2 other than to confirm that it is significantly older than 54 ka.

At the single intersection with a RES profile in which Beem and others (2021) traced their 37.6 ka IRH, our H1 intersected their profile 40 m higher in the ice column compared to their IRH. Because of the similarities in age assignment and apparent brightness of IRH it is probable that this layer is the same, and differences in the radar systems used have resulted in the offset of the IRHs (Beem and others, 2021).

Based on a combination of ice-core chronology age association and the intersections with previously dated IRHs, we therefore assign final ages and errors to our IRHs of 38.5 ± 2.2 ka for H1, 90.4 ± 3.57 ka for H2, and 161.9 ± 6.76 ka for H3 (further details in supplementary material, section S4 and S5). The extent and depths of each of these IRHs across our analysis, combined with Winter and other's (2019) connections to Vostok and EDC, and the IRH traced around South Pole (~H1) are shown in Figure 4.

Figure 4. Extent and depth of traced IRHs across East Antarctica including those traced for this study, by Winter and others (2019) and Beem and others (2021). Panels a–c reveal spatial changes in the depth of H1-H3 (below the ice surface), underlaid by ice thickness (BedMachine v2; Morlighem and others, 2020). Red circles denote ice-core sites. Ice-core data used to date layers traced here and from Winter and others (2019) are labelled on panel a (red star): South Pole Ice Core (SP); Vostok Ice Core (Vostok); Dome C Ice Core (EDC).

3.3. Age-depth modelling

Table 1 provides the age estimates for each IRH that were estimated from our 1-D modelling as an independent check of the ‘final’ age estimates derived above; fuller details can be found in Supplementary Table 4, 5 and 6. For H1, the modelled age range most consistent with that indicated by SPICEcore (Section 3.2) – ranging from 36.9 to 38.5 ka – was reached when assuming an average modern accumulation rate (~0.021 m w.e. yr⁻¹(Minghu and others, 2011)). Modelled average accumulation rates (Siegert, 2003; Cavitte and others, 2018; Wolovick and others, 2021) produced older ages for H1 and therefore suggest possible higher accumulation rates since ~39 ka. Wolovick and others (2021) modelled an average accumulation rate of 0.014 m w.e yr⁻¹; this value produced overestimates of the ages of all three IRHs here, at least compared to the ages derived in section 3.2. The accumulation rate most consistent with the ages we assigned to H2 and H3 (section 3.2) was 0.016 m w.e. yr⁻¹ (Siegert, 2003; Cavitte and others, 2018). Applying this accumulation rate within the model suggests that the age range for H2 is 88.3–95.9 ka, and for H3 is 163.9–219.3 ka.

Table 1. Age-depth modelling results (ka) for H1-H3 as an average between modelling sites M1 and M2

Accumulations based on modelled estimates from ^aWolovick and others (2021), ^bSiegert (2003) and ^bCavitte and others (2018) and present-day measurements (an average from three direct measurements it included here) from ^cMinghu and others (2011). For each IRH, a range of basal shear layer thickness has been modelled.

The 1-D modelling can determine the most appropriate basal shear thickness from the outputs that are the closest match to the results generated in Section 3.2. For H1, where the age determined by the SPICEcore is 38.5 ka, a basal shear layer thickness of 1000 to 1200 m is consistent with the most appropriate age range – 37.9 to 38.5. For H2 however, a basal shear layer thickness of 800 m produces the most consistent output based on the age estimation of 90.4 ka. Likewise, H3 requires a shear layer value closer to 800 m to achieve results consistent with the age estimation determined in Section 3.2 when applying an average accumulation rate of 0.016 m w.e. yr⁻¹ (Siegert, 2003; Cavitte and others, 2018). By applying a higher average accumulation rate closer to modern rates of (~0.021 m w.e. yr⁻¹(Minghu and others, 2011)), a basal shear layer thickness of 1000 to 1200 m produces a similar result to the age estimation for H3.

3.4. Old ice

The mean fractional depth of H3 (162 ka) was 0.65 across the entire traced region (Fig. 5). In parts of the study area, however, the minimum fractional depth for H3 was 0.45. This means that, in places, 55% of the ice column is below the deepest traced layer and must therefore encompass ice considerably older than 162 ka (Fig. 2f).

Figure 5. Histogram of fractional depth for H3 (162 ka) including the mean of 0.66 (where there is an average of 34.7% ice thickness below the IRH).

In Figure 6, we map the fractional depth of H3 and show the detailed spatial distribution of where the ice older than 162 ka ranges from highest to lowest along the flightlines. In this exercise, we emulate Winter and other's (2019) mapping around Dome C, where at EDC the fractional depth of H3 (162 ka ice) is 0.58, and the maximum age of undisturbed ice a few metres above the bed is ~800 ka. Applying the same 0.58 fractional-depth threshold, we identify (Fig. 6) flightline sections and regions that are most likely to be suitable for recovering old ice, before factors such as ice-flow dynamics and thermodynamics are taken into consideration. The two regions identified are the upper Byrd Glacier catchment and parts of the Gamburtsev Subglacial Mountains (Fig. 6).

Figure 6. Areas where the fraction depth for 162 ka IRH <0.58 (yellow) in relation to the fractional depth for 162 ka traced across the region (same as Fig. 2f). Most suitable areas for old-ice exploration are marked with a yellow circle. Fractional depth is underlain by ice thickness from BedMachine v2 (generated using the Quantarctica package in QGIS) for comparison (Morlighem and others, 2020). An area of high geothermal heat flux detected by Jordan and other (2018) is highlighted with a yellow dashed circle. Regional ice divides (Zwally and others, 2012) are noted in light grey.

4. Discussion

4.1. A coherent, mappable East Antarctic stratigraphy

We have mapped and dated key englacial features across a significant area of East Antarctica and, in doing so, established the first IRH links between South Pole and Dome A, and hence, by extension, Vostok and Dome C. We have traced three distinct, bright IRHs throughout the ice column, from depths of 313 m–2957 m, within the age ranges of 38.5 ka to 161.9 ka. (Fig. 4). We have shown that the uppermost layer, H1, can be dated consistently to 38.5 ± 2.2 ka in two ice-core chronologies, SPICEcore and Vostok/Dome C, that are >1000 km apart, providing extra confidence in the respective ice-core dating techniques and our treatment of IRHs as isochrones. Our study also provides further evidence (c.f., Winter and others, 2017) that it is eminently possible to combine radar data interpretations from some of the major different RES systems that have been used to survey the EAIS with little error (Figs 3c, d).

The widespread traceability of H1-H3, and the clear potential to have traced many more IRHs, is important to stress in the context that there are potentially significant mitigating factors to such an exercise across the EAIS. For example, snow ‘megadunes’ that satellite-imagery have shown to be pervasive in the EAIS interior (Fahnestock and others, 2000; Traversa and others, 2023) have been demonstrated to cause stratigraphic disruption of englacial layers (Welch and others, 2009). Although megadunes have been reported in the region east of the South Pole (Welch and others, 2009; Traversa and others, 2023), they have not prohibited our tracing of continuous IRHs through the majority of the study area. We therefore assume that megadune formation grid east of South Pole is likely to represent a relatively modern phenomenon, predating the formation of H1 (38.5 ka), and that such processes were not occurring when the layers were deposited or that the dunes did not eradicate older layers. We have noted that discontinuous IRHs are found throughout the ice column in RES data collected over the Gamburtsev Subglacial Mountains (Figs 3, 7) potentially as a result of recent or former surface erosion caused by wind scour or sublimation, perhaps linked to megadune formation and evolution (Siegert, 2003; Arcone and others, 2012; Scambos and others, 2012; Das and others, 2013; Winter and others, 2019).

Figure 7. (a) Evidence of the impact of megadunes (black box) in an example RES dataset (A11b) collected across the Gamburtsev Subglacial Mountains (where traced IRH layers are shown in blue (H1), red (H2) and yellow (H3)) (b) Little to no suggestion of megadunes despite evidence from RADARSAT of the surface features, again the traced IRHs H1 (blue), H2 (red) and H3 (yellow) from the AGAP N survey (flight line A67). (c) Site map showing the extent of horizon tracking in AGAP survey (yellow) and PolarGap (red) data and inset maps highlighting the radargrams shown in panels a and b. The background is high resolution RADARSAT-1 radar imagery of the region and we have included the multidirectional hillshade (lighting angle of 45° and azimuth of 315°) of the REMA 2 m data product (Howat and others, 2019) in the inset maps to clearly show megadunes and other distinctive features on the ice-sheet surface (Jezek, 1999).

A second potentially mitigating factor against IRH tracing that exists in our region is its complex subglacial topography; a phenomenon which in other regions of Antarctica has given rise to significant IRH discontinuities (e.g., Bingham and others, 2015). Variations in subglacial topography can influence the basal heat flux, ice-flow regime and accumulation deposition pattern, causing physical disruption to IRHs (Holschuh and others, 2017). In this survey, however, despite the complex subglacial topography, we were able to identify and map bright H1-H3 reflectors in most areas. This exercise has shown that, although they are time-intensive, manual and semi-automated layer tracing methods as applied here, can extract isochronous information across radargrams of poorer quality. They are appropriate where layers are often more difficult to trace, a scenario where fully automated methods fall short (Delf and others, 2020). However, across the Gamburtsev Subglacial Mountains, While H1- H3 are visible and traceable above the peaks in the topography and over steep mountainous slopes in the upper ice, deep reflections close to the bed are often untraceable (e.g. Figures 3, 7) This is because of the phase shift of the reflection where IRHs dip sharply from the horizontal (Holschuh and others, 2014), often hindering the tracing and dating of layers older than H3 (~162 ka) in the AGAP dataset, at least as it is currently processed.

4.2 Considerations on dating control and age-depth modelling

We have extended the EDC-Vostok-Dome A stratigraphy (Winter and others, 2019) to South Pole. While tracing over such large regions is laborious, the benefits of manually tracing englacial layers across dynamic ice and complex bed topography have outweighed any currently automated techniques (Delf and others, 2020). Uncertainties in the IRH depths arise as a result of the speed of electromagnetic wave variation through ice (Fujita and others, 1999), a firn correction requirement (Cavitte and others, 2016), and uncertainties in the range-resolution for the radar system (Cavitte and others, 2016). Combined, these factors led to a conservative uncertainty of ±17 m, ±21 m and ±27 m for H1 to H3 respectively. The larger IRH depth uncertainties for the deeper IRHs lead to greater uncertainties in age association. However, here these accounted for less than 5% of the ages between intersecting associated IRHs and we are therefore confident in the assignment of the ages between surveys. The primary uncertainties originate from the ice-core age-scales from Vostok/EDC and SPICEcore, and the connecting IRHs with other studies. Further crossover analysis with IRHs traced by Winter and others (2019) highlights a RMS error that falls within the age uncertainties of the ice cores used to date the IRHs. This demonstrates transferability across different studies that use RES data with similar vertical range resolutions.

By applying 1-D modelling as an independent method for gauging the ages of our IRHs, we have also generated information that is useful for assessing how accumulation patterns may have changed spatially and temporally over the study region, and determined an appropriate range of basal shear layer thicknesses that are suitable for future Dansgaard/Johnsen modelling of IRH ages. This latter consideration is perhaps the more important finding for future applications where IRHs may be traced across some regions of EAIS without a direct link able to be made, as here, to one or more ice cores for dating control. Age-depth modelling for IRHs H1-H3 demonstrated a higher sensitivity to the selected accumulation rate value compared to the choice of basal shear layer thickness in the model. We determined that a value between 800 and 1200 m of the basal shear layer thickness was the most appropriate for age depth modelling of the IRHs in this study. To obtain ages similar to those dated through our other methods, an accumulation rate of 0.016 m w.e. yr⁻¹ was most appropriate. Overall, the accumulation rate was the primary influence on the age-depth model (supplementary material, section S6).

4.3. Old-ice identification

The spatial extent of the deepest IRHs can inform the selection of ice-core sites for the recovery of old ice, and we have shown that in our study region candidate locations comprise the upper Byrd Glacier Catchment and the Gamburtsev Subglacial Mountains region (Fig. 6).

To retrieve old ice, ice flow throughout the history of the site should be as slow as possible, ideally stagnant, and the base of the ice sheet should not have experienced melting or refreezing since formation. Areas with the potential for basal melting over the 1.5 million year period to which present ice-core drilling initiatives aspire should be rejected (Wolff and others, 2005; Fischer and others, 2013; Van Liefferinge and others, 2018). Present-day ice flow in the upper Byrd Glacier catchment exceeds 2 ma⁻¹ (Mouginot and others, 2019), and across the catchment ubiquitously thick ice >2700 m is likely. According to the criteria of Bell (2008), this is likely to induce widespread pressure melting at the bed. Collectively, these considerations render upper Byrd Glacier Catchment likely to be an unfavourable candidate for the retrieval of very old ice – although the significant thicknesses of ice deeper than H3 demonstrate that it is still an important repository of ice significantly older than 162 ka.

Low geothermal heat flux (~55 mW m⁻²; (Martos and others, 2017)), surface ice velocities below 2 m a⁻¹ (Mouginot and others, 2019) and thick ice within valleys (>3000 m) means that it is likely that the oldest ice in the study area is within the Gamburtsev Subglacial Mountains, a finding consistent with previous research (Creyts and others, 2014; Van Liefferinge and others, 2018; Zhao and others, 2018). Despite the Gamburtsev Subglacial Mountains potentially being the most suitable region identified by this study, careful attention should be taken to consider the complex processes of ice accretion and folding due to ice flow over the rough topography (Bell and others, 2011), as well as the impact of variable surface accumulation (Fig. 7) and basal melting (Livingstone and others, 2022).

A potential motivation for drilling into deep ice in the Gamburtsev region, complementary to plans for further ice cores in the vicinity of Dome C (Chung and others, 2023), is the likely different resolutions of ice at different ages with depth between dome and dome-flank sites. At dome sites such as Dome C, the age-depth profile progresses rapidly with depth then slows near to the bed, giving greater resolution to the oldest layers. However, at dome-flank sites (as represented by the Gamburtsev region, where the overlying ice is on the flank of Dome A) age-depth profiles are typically more linear with depth, giving rise to greater resolution for intermediate-age ice (Fudge and others, 2014). We therefore posit that the Gamburtsev region may possess an important climate archive for intermediate-age ice >162 ka in East Antarctica.

4.4. Why do internal reflection horizons draw down near South Pole?

The consistent increase in IRH fractional depth observed on all flightlines within ~300 km of the South Pole (Fig. 6) manifest a drawdown of IRHs that could be due to three possible factors and associated processes. Firstly, IRH drawdown could be the result of high geothermal heat flux and a loss of ice through increased basal melt (Jordan and others, 2018). Numerous stable subglacial lakes in this region (Willis and others, 2016; Livingstone and others, 2022), combined with evidence of IRH drawdown, suggest that there is currently, or has previously been, high basal melt. Increased basal melt is often due to a combination of factors including thick ice leading to the pressure-melting point being reached at the bed (Livingstone and others, 2022) (i.e. ~2800 m deep at South Pole), elevated geothermal heat flux (Martos and others, 2017) and increased friction as a result of enhanced flow velocity (Karlsson and others, 2021). IRHs near South Pole demonstrate significant drawdown towards the centre of the ice sheet which has been attributed to a local zone of high geothermal heat flux (Fig. 6) (Jordan and others, 2018; Ashmore and others, 2020). Regional Antarctic geothermal heat flux maps suggest a peak in geothermal heat flux to the west of South Pole (Shapiro and Ritzwoller, 2004; Maule and others, 2005; Martos and others, 2017). Estimates of geothermal heat flux in Antarctica typically suggest high levels in West Antarctica linking to a major Cretaceous to Cenozoic rift system (Davey and others, 2016; Jordan and others, 2018). Our results, however, demonstrate that drawdown of IRHs originates further east than the higher background heat flux as noted by Jordan and others (2018) (Figs 2, 3a), and could potentially be triggered as a result of localised heat flux anomalies. Despite this, it is unlikely that localised anomalies would cause the drawdown to extend across a 300 km radius from the South Pole, as we see here (Fig. 2).

Secondly, IRH drawdown could have been caused by past ice dynamics. Englacial folding, representing convergent flow, can be present where there is presently slow ice flow (Bingham and others, 2015). Likewise, a drawdown of englacial stratigraphy could suggest higher latent heat caused by frictional process when sliding is enhanced, leading to increased melt at the bed (Beem and others, 2018). Although ice flow in the South Pole region is currently relatively slow (~10 ma⁻¹) (Mouginot and others, 2019), it has been hypothesised that previous enhanced flow and consequential significant frictional heat from sliding produced conditions suitable to form South Pole Lake (Beem and others, 2018). Evidence of disrupted englacial stratigraphy is commonly associated with areas of faster ice flow (Rippin and others, 2003; Siegert and others, 2003; Bingham and others, 2007; Karlsson and others, 2014), and this has been presented as supporting evidence for previous fast flow in the region (Bingham and others, 2007). The drawdown of englacial layers we evidence in this study is consistent with this hypothesis of formerly enhanced flow around South Pole and shows that the area of drawdown is much more extensive region than previously documented (i.e., extending ~300 km east of South Pole) (Fig. 2).

Thirdly, drawdown of IRHs may have resulted from increased surface accumulation near South Pole in comparison to parts of the study area further grid east. Given the low mean-annual air temperatures at South Pole (−49°C) (Lazzara and others, 2012), South Pole has a relatively high annual accumulation rate (0.08 m w.e. yr⁻¹) (Casey and others, 2014; Kahle and others, 2021). This compares to rates of ~0.025 m w.e. yr⁻¹ across the East Antarctic Plateau (Cavitte and others, 2018). Higher accumulation rates at South Pole lead to layer thickening (Holschuh and others, 2017; Born and Robinson, 2021) and therefore could lead to the drawdown of older IRHs.

Current observations and evidence do not permit us to disentangle which of these three processes is the one with the most important influence on IRH form and position in the ice column in the vicinity of South Pole. However, we note that high surface accumulation is evidenced by very robust observational evidence, while the other two influences are less well constrained or evidence for them is geographically restricted. Despite this, only a physically unrealistic accumulation rate could explain the magnitude of thickening and drawdown in this area (Leysinger Vieli and others, 2011; Beem and others, 2018). It is therefore likely that the drawdown of IRHs at South Pole results primarily from higher basal melt as a result of thick ice causing a warm ice-sheet bed, higher geothermal heat or past enhanced ice flow.

5. Conclusion

We have identified three spatially extensive IRHs and traced them through multiple RES datasets to make the first direct englacial stratigraphic connections between South Pole and Dome A, East Antarctica. Building on the work of Winter and others (2019), we have used ice-core chronologies from Dome C (that previously connected to Dome A) and South Pole to date three IRHs to 38.5 ± 2.2, 90.4 ± 3.57 and 161.9 ± 6.76 ka. We have used a 1-D ice-flow model to independently verify the IRH ages and develop our understanding of the long-term average accumulation rates and likely thickness of basal shear in this region. While the three IRHs are widely traceable across the region, complex stratigraphy in places, such as in the Gamburtsev Subglacial Mountains, mitigated some IRH tracing. However, the existence of surface megadunes did not typically impact IRH tracing at depth. Using the deepest IRH as a proxy for the age of the basal ice, we have mapped the potential for old ice and concluded that the deep subglacial valleys of the Gamburtsev Subglacial Mountains hold the most promise for its oldest ice. We have observed a drawdown of IRHs at South Pole but over a much larger spatial extent (300 km from South Pole) than identified by previous studies. We suggest that previous enhanced ice flow combined with thick ice and higher geothermal heat flux around the South Pole produce increased basal melt, leading to the drawdown of IRHs.

The dated IRHs generated here provide a valuable addition to a widespread database for dated IRHs which contributes to our knowledge of Antarctic ice-sheet architecture (e.g., https://www.scar.org/science/antarchitecture/home/). This study has notably demonstrated the existence of traceable IRHs covering a region of East Antarctica close to the Transantarctic Mountains; together with the previous work of Ashmore and others (2020) on the opposing flanks in West Antarctica, sets the stage for identifying IRHs that connect the South Pole region with West Antarctica. This would provide a valuable resource for unravelling the combined ice-sheet histories of both East and West Antarctica.

Meanwhile, such IRHs provide a growing resource for identifying candidate sites for drilling into different ages (and hence climate archives) of Antarctica's ice. Already in Antarctica's Gamburtsev Subglacial Mountains region a project under the auspices of the U.S. Center for Oldest Ice Exploration (www.coldex.org) is engaging in detailed airborne mapping using a similar radar system to that used by Beem and others (2021). Follow up ground campaigns will deploy dust logging melt probes that could further constrain local age depth profiles, testing the age structure mapped out here.

If incorporated into ice-sheet models we anticipate that our dated IRHs will provide constraints on past accumulation rates and patterns and could therefore improve our understanding of past ice-sheet evolution.