A thinner-than-present West Antarctic Ice Sheet in the southern Weddell Sea Embayment during the Holocene

Making accurate measurements and predictions of the West Antarctic Ice Sheet's (WAIS) contribution to present and future sea-level rise fundamentally depends on knowing its trajectory over the last few thousand years. We present new in situ ¹⁴C concentrations from subglacial bedrock cores collected from the southern Weddell Sea sector of the WAIS. Critically, these concentrations are above levels that can be produced under present-day ice thicknesses at the core sites. The cosmogenic nuclide inventories provide clear evidence for the ice sheet being thinner-than present at some point during the Holocene following initial thinning from its Last Glacial Maximum configuration. Forward modelling of nuclide concentrations indicates that our results are best explained by ice-surface lowering of at least 20 m. This period of thinner ice persisted for 300–3800 years and occurred after 6–4 ka. We suggest that thinning at our core sites is most likely to reflect a regional, dynamic response to grounding-line retreat rather than a localised change in ice-surface elevation. Our data are the first direct geological evidence for a thinner-than-present WAIS in the Weddell Sea sector and are consistent with Holocene retreat that culminated inboard of present-day limits. Glacio-isostatic adjustment has been inferred as a driving mechanism, causing re-grounding of floating ice and increased buttressing allowing the grounding line to stabilise and readvance. These data allow dynamic retreat-readvance behaviour of this nature to be tested in ice-sheet models, improving predictions of future sea-level rise in this critical sector of West Antarctica.

1 Introduction

Constraining the recent (i.e. Holocene) trajectory of the West Antarctic Ice Sheet (WAIS) is critical for making accurate measurements of current ice loss and for making predictions of its future evolution. In the first case, glacio-isostatic adjustment (GIA) is a key correction for current estimates of mass loss by satellite gravimetry (e.g. Caron and Ivins, 2020). Depending on the recent ice-sheet history the magnitude and/or sign of the GIA signal (i.e. uplift vs. subsidence) can be different (Bradley et al., 2015; Gomez et al., 2018; King et al., 2022). With sparse observations recent ice-loading histories are required for the GIA models that are used to correct satellite gravimetric estimates of current ice mass loss. In the second case numerical ice-sheet models used to make predictions of future sea-level change are tested (and tuned) against their ability to accurately reproduce known scenarios of past ice-sheet change (DeConto and Pollard, 2016; Albrecht et al., 2020; Pittard et al., 2022).

There is increasing evidence that parts of the WAIS retreated inboard of their pre-industrial limits before readvancing to these limits in the mid-late Holocene (for a review see Johnson et al., 2022; Jones et al., 2022). Direct geological evidence for retreat-readvance, in the form of measurable amounts of radiocarbon in subglacial sediments and inventories of cosmogenic nuclides in subglacial bedrock (Fig. 1), has been observed in both the Ross and Amundsen Sea sectors of the WAIS (Kingslake et al., 2018; Venturelli et al., 2023; Balco et al., 2023). Several numerical models also retrodict major retreat-readvance in the southern Weddell Sea Embayment (WSE) of West Antarctica (Kingslake et al., 2018; Albrecht et al., 2020 Pittard et al., 2022) where ice is only lightly grounded over the deep Robin Subglacial Basin (Ross et al., 2012). Indirect evidence, such as low present-day uplift rates (Bradley et al., 2015) and (sub-) surface morphology of ice rises (Siegert et al., 2013; Kingslake et al., 2018; Brisbourne et al., 2019), is consistent with retreat-readvance but could also be explained by major flow re-organisation without grounding line retreat (Siegert et al., 2019). Consequently, the Holocene trajectory of the WAIS in the WSE is yet to be unequivocally constrained. This is of critical importance as the Weddell Sea is a key sector for both the WAIS and EAIS, draining ∼20 % of their total grounded ice area with a SLE of >12 m (Rignot et al., 2019), and a key locale for global oceanic circulation with up to ∼50 % of Antarctic Bottom Water formation occurring there (Zhou et al., 2023).

Figure 1 Overview map of Antarctica showing locations where geological evidence of Holoence retreat inboard of the present day grounding line has been recovered. Modelled zones of retreat in the Ross and Weddell Sea embayments are taken from Kingslake et al. (2018), other models show similar patterns. Black box shows location of Fig. 2a. RSB: Robin subglacial basin. Base image is BEDMAP3 subglacial topography (Pritchard et al., 2025).

Here we present the first analyses from subglacial bedrock cores collected from the southern WSE to directly test for Holocene retreat-readvance of the WAIS in this sector. The basic premise (cf. Balco et al., 2023) of our approach is that dynamic thinning associated with grounding line retreat would expose bedrock that is presently buried by several tens of metres of glacier ice, at locations upstream from the present-day grounding line. Any Holocene thinning would result in significantly increased production rates of cosmogenic nuclides and potential exposure of bedrock surfaces to sunlight. In situ ¹⁴C (produced in rocks) has a half-life of 5730 years which means that any ¹⁴C present at the onset of the Last Glacial Maximum (LGM) would have been removed by radioactive decay during complete shielding by a thick LGM ice sheet. Thus, measurable ¹⁴C in subglacial bedrock requires production to have occurred following thinning of the LGM ice sheet (i.e. during the Holocene). Complementing this, luminescence measurements in rock surfaces are sensitive to exposure to light (Sohbati et al., 2012), thus measurements can constrain ice thinning beyond present. We present results of in situ ¹⁴C and luminescence analyses from four bedrock cores collected during the 2022–2023 Austral summer in the Enterprise Hills, a constituent massif of the Heritage Range within the Ellsworth Mountains, southern Weddell Sea Embayment.

