Two-million-year-old snapshots of atmospheric gases from Antarctic ice

Over the past eight hundred thousand years, glacial–interglacial cycles oscillated with a period of one hundred thousand years (‘100k world’1). Ice core and ocean sediment data have shown that atmospheric carbon dioxide, Antarctic temperature, deep ocean temperature, and global ice volume correlated strongly with each other in the 100k world2–6. Between about 2.8 and 1.2 million years ago, glacial cycles were smaller in magnitude and shorter in duration (‘40k world’7). Proxy data from deep-sea sediments suggest that the variability of atmospheric carbon dioxide in the 40k world was also lower than in the 100k world8–10, but we do not have direct observations of atmospheric greenhouse gases from this period. Here we report the recovery of stratigraphically discontinuous ice more than two million years old from the Allan Hills Blue Ice Area, East Antarctica. Concentrations of carbon dioxide and methane in ice core samples older than two million years have been altered by respiration, but some younger samples are pristine. The recovered ice cores extend direct observations of atmospheric carbon dioxide, methane, and Antarctic temperature (based on the deuterium/hydrogen isotope ratio δDice, a proxy for regional temperature) into the 40k world. All climate properties before eight hundred thousand years ago fall within the envelope of observations from continuous deep Antarctic ice cores that characterize the 100k world. However, the lowest measured carbon dioxide and methane concentrations and Antarctic temperature in the 40k world are well above glacial values from the past eight hundred thousand years. Our results confirm that the amplitudes of glacial–interglacial variations in atmospheric greenhouse gases and Antarctic climate were reduced in the 40k world, and that the transition from the 40k to the 100k world was accompanied by a decline in minimum carbon dioxide concentrations during glacial maxima. Analysis of two-million-year-old ice from Antarctica provides a direct comparison of atmospheric gas levels before and after the shift from glacial cycles of 100 thousand years to 40-thousand-year cycles around one million years ago.


Article
strong katabatic winds lead to the exhumation of old ice from depth to the surface 16 .In 2015-2016, two additional cores were drilled.The first extended the drilling to a depth of 147 m at the BIT-58 borehole (named ALHIC1503 hereafter), where approximately 1-Myr-old ice had previously been discovered.The second (ALHIC1502), drilled 56 m to the southwest, reached bedrock at 191 m.These ice cores are dated by measuring the 40 Ar deficit in trapped air relative to the modern atmosphere 17 ( 40 Ar atm ), and represent an archive of the Earth's climate extending well into the 40k world.A nearby blue-ice core (Allan Hills Site 27; S27) gives a continuous record from about 100 to 250 ka, and serves as a reference section for comparison with ALHIC1502 and ALHIC1503 18 .
Age-depth profiles of the ALHIC1502 and ALHIC1503 cores (Fig. 1) show that both sites contain 300-500-kyr-old ice, overlying a basal unit confined to approximately the bottom 30 m, with ages that range from more than 1 to 2.7 Myr.At ALHIC1503, the deepest sample lies within about 1 m of bedrock and has the oldest 40 Ar atm age, 2.7 ± 0.3 (1σ) Myr.However, stratigraphic disturbance in both ALHIC1502 and ALHIC1503 is evident from abrupt age discontinuities (Fig. 1), accompanied by large swings in δ 18 O of O 2 (δ 18 O atm ) and δD ice values 15 .In light of these complications, we interpret the paleoclimate records from the Allan Hills BIA as discrete 'snapshots'-discontinuous but accurate portraits-of the climate state, rather than continuous time series.

Climate properties in the old ice
Paleoclimate reconstructions from stratigraphically discontinuous ice sections face several challenges.First, there is the potential for sampling bias in the preserved climate states due to possible differences in accumulation rate and/or ice rheology.Second, diffusion in million-year-old ice may lead to smoothing of paleoclimate properties or proxies 19,20 , though this effect is likely to be minor in Allan Hills ice (see Methods).Third, discrete sampling might not capture the complete glacialinterglacial range of properties of interest.Overall, we estimate that 76 ± 14% (2σ) of the total atmospheric CO 2 variability in the 40k world is recovered by our CO 2 samples.This estimate is based upon our deductions that 90% of the true glacial-interglacial climate variability is preserved in the Allan Hills ice, and that 84 ± 15% (2σ) of CO 2 variability within the ice is captured by discrete samples, assuming a linear relationship between CO 2 and benthic δ 18 O (see Methods).
Measured climate proxies and ice core properties across time (Fig. 2) include proxy records of regional temperature (δD ice ), CH 4 , CO 2 , and global O 2 fractionation (δ 18 O atm ) 21 .Samples older than 800 kyr are assigned the age of the nearest 40 Ar atm -dated sample and then binned into two age groups: the MPT (800-1,200 kyr), and the 40k world (1.2-2.7 Myr).Within the 40k world bin, samples are further subdivided into three groups with average ages of 1.5, 2.0, and 2.7 Myr (see Methods and Fig. 2).
Improbably high CO 2 concentrations or δ 13 C values of CO 2 indicating contributions from respiratory CO 2 were present in all 2.7-Myr-old samples, common among 2.0-Myr-old samples, but absent from 1.5-Myrold samples (see Methods).We chose −7‰ as the cut-off value for δ 13 C and rejected all CO 2 samples below the shallowest δ 13 C sample with an isotopic value of lighter than −7‰.The application of these criteria excluded 21 CO 2 samples from our analysis.CH 4 and δ 18 O atm samples were also rejected when δ 13 C or CO 2 measurements indicated that these properties were compromised.Ar atm ages plotted against distance above bedrock with the same colour coding as in a, b.Data in circles were measured at Princeton University 15 .New measurements made at Scripps Institution of Oceanography are plotted as squares.Error bars represent the external reproducibility (1σ; 110 kyr and 220 kyr for samples measured at Scripps and Princeton, respectively) of the measurement or 10% of the sample age, whichever is greater.Shading highlights sections in which ice that is more than 1 Myr old is present.Sections affected by respiration are marked by the vertical arrows (blue, ALHIC1502; red, ALHIC1503). .c, The 800-kyr ice core CO 2 record 3,5,6,28 (purple line) and the binned Allan Hills CO 2 data (purple circles).Note that there are no reliable CO 2 and CH 4 analyses in the 2.7-Ma bin (see Methods).d, The 800-kyr ice core δ 18 O atm record 29,30 (brown line) and binned values of Allan Hills δ 18 O atm samples (brown circles).All δ 18 O atm values are normalized to the modern atmosphere.e, The globally distributed benthic oxygen isotope stack over the Plio-Pleistocene (LR04) 2 showing decreased glacial-interglacial variability and less extreme glacials before 800 ka.
Ice dating back to about 100−250 ka, sampled at a nearby site 18 , records a δD ice range of −280 to −360‰ in the 100k world.δD ice values ranged from −284 to −316‰ in 40k-world ice and from −288 to −332‰ in MPT ice.Maximum values of δD ice were similar during all three time periods, lying within an 8‰ range.However, minimum values were higher in the MPT and the 40k world by at least 28‰.This implies that interglacial periods during the 40k world in Antarctica were not significantly warmer than those in the 100k world, but that glacial maxima were lower after the MPT.
CO 2 concentrations ranged from 221 to 277 ppm (n = 23) 15 in MPT ice, and from 214 to 279 ppm (n = 37) in ice from the 40k world.Given that 76 ± 14% (2σ) of total CO 2 variability was captured by 37 discrete samples, and assuming that glacial and interglacial extremes are equally absent, the expected true range of 40k world CO 2 is 86 −13 +19 ppm (2σ; Fig. 3).Our best estimate of 40k-world CO 2 concentrations is therefore 204-289 ppm.
These ranges indicate that interglacial CO 2 concentrations in the 40k and 100k worlds were similar, whereas glacial CO 2 concentrations are likely to have been 24 ppm higher in the 40k world than in the 100k world.By comparison, reconstructions of atmospheric CO 2 from individual boron isotope measurements in foraminifera ranged from 190 to 320 ppm (Fig. 3), with the average extreme values being 234 and 277 ppm 9,10 .