2 Site description and regional ice-sheet history

The Enterprise Hills are located within the Heritage Range, ∼30 km from the present-day grounding line at Hercules Inlet (Fig. 2a and b). While further rock outcrops closer to the grounding line, available satellite imagery suggests extensive snow/firn at these sites which prevents deployment of the Winkie Drill as modified at Durham University (cf. Boeckmann et al., 2021; Braddock et al., 2025). The Enterprise Hills delimit the northern flank of Horseshoe Valley within which ice flows in a general north-westerly direction towards the grounding line. Several glaciers breach the northern bounding escarpment of Horseshoe Valley and merge with extensive areas of blue ice found on the leeward side of the Enterprise Hills (Fig. 2b). At the mouth of one such glacier (Plummer Glacier) is a small nunatak, informally named Plummer Nunatak, that reaches ∼60 m above present-day ice. Plummer Nunatak is flanked by a small outlier (∼10 m high) on its northern side. A further ∼250 m from this outlier is a small patch of exposed bedrock (Fig. 2c). This outcrop (informally referred to as Patch Rock) is surrounded by a small snow field that, where probed was always <2 m thick and directly/unconformably overlies blue ice (i.e. there is no snow-firn-ice transition). The Enterprise Hills are composed of quartzose sandstone (sometimes referred to as quartzites) of the Crashsite Formation (Spörli, 1992) and thus contain quartz suitable for cosmogenic nuclide analyses.

Figure 2 (a) Map of the southern WSE showing major ice streams and geographical features mentioned in the text. The dashed line is zone of modelled retreat (Kingslake et al., 2018). Red star is location of our drill site. (b) Location map of Horseshoe Valley and Enterprise Hills. Location of Plummer glacier is shown. Red star is location of our drill site. $\text{MH}/\text{PH}/\text{IH}=\text{Marble}/\text{Patriot}/\text{Independence}$ Hills where cosmogenic nuclide data evidences post LGM thinning. (c) Google Earth image (Maps data © Google, © 2025 Maxar Technologies) of Plummer Nunatak. Patch Rock is visible and locations of boreholes shown. Figure generated in Quantarctica version 3 (Matsuoka et al., 2021).

Following the LGM the WAIS in the region of the southern Heritage Range thinned rapidly by ca. 400 m during the mid-Holocene with ice reaching present day elevations by ca. 6–4 ka as evidenced by exposure age data from the Marble, Patriot, and Independence Hills (Bentley et al., 2010; Hein et al., 2016). Additional evidence from an ice core at Skytrain Ice Rise also suggests rapid surface lowering of ca. 450 m at 8 ka (Griemann et al., 2024). Given the consistency between these records, it seems reasonable to infer a similar and consistent ice-sheet history across the Heritage Range, including at our site in the Enterprise Hills. Specifically, this history is inferred to be a significantly (ca. 400 m) thicker ice sheet at the LGM which subsequently thinned rapidly during the Holocene to reach present day ice levels by ca. 6–4 ka. The two alternate scenarios for the subsequent ice history are (1) ice remained at that level for the remainder of the Holocene, or (2) ice thinned below present post 6–4 ka (cf. Johnson et al., 2022; Jones et al., 2022; Balco et al., 2023).

3 Methods

3.1 Subglacial bedrock sampling

To aid site selection we used a PulseEkko Pro Ground Penetrating Radar (GPR) system equipped with 100 MHz antennas, arranged in a perpendicular broadside configuration. The antennas, transmitter, and receiver were attached to a plastic sledge and hand-towed across the survey area at $\sim \mathrm{1}\mathrm{m}{\mathrm{s}}^{-\mathrm{1}}$. Due to equipment failure, we could not co-record locational information, so survey lines were marked out by flags at 25 m intervals to allow the data to be rubber-banded in post-processing (see Supplement). Following initial survey immediately prior to commencing drilling activities, we identified suitably shallow bedrock extending north from a small exposure of bedrock close to the main outcrop of Plummer Nunatak (Fig. 3). The ridge crest maintained a consistent depth beneath the ice surface (ca. 7–10 m) before it passed below a steep ramp of blue ice. As a result, and to obtain cores from a range of depths, we decided to drill at two locations offset from the main ridge crest (BH02 and BH03). This introduced additional uncertainty in our estimates of ice depths due to off-nadir reflectors and restricted time for post-processing. That said, as the depths involved were comfortably within range of the Winkie Drill this was not considered to be a critical issue.

Figure 3 (a) Patch Rock viewed from the top of Plummer Nunatak. The locations of the bedrock cores are shown. The black dashed line denotes the approximate trace of the subglacial ridge as surveyed before it passes underneath the blue ice ramp. Vertical relief of the ramp is ca. 10 m. GPR Line 5 is shown as a red dashed line in (a). (b) Interpreted radargram of GPR Line 5 showing locations of cores BH01 and BH03 along with measured depths to bedrock (n.b. for BH03 this is greater than depth extrapolated from two-way-travel time). The numbers 5.2 etc. denote locations of flags used as fiducial markers. GPR Lines 1–4 were located at 10 m intervals in increasing proximity to Patch Rock (see Supplement).

We used a modified Winkie Drill to collect subglacial bedrock samples from depths up to 29.5 m below present-day ice levels. Our drill set-up is similar to the US Ice Drilling Program Winkie Drill deployed to the Ohio Range in 2016/2017 (Boeckmann et al., 2021; Braddock et al., 2025) being powered by a 2-stroke internal combustion engine rather than an electric powerhead. As all of our boreholes were drilled directly into blue ice we used modified Kovacs ice augers attached to the Winkie Drill transmission to make our access holes. The augers performed well at depths up to 29.5 m, clearing>95 % of chippings from the boreholes. Unlike some previous drilling campaigns we did not encounter any sediment-laden ice at the base of our boreholes (cf. Braddock et al., 2025). We attribute this to some combination of persistent cold-based conditions in the blue ice area (i.e., ice thickness during the LGM was not sufficient to produce erosive warm-based ice) and the nature of the quartzite bedrock which is highly cemented and thus unlikely to produce fine grained weathering products. Once bedrock was reached, we replaced the augers with an AW34 drill rod string and cored bedrock using an IAWS core barrel using D60 drill fluid in a forward fluid circulation configuration. As all our boreholes were created directly into blue ice there was no requirement to case the hole. We would also note that working in a blue ice area almost entirely negates hazards from unidentified crevasses.

In a total of 11 d of drilling activity we drilled six access holes to depths of 7.0–29.5 m below present-day ice (Table 1). We successfully recovered bedrock cores from four of these boreholes, BH01, BH02, BH03, and BH06 (Table 1 and Fig. S2 in the Supplement). All cores were composed of medium-fine grained quartz arenites of the Crashsite formation. In BH04 we were unable to establish fluid circulation, presumably due to fluid escaping in through a fracture in the bedrock. At BH05 we hit bedrock at an unexpectedly shallow depth (7 m) and decided against sampling. At BH01 we recovered 10 cm of bedrock in our first core run but during the second core run the core barrel became irretrievably stuck. This was due to inefficient clearance of chippings following premature cessation of fluid circulation (pack off).

Table 1 Borehole locations, elevations, and core recovery. n/a –not applicable

3.2 Luminescence measurements of BH01

BH01 was recovered in the dark to facilitate analysis using rock surface luminescence (Sohbati et al., 2011). The sample was prepared at the University of Liverpool luminescence laboratory under subdued-lighting conditions to prevent sunlight contamination of the luminescence signal. Due to the irregular core-top morphology of BH01, only a single luminescence-sample core could be drilled in the luminescence laboratory (i.e. a sub-sample of the larger core collected by the Winkie Drill). This subsample was sliced to a thickness of ∼0.7 mm and slices mounted in stainless steel cups for luminescence measurements. In hand specimen our samples appeared to be dominated by quartz thus luminescence measurements followed standard protocols employed at the University of Liverpool luminescence laboratory utilising the blue signals of quartz detected in UV wavelengths (see Supplement). Luminescence depth profiles were determined by measuring the natural signal (L_n) normalised using the signal measured in response to a test-dose (T_n), termed the $L_{\mathrm{n}}/T_{\mathrm{n}}$ signal. Dose-recovery experiments were consistent with a suitable analysis protocol (see extended discussion in the Supplement).

3.3 In situ cosmogenic ¹⁴C analyses

Of the remaining cores, BH02 and BH03 were sliced at 5 cm intervals while BH06 was processed as a single sample. After crushing and milling at the University of Durham quartz was purified at the Scottish Universities Environmental Research Centre (SUERC) Cosmogenic Nuclide Facility. Given the assumed high quartz content and low feldspar content of the samples froth flotation was not performed. Samples were subject to a total of 3–5 successive 48 h leaches in 1 % $\mathrm{HF}/{\mathrm{HNO}}_{\mathrm{3}}$. Sample purity was assessed using ICP-OES with Al contents<200 ppm considered sufficiently pure. A total of ten samples were sent to the ANSTO-University of Wollongong (ANSTO-UoW) in situ ¹⁴C extraction laboratory (Fülöp et al., 2019) for carbon extraction. Prior to extraction samples underwent an additional pre-cleaning using a protocol modified from Fülöp et al. (2019). The ¹⁴C extraction scheme at ANSTO exploits the high temperature phase transformation of quartz to cristobalite to quantitatively extract the carbon as CO₂. As part of this process samples are pre-heated to 600 °C for 1 h to remove meteoric ¹⁴C. Impurities of chlorite and/or glauconite produced reactions in two of our samples (BH02 0–5 cm and BH02 10–15 cm) presumably due to presence of Cl and/or K. As a result two sample tubes cracked during in vacuo cleaning. Following recovery of these samples from the coil-heaters and subsequent carbon extraction distinct discolouration of the residual material suggested contamination by the coil-heater itself and the initial samples were disregarded (Fig. S10). For one sample (BH02 0–5 cm; OZ12602) additional material was recovered from the fine (75–125 µm) fraction which was subsequently extracted. All other extractions were carried out on the 125–500 µm fraction. We thus present results from a total of nine samples (Table 2). Following extraction CO₂ gas is converted to graphite targets and analysed using ANSTO's ANTARES 10MV tandem accelerator (Fink et al., 2004; Smith et al., 2010). Samples are normalised to small OX-II standards, corrected for a global graphitisation blank by 1/m analysis and data subsequently reduced according to Hippe and Lifton (2014).

Table 2 Results of in situ ¹⁴C analyses of samples from bedrock cores. Total inventories are not blank corrected; concentrations include global blank correction of 15 500±9500 atoms.

* 75–125 µm size fraction, all other samples utilised 125–500 µm.

3.4 Forward modelling of in situ ¹⁴C concentrations