Implications for climate evolution
Our findings appear to be inconsistent with hypotheses that attribute the transition into the 100k world to a long-term decline in both interglacial and glacial atmospheric CO 2 22,23 . Our data instead support hypotheses that link the MPT to greater ice sheet size and CO 2 drawdown during glacial maxima.Those hypotheses include, but are not limited to, enhanced dust delivery to the southern ocean 24 and changes in global ocean circulation 25 that resulted in an additional drawdown of atmospheric CO 2 , thereby enhancing the growth of continental ice sheets and leading to the emergence of 100-kyr glacial cycles.
Compared to CO 2 concentrations and δD ice , CH 4 concentrations and δ 18 O atm values in ice from the MPT and 40k world exhibited even less glacial-interglacial variability.CH 4 concentrations ranged from 405 to 569 parts per billion (ppb) in MPT ice (n = 31) and 457 to 592 ppb in 40k-world ice (n = 16).These ranges are only about 40% of the range in the 100k world (approximately 350 to 800 ppb) 26 .Undersampling may play a role in this reduction, as CH 4 maxima in the 100k world have short durations 26 .In addition, discrete CH 4 measurements require more ice than CO 2 measurements (60 g versus 10 g).They therefore average over a longer interval of time than CO 2 samples.δ 18 O atm values ranged from 0.01‰ to 0.92‰ in MPT ice (n = 63) and from 0.00‰ to 0.33‰ in 40k-world ice (n = 29).In comparison, the range is from −0.30‰ to +1.48‰ in ice from the 100k world 21 .The δ 18 O atm values in ice from the 40k world are remarkable in that they lack values higher than +0.5‰ and the total range of δ 18 O atm values is only about 20% of that observed in ice from the 100k world.Smaller changes in seawater δ 18 O values (Fig. 2) help to explain the absence of δ 18 O atm values higher than +0.5‰.However, additional mechanisms are likely to be needed to explain the lack of precession-modulated variability in δ 18 O atm associated with biosphere productivity and changes in the low-latitude water cycle 21 .
Plots of atmospheric CO 2 against δD ice and CH 4 against δD ice (Fig. 4) show that samples from the MPT and 40k world fall within the envelope of the 100k world but exhibit only a fraction of its range.Samples from the MPT and 40k world do not have the isotopically light δD ice and low CO 2 values characteristic of glacial samples in the 100k world.Both CO 2 and CH 4 are positively correlated with δD ice in ice from the 100k world (S27).By contrast, in the MPT and 40k world, only δD ice and CO 2 are significantly correlated (P = 0.002 and P = 0.003, respectively).For all three sample sets, the slopes of the correlation between CO 2 and δD ice are statistically indistinguishable (Fig. 4), in agreement with previous reconstructions from boron-based proxies and benthic δ 18 O values.Together, these observations suggest a consistent and persistent coupling between Antarctic climate and atmospheric CO 2 throughout the Pleistocene 8,9 , and that glacial cycles in the 40k world are truncated versions of glacial cycles in the 100k world.dashed lines represent the 95% confidence intervals) 9,10 .Shading represents the expanded range estimate, given the possibility that discrete Allan Hills CO 2 samples may not fully capture the true range of CO 2 variability.Black dashed lines represent the glacial-interglacial range of CO 2 in the 100k world 3 .

Conclusions
Shallow ice cores drilled in the Allan Hills BIA provide a direct observation of the variability in atmospheric CO 2 during the 40k world, confirming previous findings that minimum CO 2 concentrations declined after the MPT 8,9 .These results complement plans to drill a continuous ice core record back to around 1.5 Ma 27 .These snapshots of the Earth's climate system from shallow ice cores in BIAs provide an incomplete picture, especially with regards to the dynamics.Nevertheless, the oldest ice will probably be found in discontinuous sections.Our work demonstrates that BIAs can be exploited to extend climate records, including atmospheric CO 2 concentrations, well back into the early Pleistocene and perhaps even the Pliocene.

Sample location and description
Allan Hills Blue Ice Area (BIA).The Allan Hills BIA is located ~100 km to the northwest of the McMurdo Dry Valleys.Drilling sites are located to the west of the Allan Hills nunatak in the Main Ice Field (MIF) of the Allan Hills BIAs (Extended Data Fig. 1).Here, net ablation leads to the exposure of ancient glacial ice at the surface of an ice sheet 31 .Elevation, bedrock topography, and ice velocities in the MIF have been previously documented 32,33 .Site meteorology and surface snow properties are discussed in more detail by Dadic et al. 34 and ice and gas properties described in Higgins et al. 15 .
Site ALHIC1503 (previously named BIT-58).Site ALHIC1503 (76.73243° S, 159.3562°E) is located near a crest of a northwest-southeast trending ice ridge off the MIF.The precise direction and the velocity of ice flow are not known at this location.Measurements from the nearby ice field are consistent with primary flow to the east or northeast.Complex dust bands are visible from high-resolution satellite imagery (Extended Data Fig. 1, left), implying a complex history of glacial flow.Typical ice flow velocities near ALHIC1503 are low (0.015 m yr −1 ) 33 .Bedrock topography, determined using ground penetrating radar (GPR), indicates that ALHIC1503 is located on a steep (~45°) slope leading to a local bedrock high with an overlying ice thickness of ~150 m (Extended Data Fig. 2).In the 2015-2016 field season, 23 m of ice was drilled, extending the earlier 124-m-long ice core 15 to the bedrock at 147 m below the ice surface.
Site ALHIC1502.Site ALHIC1502 (76.73286° S, 159.35507°E) is located 56 m upflow from ALHIC1503.The ice thickness is ~204 m, consistent with a steep bedrock slope (Extended Data Fig. 2).Drilling at site ALHIC1502 recovered 197 m of ice.
Site 27 (S27).S27 (76.70° S, 159.31°E) is located along the main ice flow line of the Allan Hills region.Two hundred and twenty-four metres of ice was retrieved, representing ~70% of the estimated column thickness.Spaulding et al. 18 provide a comprehensive suite of stable water isotope and gas analyses for S27.The S27 ice chronology has been established by matching δD ice features to those of the Vostok ice core.This timescale was validated by correlating variations in δ 18 O atm (δ 18 O of paleoatmospheric O 2 ) into the record of δ 18 O atm in Vostok ice cores 35 .