To explore scenarios of ice thickness change that are consistent with the measured concentrations we use a simple forward model for cosmogenic nuclide production under 10 000 random three-stage ice histories (Fig. 4). Stage 1 has a random start time ($T_{\mathrm{1}}=\mathrm{10}-\mathrm{6}\mathrm{ka}$) and represents the early Holocene when the regional ice-sheet history mandates that our site was covered by an ice that was thicker-than-present. During Stage 1 ice thickness (H_S1) is constrained to be no thinner than 50 m to allow for some production by deeply penetrating muons to assess sensitivity of the model to initial ice thickness. Stage 2 (T₂−T₃) is a period of thinner-than-present ice that begins at a random time after T₁. In all boreholes we observed a clean ice-rock transition at the bed. Balco et al. (2023) used the presence of a dirty ice layer at the base of their boreholes to infer that ice had not thinned beyond their core top elevations during the Holocene. The absence of such a layer at our sites would, by the same logic, imply that ice could have thinned below our core top elevations but does not mandate that it did so. Consequently, during Stage 2 ice thickness (H_S2) is allowed to randomly vary between 0 m (i.e. total exposure of our cores) and its present-day thickness with 10 % of model iterations set to zero thickness. This allows our model to explore what duration of complete exposure is consistent with our measured concentration and what thickness of residual ice cover is possible if the assumption of complete exposure is incorrect. Finally, Stage 3 (T₃−T₄) represents the current period of burial by ice with its present-day thickness (H_S3). The start-time of Stage 3 is also random but with a minimum value of 60 a as available TMA aerial imagery shows ice at our site was in its present-day configuration in 1962 (Fig. S3).

Figure 4 Schematic of the forward model we use to predict downcore nuclide concentrations under varying ice thickness histories. T₁ is when the model begins, T₂ is the onset of hypothesised Holocene low stand, T₃ is end of low stand, and T₄ is the present day. H_S1, H_S2, and H_S3 refer to the thickness of the ice column during Stages 1, 2, and 3 (S1–3) respectively.

For each stage (S_x) we calculate total mass depths in the core ($z_{\text{sum}}=z_{\left(i\right)}+z_{\left(c\right)}$) where z_(i) is the mass depth (g cm⁻²) at the bottom of an ice column H (cm) thick with a density (ρ_ice) of 0.917 g cm⁻³:

$$\begin{array}{}\text{(1)}& z\left(i\right)=H\left(S_x\right)\cdot {\mathit{\rho }}_{\text{ice}}\end{array}$$

the mass depths in the core (z_(c)) are constant for each stage of the model and defined by:

$$\begin{array}{}\text{(2)}& z\left(c_x\right)=D\left(x\right)\cdot {\mathit{\rho }}_{\text{rock}}\end{array}$$

where D(x) is the depth in the core (cm) and ρ_rock is an assumed rock density (2.68 g cm⁻³). The total mass depths (z_sum) are used alongside the surface elevation of the ice sheet during the corresponding stage to calculate downcore nuclide production rates during each stage of our model (P_S1, P_S2, P_S3 in $\mathrm{atoms}{\mathrm{g}}^{-\mathrm{1}}{\mathrm{yr}}^{-\mathrm{1}}$). Production rates are calculated as the sum of spallogenic and muonogenic production as implemented by Balco et al. (2023) using the scaling method of Stone (2000).

The resulting nuclide concentration for each stage of our model (N_Sx) are calculated at each depth in our core as:

$$\begin{array}{}\text{(3)}& \begin{array}{rl}{\displaystyle }{N}_{\mathrm{S}x}=& {\displaystyle }{N}_{\mathrm{S}x-\mathrm{1}}\exp (-\mathit{\lambda }\cdot T_{\mathrm{S}x})\\ {\displaystyle }& {\displaystyle }+\left({\displaystyle \frac{P_{\mathrm{S}x}}{\mathit{\lambda }}}\cdot (\mathrm{1}-\exp (-\mathit{\lambda }\cdot T_{\mathrm{S}x}\left)\right)\right)\end{array}\end{array}$$

Where N_Sx−1 is the nuclide concentration at the end of the preceding stage (assumed to be zero at the start of the model), λ is the decay constant of ¹⁴C ($\mathrm{1.2097}\times {\mathrm{10}}^{-\mathrm{4}}{\mathrm{yr}}^{-\mathrm{1}}$), T_Sx is the duration of the stage (S_x) in years (e.g., $T_{\mathrm{S}\mathrm{1}}=T_{\mathrm{1}}-T_{\mathrm{2}}$), and P_Sx is the depth specific production rate for the relevant stage of the model ($\mathrm{atoms}{\mathrm{g}}^{-\mathrm{1}}{\mathrm{yr}}^{-\mathrm{1}}$).

In all boreholes we encountered ice that was completely frozen to its bed. Considering low mean annual temperatures at our site ($<-\mathrm{20}\mathrm{°}\mathrm{C}$; cf. Wang et al., 2023), basal melting under present-day ice thickness is glaciologically implausible making the occurrence of significant erosion during the most recent period of (thin) ice cover (Stage 3) unlikely. Similarly, while sub-aerial erosion does occur in Antarctica quantified rates are generally very low ($<\mathrm{1}\mathrm{mm}{\mathrm{ka}}^{-\mathrm{1}}$), including in quartzites (Marrero et al., 2018). Combined with observations of striations on polished bedrock surfaces above ice we consider it likely that sub-aerial erosion during any period of Holocene exposure (Stage 2) was negligible. Finally, any erosion during thick ice cover (Stage 1) is not relevant to our model as the short half-life of ¹⁴C (5730 years) requires any pre-LGM inventory to have decayed regardless of whether erosion occurred. Consequently we assume zero erosion in all stages of our model.

4 Results

4.1 Luminescence depth profile (BH01)

The results of the luminescence measurements (Fig. 5a) are inconsistent with the regional ice-sheet history. Specifically, we calculate an apparent burial age of ca. 21.4±1.6 ka (Fig. 5b) for slices 2 and 3. Such an age would imply exposure to light during the LGM when the ice-sheet surface across the WSE was several hundred metres thicker than present and the grounding line was >300 km outboard of its present-day limits (Bentley et al., 2010; Hillenbrand et al., 2014; Hein et al., 2016; Nicholls et al., 2019; Johnson et al., 2019). Although the sample passed dose-recovery experiments, we conclude that the luminescence signals of these rock slices were not dominated by quartz and thus were not suitable for measuring burial ages using the applied protocol. Previous studies have reported dim luminescence signals measured for quartz from rock slices (e.g. Sohbati et al., 2011) similar to the signals reported here for sample BH01 (Fig. S8); this is likely due to the lack of sensitisation of the signal prior to burial within a rock. The SEM-EDS measurements of sample BH01 (Fig. S6) show that K is present in the rock slices which is potentially related to the presence of K-feldspar. The presence of this mineral would be consistent with our measurements of high fading rates ($\mathrm{9.2}\pm \mathrm{2.3}\mathit{\%}/\text{per decade}$). It has previously been shown that feldspars emit detectable signals when stimulated with blue wavelengths and detected in UV (e.g. Thomsen et al., 2008). In such a scenario, the dim quartz signals give greater opportunity for minerals with less suitable luminescence properties (e.g. feldspar when stimulated and detected in blue and UV wavelengths, respectively) to dominate the signal used for dating.

Figure 5 Luminescence depth plots for sample BH01 showing D_e values (a) and apparent burial ages (b) measured using the OSL signal of quartz. Note that only Discs 2 and 3 (depths 1.5 and 2.5 mm, respectively) were used for age calculation.

Our measured signal it is not a quartz signal or a K-feldspar signal stimulated by infra-red wavelengths typically used for dating. The unsuitable luminescence properties could thus manifest in both age underestimation (because of the high fading rate), but also in age overestimation as the signal may not deplete rapidly in response to sunlight exposure (i.e., for these minerals, at this wavelength) While the rock appeared to be dominated by quartz in hand specimen it is notably fine grained and other K bearing minerals are clearly present. In future, the mineralogy should be more confidently quantified prior to choice of measurement protocol given the difficulty in collecting such samples. Ultimately, and unfortunately, the luminescence results from BH01 are not interpretable.

4.2 In situ ¹⁴C concentrations in bedrock cores

The nine samples analysed for in situ ¹⁴C yielded total inventories of $\mathrm{4.5}\times {\mathrm{10}}^{\mathrm{4}}-\mathrm{9.1}\times {\mathrm{10}}^{\mathrm{4}}\mathrm{atoms}$ ¹⁴C (Table 2). The reported total inventories of ¹⁴C are $\sim \mathrm{3}-\mathrm{6}$ times higher than the long-term ANSTO-UoW global blank value of $\mathrm{1.55}\pm \mathrm{0.94}\times {\mathrm{10}}^{\mathrm{4}}\mathrm{atoms}$ ¹⁴C which represents the average and standard deviation of all process blanks measured at ANSTO over the last 6 years. We processed and measured 15 individual process blanks alongside the samples presented here. The values from these 15 blanks are indistinguishable from the global blank value and show no time-dependence (Figs. 6 and S11). None of the samples overlap with any individual blank or with the long-term blank value at 68 % confidence. Only one sample overlaps with the global blank at 95 % confidence. As argued by Balco et al. (2022), a key point is that if upper core samples contain inventories of ¹⁴C that are above background then, given production systematics, ¹⁴C must be present in down core samples even if those samples have inventories that are close to blank levels. The uppermost samples from BH02 (5–10 cm), BH03 (Core top), and the single sample from BH06 are 2.9, 4.3, and 3.4 times above blank respectively, and do not overlap with the blank at 95 % confidence. We thus conclude that in situ ¹⁴C is present in our cores. The generally observed downcore decrease in concentrations of ¹⁴C (atoms g⁻¹) in both BH02 and BH03 is consistent with production systematics of cosmogenic nuclides and further support this conclusion.

Figure 6 Probability distribution functions showing total measured inventories (total number of atoms) of ¹⁴C in bedrock samples (red) alongside the 15 individual process blanks and the six-year global blank value of 15 500±9500 atoms (blue). Upper core samples are highlighted.

5 Discussion

5.1 Was ice thinner-than-present during the Holocene?

In situ ¹⁴C can be produced at depths of tens of metres by deeply penetrating muons (cf. Lupker et al., 2015). Thus, it is possible that at least some of our measured concentrations could be produced without ice thinning beyond present if the duration of any quasi-constant ice cover was sufficiently long. To test this scenario we modelled expected downcore nuclide concentrations under constant (present-day) ice thickness for periods of 2, 6 and 8 ka. This represents a conservative time range as defined by existing constraints on ice-sheet history in the region (Bentley et al., 2010; Hein et al., 2016; Griemann et al., 2024). Notably, the predicted concentrations are significantly lower than the measured concentrations of the upper samples in both BH02 and BH03 (Fig. 7). In addition, the “straight” shape of the modelled profiles reflects the dominance of muonogenic production under $\sim \mathrm{25}-\mathrm{29}\mathrm{m}$ of ice (i.e. the present thickness at the sites of BH02/BH03). In contrast, the data from BH02 and, in particular, BH03 define a production profile with a spallogenic contribution. These two observations support an inference that the nuclide concentrations in our bedrock cores result from production during a period when ice at our core sites was thinner-than-present allowing for increased production to occur.