Analytical procedures 40
Ar atm and δXe/Kr. 40Ar is produced in the solid earth by the radioactive decay of 40 K.As 40 Ar slowly leaks into the atmosphere, its concentration increases with time. 38Ar and 36 Ar, on the other hand, are stable, primordial, and non-radiogenic, so their atmospheric concentrations are treated as constant.The ratio of 40 Ar/ 38 Ar thus rises towards the future and decreases towards the past.The term of merit here is the gravitationally corrected paleoatmospheric 40 Ar/ 38 Ar ratio: 40 Ar atm = δ 40 Ar/ 38 Ar − δ 38 Ar/ 36 Ar.
The latter term corrects for the gravitational fractionation in the firn column 36 .Studies of 40 Ar atm as a function of age in the Dome C and Vostok ice cores constrained the rate of change to be 0.066 ± 0.006‰ Myr −1 , allowing us to date the air trapped in the ice without the prerequisite for stratigraphic continuity 17 .
A potential complication of 40 Ar atm dating comes from the in situ production of 40 Ar from the bedrock and/or from dust, a tiny amount of which is present in polar ice cores.The deepest sample at ALHIC1502 has a much younger age (0.8 Myr) than the overlying ice (2.2 Myr).This inverted age-depth relationship may result either from folding or from radiogenic input of 40 Ar from the 40 K in the underlying bedrock, a phenomenon observed for basal ice at GISP2 in Greenland 37 .The 40 K from dust in the ice would produce a negligible amount of 40 Ar.
Ar isotope analyses on the original BIT-58 core between the surface and 126-m depth were carried out at Princeton University with a procedure modified from Bender et al. 17 and described by Higgins et al. 15 .The standard deviation of 40 Ar atm measured across 68 external standards is ±0.014‰, corresponding to an age uncertainty of ±210 kyr (1σ).For Allan Hills samples drilled in 2015-2016, 40 Ar atm was measured at Scripps Institution of Oceanography with slight modification.First, sample size was increased from 500 to 800 g for better precision.Second, pumping time for evacuating the ice-bearing vessel before gas extraction was reduced to 15 min to prevent gas loss.Third, extracted and purified gases were measured on a newer mass spectrometer, which offers better signal sensitivity and stability.Finally, following published procedures 38 , we measured the Xe/Kr ratios, expressed as δXe/Kr, an indicator of mean ocean temperature 39 , in the same aliquot of the extracted and getterred gases.A typical 40 Ar atm and δXe/Kr measurement takes ~3.5 h on the mass spectrometer.
Two measurement campaigns at Scripps Institution of Oceanography were carried out and each was associated with slightly different analytical procedures.In the first campaign, a fixed amount of reference gas was introduced into the bellow of the mass spectrometer, leaving different amount of gas in sample and reference sides.As the sample and reference gases deplete during analysis, the pressure and ion current fall more rapidly for the side containing less gas.As a result, 40 Ar atm was observed to depend on the volume of the analyte relative to the reference.All measured 40 Ar atm values were subsequently corrected for volume differences, using the voltage readings of the mass-40 ( 40 Ar) beam when the gases were introduced to the fully expanded bellow on the mass spectrometer.The magnitude of this correction is generally less than 0.02‰, or 300 kyr, and the uncertainty in the correction is much smaller.In the second campaign, the amount of gas in the standard-side bellows was matched to that of the reference-side bellows.No volume correction was applied in the second campaign.Based on replicate analyses, the external reproducibility was ±0.006‰ (1σ) and ±0.007‰ (1σ) for dry La Jolla air samples and Taylor Glacier blue-ice samples, respectively.The samples measured at Scripps therefore have an age uncertainty of ±110 kyr (1σ) owing to the analytical uncertainty, or 10% of the sample age owing to the uncertainty in the 40 Ar outgassing rate, whichever is greater.Individual sample uncertainties are reported along with the measured values.
Xenon-to-krypton ratios (δXe/Kr) were measured immediately after the 40 Ar atm measurements on the same mass spectrometer by 'peak hopping' between 84 Kr and 132 Xe.The δXe/Kr measurement takes about half an hour to complete.The final reproducibility is ±0.21‰ and ±0.27‰ for dry La Jolla air samples and Taylor Glacier ice replicates, respectively.In order to calculate mean ocean temperature from δXe/Kr, which is normalized to La Jolla air, we need to make an assumption about the volume of the ocean 39 .For this purpose, we use a sensitivity of 1.1 ± 0.2 °C/‰.We obtain this value from the observed last glacial maximum (LGM)-present δXe/Kr change 39 with minor modifications.We adopted all the same sealevel-related corrections as in ref. 39 except one: the assumption that Kr and Xe saturation state was higher during the LGM than at present.This assumption lowers the sensitivity of mean ocean temperature (MOT) to δXe/Kr from 1.2 to 1.0 °C/‰ (for example, a MOT change of 2.57 ± 0.24 °C for a 2.5‰ change of δXe/Kr 39 ).To account for this uncertainty we adopt an intermediate value of 1.1 and a larger error of ±0.2.δD ice and δ 18 O ice .Stable water isotope measurements (δ 18 O ice , δD ice ) were taken from a slab of ice, 2.5-cm wide, cut (parallel to the vertical axis) from the outside of each ice core.Because some portions of the core were fractured, continuous sections were not available at all depths.These slabs were sub-sampled at 10-15-cm resolution at the Climate Change Institute, University of Maine.Each sub-sample was melted in a sealed plastic bag at room temperature, vigorously shaken, and then decanted into an 18-ml plastic scintillation vial.All vials were refrozen and stored below −10 °C until the time of analysis.Water isotope samples were measured via cavity ring-down spectroscopy (CRDS) using a Picarro Model L2130-i Ultra High-Precision Isotopic Water Analyzer coupled with a High Precision Vaporizer and liquid auto sampler module.δ 18 O ice and δD ice were measured simultaneously with an internal precision (2σ) of ±0.05‰ for δ 18 O ice , and ±0.10‰ for δD ice .The instrument drift adjustment and daily calibration were accomplished by periodically measuring three secondary standards: BBB (an average Maine freshwater from the Bear Brook watershed), ASS (Antarctic surface snow) and LAP (light Antarctic precipitation).These standards were calibrated against the internationally accepted water isotope standards: SMOW (standard mean ocean water), SLAP (standard light Antarctic precipitation) and GISP (Greenland ice sheet precipitation).δ 15 N, δ 18 O atm , and δO 2 /N 2 .Analyses for the original BIT-58 O 2 /N 2 /Ar elemental ratio, δ 18 O of O 2 , and δ 15 N of N 2 were carried out as described 15,30 .Samples from ALHIC1502 and ALHIC1503 were measured using a slightly modified procedure.In brief, ~20 g of ice was placed into a glass flask chilled in a dry-ice-isopropanol cold bath (−78.5 °C) and the vessel was turbo-pumped for 3 min.The final pressure of the residual gas inside the vial was no higher than 1.5 mTorr.The pressure of the air trapped in ice, when expanded to the same space, was ~3 Torr.The reduced pumping time thus introduces, at most, about 0.05% of modern atmosphere into the sample gas.
After melting, the gas-meltwater mixture was equilibrated for at least 4 h.The majority of the meltwater was drained and leftover water inside the glass vial refrozen at −30 °C.Once all meltwater had refrozen, headspace gases were collected by condensation at liquid helium temperatures.The sample was admitted to a mass spectrometer (Thermo Finnegan Delta Plus XP) for elemental and isotope ratio analysis after homogenizing at room temperature for 30 min.The reproducibility of δ 15 N, δO 2 /N 2 and δ 18 O atm values from ALHIC1502 and ALHIC1503 samples analysed using this procedure, calculated as pooled standard deviation, is ±0.008‰, ±3.164‰, and ±0.024‰, respectively.
Overall, the reduced pumping time leads to a reduction in gas loss and an improvement in precision.We note that the gas loss processes that alter δO 2 /N 2 did not appear to affect the isotopic composition of O 2 and N 2 in our samples (Extended Data Fig. 3).As a result, no gas loss correction was applied to the gas isotopes, in contrast to some previous studies in heavily fractured ice (for example, Siple Dome 40 ).All δ 15 N, δO 2 / N 2 and δ 18 O atm values were normalized to modern atmosphere.CH 4 , CO 2 and δ 13 C of CO 2 .CH 4 concentrations were analysed at Oregon State University using a melt refreeze technique 41 .Samples (~60-70 g of ice) were trimmed, melted under vacuum, and refrozen at about −70 °C.CH 4 concentrations in released air were measured by gas chromatography and referenced to air standards calibrated by National Oceanic and Atmospheric Administration Global Monitoring Division (NOAA GMD) on the NOAA04 scale.Precision was generally better than ±4 ppb.Individual sample uncertainties, where available, are reported along with the measured values.
CO 2 concentrations, referenced to air standards calibrated by NOAA GMD on the World Meteorological Organization (WMO) scale, were measured using a dry extraction (crushing) method 42 .Wherever possible, ice samples were processed and analysed in replicate for each depth and results averaged to obtain final CO 2 concentrations.However, four pairs of replicates did not come from exactly the same depth and their reproducibility showed substantial deterioration beyond usual analytical uncertainties.In this case, we treated those four pairs of replicates as eight individual, unreplicated data points.For greenhouse gas concentrations, no gravitational fractionation correction was applied; the correction would be less than 0.8 ppb for CH 4 and 0.5 ppm for CO 2 .
The measurement protocol for δ 13 C of CO 2 was as described 43 with improvements in the efficiency of air extraction.Approximately 200 g of ice was dry-extracted with a custom-made 'ice grater'.The reproducibility of the method, based on replicate analysis of Taylor Glacier blue-ice samples, is 0.02‰ for δ 13 C-CO 2 .The method also provides a measurement of CO 2 concentrations.The uncertainties are ±2 ppm (1σ of the pooled standard deviation).The final δ 13 C values, normalized to the Vienna-Pee Dee Belemnite (VPDB) standard, were corrected for blanks, N 2 O + mass interference from N 2 O, and gravitational fractionation using the nearest δ 15 N value.
Age assignment of CO 2 /CH 4 and O 2 /N 2 /Ar samples The depths of samples analysed for CO 2 /CH 4 and O 2 /N 2 /Ar were different from depths of the 40 Ar atm samples that provided the chronology.CO 2 /CH 4 samples, O 2 /N 2 /Ar samples, and ice samples were assigned ages by bracketing 40 Ar atm values.Here we discuss two different models for assigning ages and show that different criteria for age assignment lead to similar conclusions.For the purpose of investigating the relationship between greenhouse gases (CO 2 and CH 4 ) and Antarctic climate (δD ice ), we binned samples into three units: MPT, 800-1,200 ka; pre-MPT, >1,200 ka; and young ice, <800 ka.We focus on the robustness of ages for CO 2 and CH 4 samples.All CO 2 and CH 4 data discussed below exclude samples affected by respiration.
Conservative approach.In the most conservative approach, an age was assigned only to samples bracketed by two 40 Ar atm ages lying within a single age bin.Otherwise, the sample was left 'unassigned'.In this binning scenario, 19 CO 2 samples were classified as MPT, 15 as pre-MPT, and 34 as unassigned.We classified 28 and 19 CH 4 samples as MPT and unassigned, respectively.Eight CH 4 samples were binned into pre-MPT.
Proximity approach.In this approach, each CO 2 /CH 4 datum was assigned the 40 Ar atm age of the closest Ar isotope sample.In this age model, 23 and 37 CO 2 samples were binned into the MPT and pre-MPT units, respectively.For CH 4 , 31 samples were classified as MPT, and 16 samples as pre-MPT.Finally, eight CO 2 and eight CH 4 samples were binned into the younger (<800 ka) age unit.
Extended Data Fig. 4 summarizes the CO 2 -δD ice and CH 4 -δD ice relationships using the two age models.Our conclusions regarding the greenhouse gas-δD ice relationship are not dependent on which model is used.As the proximity approach assigns an age to every depth, we adopted it for our interpretations.

Assessing the fraction of climate variability recovered by discrete sampling
When interpreting data from stratigraphically disturbed ice, the question arises as to what fraction of the amplitude of climate variability is accessed in the suite of samples analysed.Here we consider two factors: (1) the fraction of climate variability preserved in the ice cores compared to the true natural variability; and (2) to what extent discrete samples capture the variability that is preserved in the ice.Below we discuss these two factors separately.
Fraction of the pre-MPT δD ice range preserved in the ice.We start by considering what fraction of the true glacial-interglacial variability in the 40k world is preserved in the ice.If we assume the absence of post-depositional alteration processes such as respiration and diffusive smoothing, this fraction should be an intrinsic property of the ice itself and insensitive to the property being measured.Given the large number (>500) and the high spatial resolution (10 cm) of δD ice measurements, we used δD ice to evaluate the fraction of true δD ice variability preserved in the ice.Next, we needed to make an assumption of the true δD ice amplitude in the 40k world.It is assumed that the true amplitude of δD ice variations in Allan Hills ice scales linearly with the amplitude of deep ocean temperature.This assumption is justified for two reasons.First, most of the bottom water formation today occurs in the southern ocean.Second, Antarctic temperature has been found to correlate tightly with mean ocean temperature during the LGM-Holocene transition 39 and with South Pacific bottom water temperature through the past 800 kyr 44 .
For deep ocean temperature, we considered two scenarios.First, we assumed a constant partitioning between ocean temperature and the influence of seawater δ 18 O on the benthic foram δ 18 O record throughout the Pleistocene.This means that the amplitude of δD ice should simply scale with the benthic foram δ 18 O values.The amplitudes of glacial cycles in the early and mid-Pleistocene (between marine isotope stages (MIS) 39 and 104), expressed in benthic foram δ 18 O, averages 0.70‰, representing 37% of the change from MIS 6 to MIS 5e (1.90‰).
The range of δD ice values between MIS 6 and MIS 5e observed in the Allan Hills ice from the 100k world was 80‰.Multiplying this 100k-world range by the scaling factor (37%) based on benthic foram δ 18 O predicts the range of δD ice in 40k-world ice to be 30‰.By comparison, the range of 40k-world δD ice observed in the Allan Hills ice is 32‰.Using this approach, we estimated that our pre-MPT ice samples represent ~100% of the average amplitude of glacial-interglacial δD ice cycles in the 40k world.
In an alternative scenario, we considered deep ocean temperature inferred from a seawater δ 18 O record constrained independently from records of carbonate microfossil δ 18 O from the Mediterranean Sea 45 .A slightly higher variability in the 40k world deep ocean temperature, which corresponds to 50% of the average variability in the 100k world, lowers our estimate of the recovered climate variability to ~80%.Consequently, although other factors might influence the recovery of the range of climate variability, given these two approaches our best estimate is that the Allan Hills ice samples presented here preserve 90% of the true range of the climate variables in the 40k world.
Fraction of the pre-MPT CO 2 and CH 4 range recovered by discrete samples.CO 2 and CH 4 concentrations were measured on a relatively small number of samples.It is therefore possible that we have not retrieved the full range of variability of these properties in the 40k world, even if their true variability is well preserved in the ice cores.We evaluated this possibility by first generating a synthetic CO 2 series based upon the assumption that the variation of CO 2 in the 40k world scales linearly with benthic foram δ 18 O values.A synthetic CO 2 time series between 1.2 and 2.0 Ma was constructed on the basis of the LR04 global benthic foram δ 18 O stack 2 .A linear regression between benthic δ 18 O and atmospheric CO 2 between 0 and 800 ka resulted in correlation parameters that were used to create an artificial atmospheric CO 2 record further back in time (Extended Data Fig. 5a, b).After resampling to interpolate the synthetic data at a temporal resolution of 2.5 kyr, a Monte Carlo simulation randomly selected 37 samples from the synthetic CO 2 record.The simulation was repeated 10,000 times.We then calculated the observed range of the 37 random CO 2 samples and compared it to the 'true' range, yielding a ratio for each simulation run.Finally, the distribution of this ratio is shown as a series of histograms (Extended Data Fig. 5c-h).Note that the term of merit for each simulation is the ratio of the sampled range to the hypothetical true range.
To account for the fact that interglacial accumulation rates are generally higher than glacial rates 46 , we multiplied the modelled occurrence of interglacial ice (CO 2 > 250 ppm) by 2, 4, 8, 16, and 32.When the presence of the interglacial ice is two times the glacial ice, 84 ± 15% (95% confidence interval) of the total CO 2 range can be captured by 37 samples (Extended Data Fig. 5d).This ratio was derived from Vostok, an inland coring site of East Antarctica 47 .
We note that 35 out of the 37 40k-world CO 2 samples involve replicates measured at the same depth.The assumption here is that replicates at the same depth have exactly the same age, which remains unevaluated in Allan Hills ice cores.It is possible that the ice flow and dip of ice layers make samples from the same depth different in age.If this is the case, each pair of replicates should be treated as two individual measurements, and the total number of samples would be 72.In this scenario, the recovered range is from 214 to 281 ppm, and 72 discrete samples capture 88 ± 13% of the CO 2 variability preserved in the ice.Treating replicates as individual samples does not significantly change the range (from 65 to 67 ppm), but it increases the likelihood that we are recovering close to the full range of CO 2 in the ice (from 84% to 88%).
Considering the discussion above, our best estimate is that the 37 CO 2 samples cover 76 ± 14% of the true CO 2 range in the 40k world (90% of true variability preserved in the ice multiplied by 84 ± 15% that is captured by discrete sampling).Furthermore, we note that the δD ice values of the ice samples analysed for CO 2 as well as CH 4 span more than 90% of the observed 40k-world δD ice range (Fig. 4).This number is higher than the estimate yielded by Monte Carlo simulations, meaning that the discrete samples are likely to capture more variability than probability alone would dictate.As a result, there is a greater chance that we have captured a significant fraction (>90%) of the glacial-interglacial CO 2 range in the 40k world.
For δ 18 O atm values the estimated recoverability should be no less than 70%.First, 29 samples were analysed for δ 18 O atm .Second, the co-depth δD ice values of δ 18 O atm samples span 97% of the 40k-world δD ice (Extended Data Fig. 6a).