Figure 7 Downcore measured nuclide concentrations for and predicted depth profiles under invariant ice cover for durations of 2 ka, 4 and 8 ka for samples from BH02 (a), BH03 (b) and the single sample from BH06 (b). Solid box and outline denote 68 % and 95 % uncertainties on measured nuclide concentrations respectively.

Unlike the samples from BH02 and BH03, the ¹⁴C concentration measured in the single sample from BH06 is consistent with, but does not mandate, production during prolonged (8–6 ka) burial by ice of present-day thickness which at this site is 9.5 m. However, given their respective locations, if ice thinning exposed the BH02/03 core sites then the BH06 core site must also have been ice free. Similarly, while the OSL results presented here cannot be interpreted with respect to past exposure of the BH01 core surface to sunlight, as with BH06 which is at a comparable depth, if ice thinning occurred at the sites of BH02 and BH03 then ice thinning must also have occurred at the site of BH01. We cannot exclude the possibility that BH06 was shielded by sediment and/or snow during any Holocene low stand, However, without observations such a scenario remains speculative. Importantly, we emphasise that the nuclide concentrations in both BH02 and BH03 show an unambiguous signal of increased production above that possible due to muonogenic production under present-day ice thicknesses thus requiring ice to be thinner at our core sites at some time during the Holocene.

5.2 Ice surface histories inferred from forward model of nuclide concentrations

Our simple forward model of nuclide accumulation allows us to predict downcore nuclide concentrations in our bedrock cores under various synthetic ice histories (Fig. 8). Regardless of the relative contributions of different production pathways in situ ¹⁴C production rates decrease with depth (cf. Gosse and Philips, 2001). Consequently, as all depths in any individual core have inherently experienced the same exposure-burial history the actual (i.e. not measured) concentrations must decrease with depth. To do an assessment of model fit we calculate a reduced chi-squared misfit between the concentrations predicted by any individual synthetic ice-sheet history and the measured concentrations in our cores. Histories that can't be rejected at 98 % confidence (i.e., with a reduced χ²<3.279 for 3 degrees of freedom) were retained for further analysis and discussion. This results in a total of 1493 and 244 fitting iterations for BH02 and BH03 respectively. We acknowledge that the 98 % threshold is a relatively conservative criterion. However, use of a narrower threshold (e.g., 95 %) reduces the number of fitting iterations for BH03 to effectively zero (n=2) primarily because, as pointed out by Balco (2025), the three uppermost measurements in BH03 have concentrations that randomly decrease with depth more rapidly than a surface spallation profile. Given that an overall pattern of decreasing concentrations with depth is observed (e.g., the upper two and lowermost samples are physically consistent with production systematics) we suggest that BH03 contains more useful information on the range of ice histories that are broadly consistent with its profile than we can meaningfully discuss from only 2 fitting iterations at the 95 % threshold.

Figure 8 (a) Predicted nuclide concentrations plotted against depth in core BH02. Measured concentrations shown as in Fig. 7. Concentrations predicted by the best (10th percentile) of fitting iterations shown in green, all fitting iterations shown in dark grey, random selection (n=1000) of remaining iterations shown in light grey. (b) Three stage ice histories used to produce modelled concentrations, colours correspond to the same sub-sets of iterations as panel (a). The histories, which all have different starting times in the model, have been aligned to the onset of the low stand for visualisation. (c, d) are the same figures but for BH03.

While some accepted scenarios for both cores are consistent with a Holocene low stand during which our drill sites were completely ice free (Fig. 9) our forward model, which is under constrained, also allows for scenarios with thin residual ice cover, but relatively longer low stand durations. Put simply, if production rates in the cores are reduced due to thin residual ice cover, then a longer time is required to accumulate the concentrations subsequently measured. For BH02 fitting iterations of our forward model predict low stand durations of 300–9200 years corresponding to residual ice thicknesses of 0–15 m. For BH03 the forward model predicts a narrower range of ice histories with low stand durations of 550–8250 years corresponding to residual ice thicknesses of 0–7 m. For BH02 the model appears to favour scenarios with intermediate residual ice thicknesses ($\sim \mathrm{2}-\mathrm{8}\mathrm{m}$) whereas for BH03 scenarios with complete exposure of the core top produce lower misfit scores (larger symbols in Fig. 9). This difference is potentially due to randomness within the measured concentrations from the limited number of samples available thus we prefer to present the results of the forward model conservatively in terms of the overall range of low stand durations and ice cover thickness represented by the fitting iterations. Considering the full population of accepted scenarios the range of possible low stand durations is $\sim \mathrm{300}->\mathrm{9200}\text{years}$ with the upper 10th percentile of this range requiring residual ice thicknesses of $\sim \mathrm{5}-\mathrm{15}\mathrm{m}$ (Fig. 9). This overall range of possible solutions can be narrowed further by considering the total time available for thinning and re-thickening as constrained by the most proximal cosmogenic nuclide data from the Heritage Range. This data suggests ice thinned to near its present-day thickness by 4 ka (cf. Bentley et al., 2010; Hein et a., 2016) which would limit the range of accepted solutions to encompass low stand durations of 300–3800 years and maximum residual ice thicknesses of ∼7 m (Fig. 10).