Potential sample mixing on various length scales
Mixing by molecular diffusion.The interval between 123 and 141 m in ALHIC1503 shows unexpectedly low variability in the elemental and isotopic ratios of the major gases (Extended Data Fig. 7).Molecular diffusion in the firn column has been shown to reduce annual variability of water isotopes 48 and gases 49 .In polycrystalline ice, diffusive mixing has also been invoked to explain the lack of sub-millennial δD ice variability preserved in the EPICA Dome C ice core, dated ~780 kyr old 20 .When annual layers are thin and temperatures are high, some gas records (for example, δO 2 /N 2 ) can also undergo extensive diffusive smoothing and are expected to lose even the glacial-interglacial variability 20 .
However, diffusive smoothing of climate records is unlikely to be significant at the Allan Hills BIA.First, current temperatures within the boreholes at the Allan Hills BIA are around −30 °C 32 .Rates of diffusion at the Allan Hills BIA should be substantially slower than at the warm basal temperatures (>−10 °C) expected at the bottom of thick (>3 km) polar ice sheets.Second, as is shown below, water isotopes (δD ice ) are expected to be much more sensitive to diffusive smoothing than gases (O 2 and N 2 ).As a result, the presence of substantial sample-to-sample variability in δD ice values measured at 10-cm resolution indicates that the effects of diffusion on this and longer length scales are likely to be minor.
Below we quantitatively estimate the characteristic diffusion time scale of water isotopes and of oxygen and nitrogen gases in the ice.Gas and water permeation coefficients in the ice have been estimated from molecular dynamics simulations 50,51 and limited observations in some 'natural experiments' 52,53 .Values obtained from different methods can differ by more than one order of magnitude.In our calculation, the fastest permeation parameters are used whenever possible to provide the most rigorous test of diffusive smoothing.
For the self-diffusion of water molecules in the ice, the diffusion time scale (τ ice ) depends solely on the diffusivity of water molecules in ice (D ice ) and the length of interest (L; 0.1 m): The diffusivity of ice is determined by 54 : where D ice 0 (9.1 × 10 −4 m 2 s −1 ) is a constant of diffusion for water molecules in ice, Q ice (5.86 × 10 4 J mol −1 ) the activation energy of diffusion, R (8.314 J mol −1 K −1 ) the gas constant, and T the temperature.The characteristic time scale of self-diffusion in ice as a function of temperature is shown in Extended Data Fig. 8. Depending on temperature, it takes 10 5 to 10 7 years for diffusive smoothing to smear the original signal on the length scale of 0.1 m.Next, we consider the molecular diffusion of O 2 and N 2 in the ice.Below is a simple model of gas permeation with a large reservoir of gases (that is, bubbles) and a medium through which diffusion takes place (that is, ice): (1) First, gases in the reservoir dissolve into the ice.The amount of dissolved gas depends on ice volume, pressure, and solubility.
(2) Next the gas will diffuse in the ice driven by the concentration gradient.This is assumed to be the rate-limiting step.(3) Diffusive gas exchange in the ice along the concentration gradient will immediately be compensated by the gas exchange between the ice phase and the bubble phase.(4) Eventually the concentration gradient in the bubbles will be balanced owing to diffusive flux.
Were there no gas reservoirs, the characteristic diffusion time scale of gas species m would simply follow the diffusivity of gas m in ice (D m ) and the length of interest (L): where L is 0.5 m because the spatial resolution of gas measurements is about 50 cm.Similar to the diffusivity of water isotopes, the diffusivity of gas m follows 50 : where D m 0 is a constant and Q m,ice is the activation energy of diffusion.
However, the presence of the reservoir would compensate for the diffusive flux.Here, we introduce the concept of partitioning function Z, the ratio of the number of gas molecules in the bubbles to the number of gases dissolved in the ice.
By definition, Z of a given gas species m is where n m,gas and n m,ice are the number of gas molecules in the gas phase and in the ice phase, respectively.It is assumed that n m,gas ≫ n m,ice .The number of gas m molecules in the gas phase can be inferred from the molar fraction of gas m (f m ) and the total number of gases in the ice (n gas ), which has indeed been measured as the gas content of the ice (V), expressed in SI units m 3 kg −1 : m m ,gas gas Given the ideal gas law, equation ( 7) can be re-written as: where m ice is the mass of the ice, p 0 is 101,325 Pa and T 0 is 273 K.The number of molecules of gas m in the ice phase depends on the solubility of gas m (S m ), the partial pressure of gas m (P m ), and the number of water molecules (n ice ): The solubility of gas m is governed by: and Q N 2 are 9,200 J mol −1 and 7,900 J mol −1 , respectively 50 .
The partial pressure of gas m is the product of the bubble pressure, which is assumed to be the hydrostatic pressure of the depth in which bubbles are located, and the molar fraction of the gas m in air: where ρ ice is the density of ice (920 kg m −3 ), g the gravitational acceleration constant (9.8 m s −2 ), and h the depth.In our case, h is assumed to be a constant 150 m.
The number of water molecules n ice can be directly calculated from the mass of the ice: ice ice water where M water is the molecular weight of H 2 O (0.018 kg mol −1 ).Therefore, the final expression for Z m is: Using equations ( 4) and ( 14) and relevant parameters, the diffusion time scale of O 2 and N 2 are computed as a function of absolute temperature (Extended Data Fig. 8).Unless at very cold temperatures (<230 K), water isotopes measured every 10 cm should be more susceptible to diffusive mixing than gases measured every 50 cm.As a result, the preservation of δD ice variability argues that flat δO 2 /N 2 and δ 18 O atm profiles are not the result of diffusive smoothing.However, if the temperature is indeed that low, the characteristic time scale of gases becomes too large (>10 7 yr) for substantial smoothing to occur in million-year-old ice.We did not calculate the permeation coefficients for CH 4 or CO 2 as these molecules have larger molecular diameters than O 2 and N 2 , and should diffuse at an even slower rate and have a longer τ (ref. 19).
Mixing within individual samples.In polar ice sheets the thickness of annual ice layers declines with depth owing to compression and lateral flow.Thinning will reduce geochemical variability if a property changes over timescales shorter than that encompassed by the physical length of a single sample.This problem is further amplified in stratigraphically disturbed ice, as the orientation of the layering is uncertain.CO 2 and O 2 /N 2 / Ar samples were cut from smaller lengths of ice than δD ice samples (which were 10-cm long).Therefore, assuming the layering is horizontal, CO 2 and O 2 /N 2 /Ar samples should average shorter time intervals than is the case for δD ice samples.On the other hand, CH 4 samples were cut in 10-cm lengths.Because δD ice samples show large variability, we consider that gas properties other than 40 Ar atm and δXe/Kr represent relatively short climate intervals without much mixing.The δD ice samples themselves may sometimes be mixed, and this possibility remains to be examined by the analysis of samples smaller than 10 cm.
Flow-induced mixing.Physical mixing of ice of different ages can result from glacial flow.Mechanical mixing of this type may create an artificial correlation between properties that are otherwise uncorrelated.Granted, many climate properties will co-vary in the absence of mixing.Therefore, property-property correlations do not necessarily indicate that mixing has occurred.Extended Data Fig. 9a shows the Pearson correlation coefficient (expressed as R 2 ) between δD ice and d within intervals of varying lengths.In normally ordered ice cores, d is found to have a complex relationship with δD ice 55 .We calculate R 2 as a function of the number of δD ice samples included.The section between 132 and 142 m shows a very strong correlation between δD ice and d, which raises the suspicion of mixing.
We next examined the gas records between 132 and 140 m.The maximum and minimum CO 2 values within this 10-m section are 253 and 217 ppm, respectively.Because the gas content in the deep ALHIC1502 and ALHIC1503 samples is very similar (~0.098 cm 3 g −1 ), a simple linear mixing curve in the gas phase is expected.If we regard these two points as end members and plot δD ice against CO 2 , the data points do not fall onto a mixing line (Extended Data Fig. 9b).Similarly, it is difficult to explain the δD ice -CH 4 and CO 2 -CH 4 plots by simple two-end-member mixing.In other words, a simple mixing scenario with two end members cannot explain our gas data, measured at every ~50 cm.We therefore consider the property-property plots (Fig. 4) to be unaffected by mechanical mixing on the length scale of 50 cm, although mixing at a much smaller length scale is still possible.

Effect of respiration on gas properties
Respiration of detrital organic matter has been shown to lead to the in situ production of CO 2 in the bottom of polar ice sheets 56,57 , resulting in elevated concentrations of greenhouse gases compared to the contemporaneous atmosphere.We observed substantially elevated CO 2 and CH 4 concentrations in the basal sections of ALHIC1503 and ALHIC1502, which we attribute to respiration and methanogenesis, respectively.Two lines of evidence support this hypothesis.First, there are extremely depleted (<−500‰) δO 2 /N 2 ratios and high (>4‰) δ 18 O atm values in samples near the bottom of ALHIC1503.For comparison, typical δO 2 /N 2 values in Antarctic ice cores vary between −5 to −20‰ 58 and δ 18 O atm ranges from −0.5 to +1.5‰ 21 on glacial-interglacial timescales.Second, measured δ 13 C of CO 2 values in some deep samples from the ALHIC1502 and ALHIC1503 cores (Extended Data Fig. 10) are outside the range (−6 to −7‰) observed for the last 160 kyr 59 .This is consistent with contributions from respiratory CO 2 with a δ 13 C value of about −25‰.The deepest measured δ 13 C sample comes from 3 m above the bedrock in ALHIC1503 and has a δ 13 C value of −22.4‰.Assuming an initial δ 13 C of CO 2 value of −6.5‰ and an δ 13 C of CO 2 value of −25‰ from respired carbon, isotope mass balance suggests that 86% of the CO 2 in this sample is derived from respiration.
In four out of the nine δ 13 C samples, δ 13 C of CO 2 was compatible with the absence of respiration, assuming a cut-off δ 13 C value of −7‰.In the other five samples, which originated within 7 m of the bottom of the two cores, however, δ 13 C of CO 2 was lighter than −7‰ and indicates respiration.We marked the shallowest δ 13 C sample with less than −7‰ isotopic as a cut-off depth and rejected all CO 2 samples below it.Admittedly it is still possible that samples immediately above those cut-off depths might also be affected, but the highest pre-MPT CO 2 value is found around 131 m in ALHIC1503, a depth that is bracketed by two δ 13 C measurements falling between −6% and −7% and is considered respiration-free.As a result, even if certain samples above the threshold depths are contaminated by respiration, our conclusions would not change.On the basis of these results, along with anomalous values of CH 4 , δO 2 /N 2 , and δ 18 O atm in the deeper sections, we excluded biogenic gas data below 185.950 m in ALHIC1502 and 139.625 m in ALHIC1503 from our evaluation.We still included δD ice from those sections in Fig. 2.   5b) by 37 randomly selected samples, given different preferential preservation of interglacial ice.Here the term 'interglacial' is operationally defined as any sample with greater than 250 ppm CO 2 .Different multipliers indicate the multiple occurrence of the 'interglacial' ice to simulate varying degrees of ice preservation biases.This analysis shows that when the presence of interglacial ice is no more than eight times greater than that of glacial ice (c-f), we could expect 37 random samples to be likely to capture 66-98% (95% confidence interval) of the CO 2 range in the ice.