Figure 9 (a, d) Low stand duration plotted against thickness of residual ice cover for fitting model results from BH02 (a) and BH03 (b). The 10th percentile of best fitting scenarios has symbol sizes scaled in order of goodness-of-fit showing that scenarios of thin ice better fit the measured concentrations in BH03. (b, e) and (c, f) are the respective histograms for the two variables, again showing that shorter, “none-or-thin” ice histories are preferred.

Figure 10 XYZ plot of durations of low stand (Stage 2) and period of present ice thickness (Stage 3) against residual ice thickness during Stage 2 for cores BH02 (a) and BH03 (b). Nearest exposure age data imply ice surface reached modern elevations at ca. 6–4 ka thus total duration of Stage 2 and Stage 3 would be shorter than this. Such constraints narrow the range of possible low stand durations and residual ice thicknesses as shown for periods of 4000 years (blue markers) and 6000 years (red markers). Remaining accepted scenarios shown with white markers.

5.3 Wider implications for Holocene ice-sheet history in the Weddell Sea Embayment

In summary, the data requires thinner-than-present ice and is most consistent with ice thinning at our core sites of at least ca. 20 m defined by the total time available for thinning and re-thickening. Our data does not exclude the possibility that our core sites were completely exposed so this value is considered a minimum constraint for mid-late Holocene (i.e., post 4 ka) ice-sheet thinning at our sites. Given that drill site choice is dominated by logistical parameters there will always be a trade-off between where we would ideally go to answer a scientific question and where drilling is likely to be successful. A question that must be addressed therefore is whether ice-sheet thinning at any given location (as represented by data such as ours) represents a broader lowering of the regional ice surface (i.e. due to dynamic thinning associated with grounding-line retreat) or a (very) localised change due to site-specific factors such as wind scoop development and/or locally enhanced ablation in the BIA due to stronger katabatic winds. Ultimately resolving this question, for example by using nested higher resolution grids around potential drill sites in ice sheet models to quantitatively link ice surface lowering to grounding line position is a challenge that remains to be addressed. However, there are multiple lines of evidence to support a scenario of thinning at our site being representative of a regional response.

Firstly, wind scoops are depressions eroded by vortexes formed when prevailing winds are obstructed (e.g., by protruding nunataks). Large examples can be eroded into glacier ice and can have floors several tens of metres below the adjacent ice surface (cf. Ackert et al., 2007). Wind scoops often form on the upwind side of obstructions, but their overall distribution reflects complex interactions of varying wind directions and local (meso-scale) topography. Presently, there is very limited wind scoop development around Plummer Nunatak, and none around the small outlier or Patch Rock. Given the lack of vertical relief around Patch Rock (i.e., <2 m above ice) it is difficult to envisage development (and subsequent disappearance) of a wind scoop of sufficient size to expose our sample sites. Additionally, BH02 and BH03 are located on opposite sides of the subglacial ridge thus even if a wind scoop did develop it is unlikely that one would form on both sides of such a low relief outcrop.

Secondly, our core sites are located in an extensive Type 1 (cf. Bintanja, 1999) blue ice area (BIA) on the leeward side of the Enterprise Hills (Fig. 2b). Type 1 BIAs form where katabatic winds, enhanced by topography, focus sublimation on the lee side of nunataks (Bintanja, 1999). The first order control on the surface elevation of BIAs is regional thinning/thickening of the ice sheet over glacial timescales (cf. Woodward et al., 2022). Once a close approximation of the current ice-sheet configuration was reached (i.e. ca. 6–4 ka) the predominant mechanism for localised lowering of the BIA surface elevation (in the absence of dynamic thinning) would be a situation where katabatic winds are enhanced, leading to the rate of surface ice loss exceeding the rate of compensating ice flow. Katabatic wind strength is the product of continuous cooling over a sloping ice-sheet surface (Van den Broeke and Bintanja, 1995; Bintanja et al., 2014) thus wind strength increases with surface gradient. Minor elevation changes (<100 m) at the WAIS divide during the Holocene (cf. Koutnik et al., 2016) would result in negligible changes in surface gradient over the >500 km between the WAIS divide and the Enterprise Hills. Additionally, modelled future changes in katabatic wind strength due to anthropogenic forcing are small (Bintanja et al., 2014). This implies that, in the absence of major ice-sheet elevation changes, variations in Holocene wind strength, when atmospheric forcing was of a similar (or smaller) magnitude, would also be minor.

Thus, we suggest that the most parsimonious explanation for a thinning, and subsequent re-thickening, of the ice surface at our sites is a dynamic change associated with fluctuations of the grounding line. If this is the case, then the surface elevation change observed at our sites is indicative of broader retreat-readvance of the WAIS in the southern WSE. The timescales involved in thinning and re-thickening of the ice surface, as constrained by the cosmogenic nuclide data from BH03, are consistent with observations of late Holocene (<4 ka) flow reorganisation over Bungenstock and Korff Ice Rises (Siegert et al., 2013; Brisbourne et al., 2019). Similarly, the modelled duration of the low stand (500–3800 years) agrees with predictions from GIA modelling constrained by GPS observations in the southern WSE (Bradley et al., 2015). Notably, Bradley et al. (2015) observe that their GIA model could potentially produce a better fit to GPS observations if there was “an increase in the spatial extent of the retreat-readvance”. The ice-loading histories used as inputs to their GIA model focus retreat over the RSB, ca. 200 km south of our site. Ice-sheet model simulations of past, and potential future retreat in the southern WSE also focus grounding line retreat in the RSB but generally not in the near vicinity of our site (Kingslake et al., 2018; Albrecht et al., 2020; Hill et al., 2021; Pittard et al., 2022). Our data therefore implies that Holocene retreat-readvance of the WAIS in the southern WSE may have been more extensive than previously hypothesised. Constraining the overall extent of this change is a future research priority (Ross et al., 2025) and data from sites near the southernmost portion of the WSE, at the margin of the Foundation Ice Stream, are forthcoming (Small et al., 2024).