Fig. 1 |
Fig.1| Age-depth profile of Allan Hills ice cores.a, ALHIC1502 data in blue; b, ALHIC1503 data in red.c,40 Ar atm ages plotted against distance above bedrock with the same colour coding as in a, b.Data in circles were measured at Princeton University15 .New measurements made at Scripps Institution of Oceanography are plotted as squares.Error bars represent the external reproducibility (1σ; 110 kyr and 220 kyr for samples measured at Scripps and Princeton, respectively) of the measurement or 10% of the sample age, whichever is greater.Shading highlights sections in which ice that is more than 1 Myr old is present.Sections affected by respiration are marked by the vertical arrows (blue, ALHIC1502; red, ALHIC1503).

Fig. 2 |
Fig.2| Climate properties over the past 2.9 Myr documented in ice core and benthic foram records.a, The continuous 800-kyr δD ice record from the European Project for Ice Coring in Antarctica, Dome C (EDC) ice core 4 (black line); the continuous S27 δD ice record covering 120-250 ka18 (red line); and the discrete Allan Hills δD ice ice samples (red circles).b, The 800-kyr ice core CH 4 record26 (orange line) and the binned Allan Hills CH 4 data (orange circle).c, The 800-kyr ice core CO 2 record 3,5,6,28 (purple line) and the binned Allan Hills CO 2 data (purple circles).Note that there are no reliable CO 2 and CH 4 analyses in the 2.7-Ma bin (see Methods).d, The 800-kyr ice core δ18 O atm record29,30 (brown line) and binned values of Allan Hills δ18 O atm samples (brown circles).All δ18 O atm values are normalized to the modern atmosphere.e, The globally distributed benthic oxygen isotope stack over the Plio-Pleistocene (LR04)2 showing decreased glacial-interglacial variability and less extreme glacials before 800 ka.

Fig. 3 |
Fig. 3 | Comparison of blue-ice CO 2 record and boron-based CO 2 reconstructions during and before the MPT.a, LR04 benthic oxygen isotope stack 2 .b, Distribution of MPT CO 2 (left, orange) and pre-MPT CO 2 (right, red) data observed in Allan Hills ice.c, Comparison between ice core CO 2 (orange and red circles) and δ 11 B-based CO 2 reconstructions (purple and blue solid lines; −1  ] is the dissolution constant, and Q m the activation energy of dissolution.10 −13 Pa −1 and 4.5 × 10 −13 Pa −1 , respectively; Q O 2