5.4 Comparison to evidence for retreat-readvance from other sectors of West Antarctica

While being mindful that different methods constrain different glaciological phenomena (i.e., grounding line retreat vs. ice-surface elevation change) it is worth comparing our data with evidence for Holocene retreat-readvance in other sectors of the WAIS. In the Ross Sea, radiocarbon dating of sub-glacial organic carbon recovered at and inboard of the Siple Coast grounding line indicates ca. 250 km of grounding line retreat in the mid- Holocene ca. 6.3–7.2 cal kyr BP, before readvance to the present-day position at an unconstrained time (Kingslake et al., 2018; Venturelli et al., 2020, 2023). In the Amundsen Sea cosmogenic nuclide data from exposed nunataks, subglacial bedrock, and a controlled moraine (cf. Evans, 2009) constrains a minimum low stand duration of 3000 years that occurred sometime after ca. 7.1 ka with readvance occurring at ca. 1.4 ka (Adams et al., 2022; Balco et al., 2023; Nichols et al., 2024).

Our new data from subglacial bedrock cores in the WSE is consistent with these previous studies although we cannot constrain the timing of ice-thinning directly. Indirectly, the nearest available above-ice cosmogenic nuclide data suggests that the ice surface initially reached it present elevation at ca. 6.0–4.0 ka (cf. Bentley et al., 2010; Hein et al., 2016). If this was the case then excess thinning in the WSE must post-date this time, broadly in line with the timings in both the Ross Sea and Amundsen Sea sectors. The most likely duration of the low stand in the WSE (300–3800 years) is, towards its lower limit, somewhat shorter than in the Amundsen Sea (>3000 years) although we note that our sites are located at higher elevations ($\sim \mathrm{500}\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$ vs. $\sim \mathrm{80}\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$ (meter above sea level)) and further from the present-day grounding line (∼30 km vs. ∼2 km) than the core sites of Balco et al. (2023). Thus, the difference in modelled low stand duration may reflect a real variation or simply a delay between grounding-line retreat-readvance and corresponding ice-surface elevation changes at our upstream core sites. Our data, therefore, does not preclude the possibility that there was near synchronous Holocene retreat-readvance of the WAIS in three of its main drainage catchments. Similarly, our data (and the other data constraining WAIS retreat-readvance) are consistent with recent global mean sea level (GMSL) reconstructions based on GIA modelling that suggest excess ice loss from Antarctica contributed to a higher than present GMSL since 7 ka (Creel et al., 2024).

6 Conclusions

The measured concentrations of in situ ¹⁴C in two of our subglacial bedrock cores, combined with the observed downcore decreases in these concentrations, provides diagnostic evidence for a thinner-than-present ice sheet at our core sites during the Holocene. Our concentrations and their depth profiles cannot be explained by an ice surface that was maintained at its present elevation following the cessation of post-LGM thinning. A forward model of nuclide concentrations produced under synthetic ice histories constrains the duration of the ice surface low stand and thickness of any residual ice cover. Low stand durations of ca. 300–9000 years, and residual ice thicknesses of up to 15 m can explain our data however the best fit is provided by a low stand lasting 300–3800 years and ice no thicker than ∼7 m. While we cannot exclude the possibility that ice surface change at our site is due to local factors, the lack of an obvious mechanism for such localised change, alongside the wider evidence for major ice-sheet change in the WSE during the Holocene, leads us to suggest that the thinning, and subsequent re-thickening demonstrated by our data reflects dynamic thinning associated with retreat of the grounding line inboard of present-day limits in the southern WSE. The areal extent of this retreat is greater than predicted by existing ice-sheet models simulations. Direct geological evidence for Holocene retreat-readvance is now observed in the Ross Sea (Kingslake et al., 2018; Venturelli et al., 2020, 2023), Amundsen Sea (Balco et al., 2023; Nichols et al., 2024) and Weddell Sea (this study) sectors of the WAIS. Retreat-readvance of the WAIS can explain a higher than present GMSL in the mid-late Holocene but the relative contributions (i.e., future vulnerabilities) of the main sectors of the WAIS have yet to be fully determined.

Finally, modelling studies have suggested that following excess retreat inboard of present limits in the WSE, ongoing glacio-isostatic adjustment led to re-grounding of ice rises, increasing buttressing and promoting grounding-line advance on a reverse bed-slope (Kingslake et al., 2018; Pittard et al., 2022). This mechanism has been proposed as a negative feedback on potential future sea-level rise from Antarctica (Gomez et al., 2015; Kachuck et al., 2020). Our data provides a testable scenario of retreat-readvance, along with potential timeframes involved, which should be used to test numerical ice-sheet models that include representations of this mechanism and are used to make predictions of future sea-level rise from Antarctica. To maximise the utility of data such as ours requires an approach to quantifying the magnitude of grounding-line retreat that explains point measurements of ice-sheet thinning at logistically constrained drill sites.