Fig. 1 |
Satellite imagery of the Allan Hills study area.Right, WorldView03 colour pan-sharpened imagery (copyright 2011, DigitalGlobe, Inc.) of the main ice field, Allan Hills Blue Ice Area (in true colour mode).Inset, Antarctica with hill shading (source: Australian Antarctic Data Centre, Map 13469; licensed under a Creative Commons Attribution 3.0 Unported License) on which our study area is marked by the black square.The black outcrop in the main image is the Allan Hills nunatak.The red arrow marks the local iceflow.Left, a magnified image of the drilling site from same source file.The imagery is enhanced (gamma-adjusted on ArcGIS 10) to highlight the colour contrast within the blue ice, which now has a greenish hue owing to the colour rendition.The brown lineations in the ice are exposed dust bands, providing a first-order tracer of surface ice stratigraphy.The locations of the cores reported in this work (see text) are marked with red circles.The location of the GPR profile in Extended Data Fig. 2 is shown as a yellow dashed line.Extended Data Fig. 2 | GPR profile proceeding from the SW to the NE, crossing (within less than 5 m) ALHIC1502 and ALHIC1503.The exact transect location is shown in Extended Data Fig. 1.The location and depths of boreholes ALHIC1502 and ALHIC1503 (red bars) and ice drilled previously as BIT58 (grey bar) are indicated.The profile was collected with an 80-MHz MLF antenna in a step-and-collect survey style with a step size of 25 cm and a stacking rate of 64 scans.Standard post-processing steps were applied, including time-zero position correction, distance normalization, 50/110-MHz finite impulse response filter, background removal, and gain adjustment.Depth estimates from two-way travel time (TWTT) and migration use a radar travel velocity of 0.165 m ns −1 .a, Non-migrated radargram showing dipping bed-parallel englacial stratigraphy and a strong dipping apparent bedrock reflection.Dashed blue line indicates the modelled location of the actual bedrock location as shown in b. b, Migrated radargram showing the 'true' location of the bedrock.Migration quality is poor at this location owing to the steeply dipping slope of the bedrock rise.Comparison between a and b suggests that the correct bedrock reflector location is translated downward and more steeply dipping than indicated in the non-migrated data.ALHIC1502 and ALHIC1503 borehole depths independently verify this interpretation.Extended Data Fig. 3 | Minimal impact of gas loss on oxygen and nitrogen isotopes.Pair differences of δ 18 O atm (Δδ 18 O atm ; red) and δ 15 N (Δδ 15 N; blue) are plotted against ΔδO 2 /N 2 in ALHIC1502 and ALHIC1503 samples.The difference between two replicate samples is attributable to gas loss.The very weak dependence of Δδ 18 O atm and Δδ 15 N on ΔδO 2 /N 2 suggests that gas isotopes are relatively insensitive to gas loss fractionation in Allan Hills ice.Extended Data Fig. 4 | Age assignment of CO 2 and CH 4 samples.Using the two age models discussed in the text-the conservative approach (top and bottom left) and the proximity approach (top and bottom right)-CO 2 and CH 4 are plotted against δD ice .Both of the age models show reduced range in the MPT and pre-MPT ice.We note, however, that the highest CO 2 and CH 4 values fall into the unassigned category in the strictest age model.Extended Data Fig. 5 | Evaluating the fraction of CO 2 range captured by 37 discrete samples.a, Ice core CO 2 record between 0 and 800 ka versus LR04 benthic foram δ 18 O stack synchronized onto the AICC2012 timescale 60,61 .The result of a linear regression is shown and used to calculate the synthetic CO 2 record between 1.2 and 2.0 Ma. b, A synthetic CO 2 record between 1.2 and 2.0 Ma reconstructed from LR04 δ 18 O and the regression parameters calculated in a.The range of CO 2 in this synthetic time series is 213-269 ppm.c-h, Fraction of the recovered CO 2 variability (observed range/true range) of the synthetic CO 2 time series (see Extended Data Fig. Data Fig.6| Evaluating the fraction of δ18 O atm and δXe/Kr range captured by discrete samples.a, b, Allan Hills δ18 O atm (a) and δXe/Kr (b) are potted against the δD ice of the same depth, colour coded according to their age units.Vertical dashed lines represent the range of δD ice observed in the 40k-world ice (−284 to −316‰).Notably, the range of δD ice associated with nine δXe/Kr values from the 40k world is only about 40% of the entire δD ice range in our 40k ice samples (−288 to −300‰); low δD ice values (less than −300‰) are missing from the δXe/Kr samples.Thus, we do not expect the nine Xe/Kr samples to fully capture the range of Xe/Kr ratios in the 40k-world ice.On the other hand, the co-depth δD ice values of δ18 O atm samples occupy 97% of the total range, implying that 29 δ 18 O atm samples are likely to cover 70% of the variability preserved in the 40k-world ice, barring any diffusive smoothing of the δ 18 O atm records.Extended DataFig.7 | Gas and ice properties in the interval between 123 and 143 m in ALHIC1503.a, δ 15 N; b, δO 2 /N 2 ; c, δ 18 O atm ; d, δD ice .The error bars associated with gas properties represent the pooled standard deviation (1σ) of all measurements: 0.010‰ for δ 15 N, 3.063‰ for δO 2 /N 2 , and 0.026‰ for δ 18 O atm .Note that the range of δ 18 O atm is less than 0.2‰.By comparison, glacial-tointerglacial δ 18 O atm variability in the 100-kyr climate cycles is about 1.5‰ 21 .Extended Data Fig. 8 | Evaluating the effect of diffusive mixing on ice and gas properties.a, Partitioning function Z (ratio of the number gas molecules in the gas phase to that in the ice phase) of O 2 (blue) and N 2 (red) as a function of temperature.b, Characteristic time scale of the self-diffusion of water molecules in ice (black), of O 2 permeation in ice (blue), and of N 2 permeation in ice (red), plotted as a function of temperature.Note that for water molecules the length scale of diffusion (L) is 0.1 m, whereas for N 2 and O 2 L = 0.5 m, reflective of their respective sampling resolutions.Extended Data Fig. 9 | Evaluating the possibility of flow-induced mixing in Allan Hills ice cores.a, A contour map shows the correlation coefficient (R 2 ) between δD ice and d in core ALHIC1503 as a function of the number of adjacent stable water isotope data points included in the calculations, from 5 to 85 with a step of 10.The sampling resolution is approximately 10 cm.The high correlation coefficient (>0.7) between 132 and 142 m raises the suspicion that mixing has taken place in this interval.However, the possibility of mixing is not supported by CO 2 and CH 4 .b-d, Cross-plots of δD ice -CO 2 (b), δD ice -CH 4 (c), and CO 2 -CH 4 (d) in intervals between 132 and 140 m in ALHIC1503.The black solid lines are mixing lines based on two end members at 135.30 m and at 139.76 m, where the minimum and maximum CO 2 values are measured, respectively.Most of the points do not fall onto the mixing line, implying that two-end-member mixing alone cannot explain the high correlation between δD ice and d.