Reconciling the apparent absence of a Last Glacial Maximum alpine glacial advance, Yukon Territory, Canada, through cosmogenic beryllium-10 and carbon-14 measurements

We present a new in situ produced cosmogenic beryllium-10 and carbon-14 nuclide chronology from two sets (outer and inner) of alpine glacier moraines from the Grey Hunter massif of southern Yukon Territory, Canada. The chronology of moraines deposited by alpine glaciers outside the limits of the Last Glacial Maximum (LGM) ice sheets potentially provides a less-ambiguous archive of mass balance, and hence climate, than can be inferred from the extents of ice sheets themselves. Results for both nuclides are inconclusive for the outer moraines, with evidence for pre-LGM deposition (beryllium-10) and Holocene deposition (carbon-14). Beryllium-10 results from the inner moraine are suggestive of canonical LGM deposition but with relatively high scatter. Conversely, in situ carbon-14 results from the inner moraines are tightly clustered and suggestive of terminal Younger Dryas deposition. We explore plausible scenarios leading to the observed differences between nuclides and find that the most parsimonious explanation for the outer moraines is that of pre-LGM deposition, but many of the sampled boulder surfaces were not exhumed from within the moraine until the Holocene. Our results thus imply that the inner and outer moraines sampled pre- and post-date the canonical LGM and that moraines dating to the LGM are lacking likely due to overriding by the subsequent Late Glacial/earliest Holocene advance.

1 Introduction

We present an in situ cosmogenic beryllium-10 and carbon-14 chronology of alpine glacier advances from interior Yukon Territory that escaped inundation by the Cordilleran Ice Sheet to infer climate conditions during the late Pleistocene. Relatively few records of Pleistocene climate exist from interior Yukon (e.g., Mahony, 2015). Alpine glaciers and ice sheets are commonly viewed as sensitive indicators of climate change, both past and present, as they integrate conditions affecting their mass balance and thus reflect regional climate (e.g., temperature and precipitation). Chronologies of ice sheet margins may differ from those of alpine glaciers because of differences in response time to climate perturbations (Jóhannesson et al., 1989), variable ice dynamics than can alter the position of ice divides (Fulton, 1991), and effects of debris cover on alpine glaciers (Scherler et al., 2011). Additionally, accumulation areas of ice sheets can be far removed from their margins, and thus part of its mass balance is set by climate different from that at the margin. Yet, the margins are where we can date deposits left behind by former ice sheets. Studies exploring such climate reconstructions in the Yukon largely have focused on ice sheet margins, notably the northern margins of the Laurentide (LIS) and Cordilleran (CIS) ice sheets (e.g., Margold et al., 2013a; Stroeven et al., 2010, 2014) and found somewhat ambiguous results largely a result of complex pre- and post-depositional processes affecting cosmogenic nuclide studies. Less studied are the current and former alpine glaciers in the region occupying the many massifs and ranges that escaped inundation by the LIS and CIS during glaciation, and that may provide a less ambiguous picture of climate forcing than records from ice sheet margins.

A common observation in cosmogenic nuclide-derived chronologies of alpine glacier and ice sheet fluctuations are a general lack of coherence from a single feature, or multiple features of the same morphostratigraphic age (e.g., Balco, 2011; Heyman et al., 2011). Both old and young outliers are observed in datasets, suggesting that pre- and post-depositional processes are operating (Applegate et al., 2012, 2010). Potential causes of old ages is inheritance of nuclides from earlier exposure periods and preservation through cover by cold-based ice. While young ages may result from post-depositional shielding of the sampled surface (e.g., exhumation of the surface, ephemeral cover). The measurement of multiple nuclides with differing half-lives, ideally one long lived and one short lived, can potentially not only elucidate the age of the moraines but also provide insight into the processes generating the observed scatter in moraine surface exposure ages.

Our objective is to develop a chronology for moraines that head from the Grey Hunter massif in Yukon Territory, which escaped inundation by the Cordilleran Ice Sheet. Such a site allows for the reconstruction of climate fluctuations during the late Pleistocene derived from simpler, and likely more responsive, alpine glacier systems relative to the adjacent Cordilleran Ice Sheet. This paper describes new measurements of in situ cosmogenic beryllium-10 (¹⁰Be) and carbon-14 (¹⁴C) of moraine boulders from the Grey Hunter massif, Yukon Territory, Canada, delineating former alpine glacier extents.

2 Study area

Grey Hunter Peak (∼ 2200 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$) and its adjacent terrain, herein referred to as Grey Hunter massif, is located within the generally low-lying MacArthur Mountains in east-central Yukon (Fig. 1). Grey Hunter massif is located adjacent to Tintina Trench and Tintina Fault, a large strike-slip fault that runs sub-parallel to the Denali Fault. The massif encompasses approximately 450 km² of terrain that preserves widespread evidence of Pleistocene glaciation as evidenced by the presence of cirques, aretes, and relatively wide U-shaped valleys. We selected the Grey Hunter massif because previous work suggests it escaped inundation by the Cordilleran Ice Sheet that flowed around the massif, but never overtopped it (Hughes, 1983). Today, the Grey Hunter massif hosts only small ice patches and rock glaciers; these features are limited to north-facing cirques. Grey Hunter massif predominantly consists of Cretaceous quartz monzonite of the Mayo Suite, which is the lithology of all collected samples. Mean annual temperature and precipitation at the nearby Stewart Crossing weather station (63.62^∘ N, 135.87^∘ W; 504 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$) are −4.3 ^∘C and approximately 300 mm yr⁻¹. The mean elevation difference between Stewart Crossing and mean elevation at Grey Hunter is approximately 1000 m. Assuming a simple moist adiabatic lapse rate of 6.5 ${}^{\circ }\mathrm{C}{\mathrm{km}}^{-\mathrm{1}}$, the mean annual temperature at Grey Hunter is approximately −10.8 ^∘C.

Figure 1 Map of the Grey Hunter region delineating the Grey Hunter massif (red box) as shown in Fig. 2 and former CIS glacial limits (shading). Also shown is the Glenlyon Range.

Figure 2 Shaded relief map of the Grey Hunter massif region. Outer moraines (orange) and inner moraines (purple) are shown. Sample locations are also shown for reference and labeled according to sample set ID. Samples comprising each same set are summarized in Table 1.

Table 1 Summary table of sample, corresponding sample set, morphostratigraphic age, exposure ages, and ${}^{\mathrm{14}}\mathrm{C}{/}^{\mathrm{10}}\mathrm{Be}$ ratios. All uncertainties are presented at the 1σ level. “Sat.” is the abbreviation for “saturated”.

3 Methods

Samples were selected based on their prominence relative to the surrounding moraine surface. We chose large boulders (> 1 m³) to minimize the chance of exhumation out of the moraine. Our sampling also targeted surfaces displaying signs of long-term exposure (e.g., lichen cover) over those with clean (e.g., unweathered and non-lichen-covered) surfaces that might be the result of recent exhumation or boulder spalling. We collected samples using a gas-powered rock saw to extract rock “brownies”, with approximate dimensions of 2 cm × 3 cm × 1 cm (L × W × H). A handheld GPS receiver was used to collect sample location and elevation, and we measured topographic shielding with a hand-held sighting clinometer by dividing the landscape into slope segments approximating the topography.

Sample preparation occurred in the Tulane University Cosmogenic Nuclide Laboratory. Samples were crushed, milled, and sieved to isolate the 250–500 µm size fraction. Quartz isolation followed standard laboratory procedures which includes rinsing with tap water, magnetic separation, froth flotation to remove feldspars, and a minimum of two 24 h 5 % HF/5 % HNO₃ etches on a shaker table and then two 24 h 1 % HF/1 % HNO₃ etches in an ultrasonic bath (Nichols and Goehring, 2019). This method consistently produces quartz with Al and Ti concentrations less than 150 ppm and ensures complete removal of dodecylamine, which may interfere with in situ ¹⁴C analysis. Preparation of quartz separates for ¹⁰Be and in situ ¹⁴C analysis is the same for each nuclide. An aliquot of material for samples where ¹⁴C was also measured was split from archived clean quartz material.

Chemical isolation of Be followed standard protocols, including sequential anion and cation columns (Ditchburn and Whitehead, 1994). All samples were spiked with ∼ 0.25 mg Be (nominal ${}^{\mathrm{10}}\mathrm{Be}/{}^{\mathrm{9}}\mathrm{Be}$ ratio ∼ 4 × 10⁻¹⁶). One process blank was prepared with each batch of eight samples. Be-isotope ratio measurements were made at the Purdue University Rate Isotope Measurement Laboratory (PRIME Lab), all normalized to dilutions as part of the standard KNSTD dilution series (Nishiizumi et al., 2007).

We extracted carbon from quartz separates using a fully automated gas extraction system optimized for the fusion of quartz using LiBO₂ and collection, purification, and graphitization of evolved carbon as gaseous CO₂ (Goehring et al., 2019). Process blanks are run every eighth sample. Additionally, the CRONUS-Earth CRONUS-A intercomparison material is run approximately every 40th sample to ensure tracking of long-term measurement stability and provide an estimate of sample reproducibility, which typically exceeds individual sample analytical precision. Carbon isotope ratios were measured at the National Ocean Science Accelerator Mass Spectrometry (NOSAMS) facility at Woods Hole Oceanographic Institution relative to the primary standard Ox-II (${}^{\mathrm{14}}\mathrm{C}/{}^{\mathrm{13}}\mathrm{C}$ = 1.4575 × 10⁻¹⁰, NIST SRM4990C). Primary and secondary (IAEA C7, ${}^{\mathrm{14}}\mathrm{C}/{}^{\mathrm{13}}\mathrm{C}$ = 5.45 × 10⁻¹¹) standards are graphitized using the same reactors used for unknown samples, ensuring fully internal standardization of results.

Exposure ages are calculated using the scaling method for neutrons and protons outlined in Lifton et al. (2014), and a simplified muon scheme presented in Balco (2017), all of which are modulated by the geomagnetic model of Lifton (2016). All ages assume zero surface erosion since deposition. Reference ¹⁰Be production rates were determined using the CRONUS-Earth “primary” calibration dataset (Borchers et al., 2016). For ¹⁴C production rates, we used the same production rate scaling method used for ¹⁰Be and calculated production rates by assuming that the CRONUS-A interlaboratory comparison standard (Jull et al., 2015; Goehring et al., 2019) displays production–decay saturation. All reported uncertainties are at the 1σ level and include full propagation of analytical errors (quartz mass, carbon mass, AMS uncertainty, blank uncertainty).

4 Results

Two groups of moraines (Fig. 2) down valley of rock glaciers fringe many of the cirques of the study area and were classified primarily by the relative weathering between the two moraine groups and the inner most rock glaciers. We refer to these two groups as the inner and outer moraines. Inner moraines display characteristics of stability and antiquity greater than the Little Ice Age (LIA). Most boulders are prominent relative to their surroundings (Fig. 3). No pronounced surface weathering of boulders is observed; most boulders have extensive, sometimes large lichen cover. Inner moraines typically extend 1.5 to 3 km down valley from cirque headwalls with a toe–headwall area ratio (THAR) equilibrium line altitude (ELA) estimate of 1565 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$

Figure 3 Example of boulders sampled from outer moraines (15-GH15; a) and inner moraines (15-GH11; b). Note in panel (a) the embedded character of the far side of the large erratic.

The outer moraines are decidedly more rounded in profile and boulders are generally more embedded, suggesting that boulder exhumation occurred. Boulder weathering characteristics are more evolved than the inner moraine. Moraines extend 5.5 to 7 km down valley from cirque headwalls with a THAR-derived ELA estimate of 1285 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$

Below we present resulting ¹⁰Be and ¹⁴C exposure ages, first from the inner moraines and then from the outer moraines (Tables S1–S3 in the Supplement). Individual moraine medians, as well as group medians (outer and inner moraines), are presented in Table 1 with uncertainties represented by the half width of the interquartile range. Where multiple ¹⁴C ages are finite (i.e., are not saturated with respect to ¹⁴C production–decay systematics), we present the median age as above.

4.1 Inner moraines

Resulting median ¹⁰Be exposure ages for the inner moraine sample sets are 18.2 ± 0.4, 15.1 ± 0.8, 16.3 ± 1.2, and 17.9 ± 1.5 ka (Fig. 4; Table S4); the group median age is 17.1 ± 1.0 ka. The quasi-standard error (ratio of half-interquartile range and median) ranges from 2.2 % to 8.5 %. Coherence for the younger moraines is better than that of the outer moraine sample sets. All measured samples from the younger moraine sample set (sample set 7) yield finite ¹⁴C exposure ages (Figs. 4 and 5) with a median ¹⁴C age of 11.2 ± 0.1 ka (n=3) and a quasi-standard error of 0.94 %.

Figure 4 Box-and-whisker plot of moraine ages for the inner moraines (purple) and the outer moraines (orange) as identified by their site number and nuclide. Median ages are indicated by black lines. Note the larger magnitude of scatter and outlier spread for the outer moraines relative to the inner moraines. The larger spread for the carbon measurements from set 5 results from the measurement of one saturated age.

Figure 5 Elevation vs. ¹⁴C concentration. Contours show the equivalent concentrations for a given exposure age as a function of elevation. The heavy black line represents saturation concentration for ¹⁴C. The measurements from the inner moraine (purple) are tightly clustered relative to those of the outer moraines (orange).

4.2 Outer moraines

The four outer moraines yield median ¹⁰Be ages of 24.7 ± 0.6, 38.0 ± 7.7, 14.1 ± 8.5, and 38.2 ± 2.1 ka (Tables S2 and S3). For a given outer moraine, we generally observe a large degree of scatter (Fig. 4); the one exception is sample set 1. The quasi-standard error for the outer moraine ranges from 2.5 % to 60.3 %. The grand group median for the outer moraines is 32 ± 8 ka.

In situ ¹⁴C results from sample set 5 yields one saturated (i.e., ≳ 30 kyr) exposure age, while two other samples yield finite ages of 9.4 ± 0.2 and 7.4 ± 0.2 ka. Similarly, all samples from sample set 8 yield finite ¹⁴C exposure ages of 7.6 ± 0.3, 5.4 ± 0.2, and 6.4 ± 0.2 ka. Results with respect to theoretical saturation concentrations for the elevation of a given sample is shown in Fig. 5 along with isochrons of exposure duration. Determination of saturation for a sample is based on a statistically significant difference at one standard deviation between the theoretical saturation concentration and the measured sample concentration.

4.3 ${}^{\mathrm{14}}\mathrm{C}{/}^{\mathrm{10}}\mathrm{Be}$ isotope ratios

A normalized paired-nuclide plot of ¹⁴C and ¹⁰Be concentrations can help identify complex exposure–burial histories (Fig. 6). Concentrations are normalized by sample specific total nuclide production rates (spallation plus muons) to yield concentrations equivalent to production rates of 1 $\mathrm{atom}{\mathrm{g}}^{-\mathrm{1}}{\mathrm{yr}}^{-\mathrm{1}}$ for each nuclide. All samples fall on or below the curve describing continuous exposure. The three inner moraine sample sets are consistent with continuous exposure and steady-state erosion or the possibility of minor amounts of ¹⁰Be inheritance (see below). Two samples from the two outer moraines display results consistent with continuous exposure, including one sample consistent with production–decay saturation (15-GH04). Most samples from the other moraines though are consistent with a history of complex exposure.

Figure 6 Paired ${}^{\mathrm{14}}\mathrm{C}{/}^{\mathrm{10}}\mathrm{Be}$ plot demonstrating that two outer moraine samples, GH05 and GH04, are consistent with a single-stage continuous exposure history with little or no surface erosion (upper bound of gray shaded region); the three inner moraines samples are consistent with continuous exposure but also could be affected by a multi-stage exposure history with ¹⁰Be inheritance; and finally four samples from the outer moraines are consistent with complex multi-stage exposure history. Note that ¹⁰Be inheritance serves to lower the ${}^{\mathrm{14}}\mathrm{C}{/}^{\mathrm{10}}\mathrm{Be}$ ratio and shift samples to the right. Outer (inner) moraine samples are shown in orange (purple).

5 Discussion

The large degree of scatter observed in the ¹⁰Be ages, as well as from any single moraine, argues for the presence of processes operating external to climatic controls on glacier extent and moraine deposition. Below, we will explore processes that could explain the observed ¹⁰Be and ¹⁴C concentrations/ages from a given moraine set. There are four possibilities affecting any given moraine dataset: (1) the apparent exposure ages are younger than the depositional ages as a result of erosion of the boulder surfaces during exposure of the sampled surfaces, (2) the apparent exposure ages are younger than the depositional ages as a result of significant exhumation from within a moraine after deposition, (3) ¹⁰Be exposure ages are older than the depositional age due to inheritance of nuclides from prior exposure, and (4) exposure ages, and therefore moraine ages represent a close minimum age for stabilization of the landform.

Weathering and erosion of boulders via the removal of mass, either physically or chemically, lowers the measured nuclide concentration from a boulder surface because, over time, these processes advect accumulated cosmogenic nuclides from a deeper depth in the rock to the surface. At Grey Hunter, where little gradient in climate is observed across the study area and lithologic variations of sampled boulders are slight, we assume negligible differences in surface erosion rates for the boulders sampled in our study. For the purposes of surface exposure dating, it is the rate of erosion that matters and is therefore presumed to be the same for both older and younger moraine samples. We stress that we do not imply uniform erosion rates between glacial and interglacial times. All samples are consistent with exposure subject to steady-state erosion (Fig. 6); due to the wide range in apparent ¹⁰Be exposure age for a given lithologically homogenous moraine and no evidence for systematic underrepresentation of exposure ages, however, we discount rock surface erosion as a significant post-depositional process affecting apparent exposure ages.

First-order observations indicate that the outer moraines have relative scatter greater than the inner moraines. However, for both the inner and outer moraines, the ¹⁰Be exposure age distributions are predominantly old biased and thus indicate the presence of inherited nuclides (Fig. 4). A common approach to interpreting such age distributions is to use the youngest sample from each moraine surface (e.g., Heyman et al., 2011). For the outer moraine sample sets (1, 3, 5, 8), this yields ¹⁰Be ages of approximately 24, 23, 9, and 32 ka. Similarly, for the inner moraine sample sets (2, 4, 6, 7), we observe youngest ages of approximately 18, 15, 14, and 16 ka. The outer moraines using the youngest age approach are generally older than the inner moraines but not in all instances. Additionally, this approach is incompatible with the morphostratigraphic requirement of the inner moraine post-dating the outer moraines within a single valley. Regardless, it is tempting to conclude that the outer moraines are a result of early retreat from the Last Glacial Maximum (LGM) and the inner moraines represent a stillstand during recession from the LGM (e.g., Margold et al., 2013a, b; Stroeven et al., 2010, 2014; Menounos et al., 2017). Applying the youngest age approach with the ¹⁰Be age is inconclusive at worst and unsatisfactory at best particularly given the observation that sample set 5 yields ¹⁰Be ages from 9 to 44 ka.

We further investigated the processes operating at Grey Hunter and more robustly determined moraine ages by measuring in situ ¹⁴C. In situ ¹⁴C should be less sensitive to inheritance given its short half-life. During stadial times, the combination of subglacial erosion and most importantly decay during ice cover reduces the number of nuclides present in the sample at the time of deposition. Further, given the differences in depth-dependent production between ¹⁰Be and ¹⁴C, we can investigate potential issues of exhumation of sampled boulders from moraine surfaces. This latter point is of potential importance since many of the sampled outer moraine boulders were less prominent relative to their surroundings than the sampled inner moraine boulders.

Before discussing paired ${}^{\mathrm{14}}\mathrm{C}{/}^{\mathrm{10}}\mathrm{Be}$ results from outer moraine samples, expected ¹⁴C concentrations and ${}^{\mathrm{14}}\mathrm{C}/{}^{\mathrm{10}}\mathrm{Be}$ concentration ratios is required. Our hypothesis with regards to moraine age implies that all samples should return finite ¹⁴C ages. Any resulting ¹⁴C ages that are saturated (Figs. 5 and 6) with respect to production decay systematics indicate that a boulder has resided at or near the surface for sufficient time (> 30 kyr) that production of ¹⁴C is balanced by decay of ¹⁴C. There is the potential for ¹⁴C inheritance, implying that a sample was delivered to a moraine with sufficient ¹⁴C to yield a saturated concentration in less than 30 kyr. We discuss this latter scenario first as we observe a single saturated sample from sample set 5.

We present two end-member model scenarios that could lead to inheritance in a moraine boulder sample, and we discuss model results in the context of our ¹⁴C measurements (Fig. 6). The first scenario involves the plucking of a pre-exposed surface from a valley bottom, and englacial incorporation of the resulting boulder and transport to its final moraine resting place (Fig. 7a). The second scenario envisions delivery of a pre-exposed paraglacial boulder to a supraglacial setting prior to burial and englacial transport to a moraine (i.e., Scherler and Egholm, 2020; Orr et al., 2021; Ward and Anderson, 2010; Fig. 7b). We use characteristic elevations for the valley bottom, cliff, and moraine settings and assume that in both scenarios, the ¹⁴C concentration is at production–decay secular equilibria and thus represents the maximum ¹⁴C concentrations that could be observed prior to delivery to the glacier surface. This means that as soon as the boulder is removed from its higher-elevation origin, there is excess ¹⁴C that begins to decay and is no longer balanced by production as the new lower elevation has an attendant lower production rate.

Figure 7 (a) Model of resulting ¹⁴C concentrations for a scenario of plucking of boulder from 1000 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$ that was saturated with respect to ¹⁴C prior to plucking. Contours show resulting ¹⁴C concentration (10⁵ atoms g⁻¹) for a range of plucking ages and moraine deposition ages. For reference, the saturation concentration for the mean older moraine elevation is 3.98 × 10^5 14C atoms g⁻¹. (b) As in panel (a), except that instead of the boulder being plucked by the advancing glacier, the boulder is delivered via rockfall from the headwall, assuming a saturated pre-transport concentration, and transported supraglacially down valley before deposition on the moraine.

Resulting theoretical ¹⁴C concentrations for the above scenarios indicate that for nearly all deposition ages, the resulting ¹⁴C concentrations are invariant with respect to the timing of plucking and entrainment or rockfall. Invariance implies that the observed ¹⁴C concentrations will be almost entirely set by the age of moraine deposition, and, due to the short half-life of ¹⁴C, only the youngest moraines (<10 ka) would have any memory of pre-exposure. The above scenarios though are unrealistic for Grey Hunter given plausible ages for the inner (e.g., ≥ Late Glacial) and outer moraines (e.g., ≥ LGM), and thus we dismiss the possibility of inheritance due to pre-exposure affecting our ¹⁴C measurements from either the inner or outer moraines, though not necessarily the paired ¹⁰Be measurements. In contrast, we conclude that the single saturated sample suggests relative antiquity (>30 ka) for sample set 5, which is further supported by the approximately 44 ka ¹⁰Be exposure age for sample 15-GH04 and a ${}^{\mathrm{14}}\mathrm{C}{/}^{\mathrm{10}}\mathrm{Be}$ ratio indicative of continuous surface exposure.

So far, we have ruled out inheritance affecting ¹⁴C measurements, thus allowing for some antiquity of the outer moraines. We also observe anomalously low apparent ¹⁴C ages from the outer moraines (Figs. 5 and 6), particularly relative to the inner moraines, along with ${}^{\mathrm{14}}\mathrm{C}{/}^{\mathrm{10}}\mathrm{Be}$ ratios indicative of burial (Fig. 6). This suggests that other processes are operating, notably exhumation of boulder surfaces from eroding moraines, as boulders from the outer moraines lack prominence relative to their surroundings and appeared embedded, in contrast to the more prominent boulders of the inner moraine samples.

The exhumation of boulders from a degrading moraine will lead to apparent exposure ages that are less than the true depositional age of the moraine since the boulders accumulated nuclides at a lower rate than at the surface while buried. The degree of scatter for a set of exhuming boulders is expected to be large because of variable nuclide inheritance at the time of deposition, but also because there is no a priori reason that the depth from which the boulder was exhumed needs be consistent from sample to sample. We can use paired ${}^{\mathrm{14}}\mathrm{C}{/}^{\mathrm{10}}\mathrm{Be}$ measurements to investigate the likelihood of exhumation processes occurring on the Grey Hunter moraines and leading to the observed scatter in apparent exposure ages and ${}^{\mathrm{14}}\mathrm{C}{/}^{\mathrm{10}}\mathrm{Be}$ ratios indicative of burial. However, because of the strong likelihood of ¹⁰Be inheritance resulting from exposure during multiple interglacial periods (Balco, 2011; Balco et al., 2014), we need to use ${}^{\mathrm{14}}\mathrm{C}{/}^{\mathrm{10}}\mathrm{Be}$ ratios with caution. A more robust approach is to exclusively look at ¹⁴C concentrations since our model above effectively rules out the presence of ¹⁴C inheritance. To explore this possibility, we model ¹⁴C concentration resulting from exhumation within an eroding moraine followed by exposure on the moraine surface. Important to our approach here is that we account for ¹⁴C production by muons, as the excitation function for production by negative muons (Heisinger et al.,  2002a, b) is large enough that significant subsurface production must be accounted for in all but the fastest of exhumation rates; to do so, we follow the simple approach outlined in Goehring et al. (2013). We assumed a density 2.4 g cm⁻³ for the eroding moraine and vary the duration of exhumation and the residence time on the surface, the sum of which represent the depositional age of the moraine and model resulting ¹⁴C concentrations for a range of exhumation depths (Fig. 8).

Figure 8 Contours of measured ¹⁴C concentrations compatible with exhumation and exposure duration pairs for outer moraine samples. The dashed black line represents theoretical and measured saturation concentration (15-GH04) and effectively sets the minimum age of the outer moraines given the unlikelihood of ¹⁴C inheritance. Orange shading represents the range of exhumation–exposure duration pairs compatible with saturated ¹⁴C concentration. The intersection of lower concentration contours with orange shading represents compatible exhumation depths and histories.

From our model, we can make some simple, yet powerful, inferences regarding the age of the outer moraines at Grey Hunter. For a boulder that has always resided at the surface of an outer moraine, the concentration is entirely set by duration of exposure and follows the expected trajectory for a sample accumulating nuclides at the surface (Fig. 6). For a given depositional age (exhumation duration plus surface exposure duration), the deeper the exhumation depth, the lower the effective or average nuclide production rate during exhumation largely due to the exponential decrease in spallation with depth. Therefore, for deeper samples a longer exposure duration is required at the surface (by consequence, exhumation rates are faster).

Contours outline exhumation–exposure duration pairs compatible with measured ¹⁴C concentrations (Fig. 8). Additionally, because ¹⁴C inheritance is unlikely (see above), we can use the single saturated sample (15-GH04) to set minimum depositional ages for an outer moraine assuming it has resided at the moraine surface for its entire history. Doing so yields a field (orange shading) of compatible exhumation–exposure duration pairs assuming all samples have a minimum depositional age of approximately 35 ka (Fig. 8). For the samples from the outer moraines with finite ¹⁴C ages, we note that one sample is potentially compatible with 200 cm of exhumation, while all other samples require at least 400 cm of exhumation over a minimum of 35 kyr. Additionally, the saturated sample was deposited at least 35 ka and older if exhumed from any depth. We can also infer that the young apparent ¹⁴C age samples must have ¹⁰Be inheritance from prior periods of exposure to yield the observed ${}^{\mathrm{14}}\mathrm{C}{/}^{\mathrm{10}}\mathrm{Be}$ ratio, and that all samples reached the moraine surface during the Holocene. These magnitudes of exhumation (erosion) are not unprecedented, particularly during Holocene warming leading to the melting of ice lenses and solifluction events as observed elsewhere in the Arctic during the Pleistocene–Holocene transition (e.g., Mann et al., 2002). An outstanding question that arises, however, is why the inner moraines that pre-date the Holocene do not require similar magnitudes of exhumation. We speculate in this case that given the higher mean elevation of the inner moraines means that they were effectively insulated from early to mid-Holocene warmth and never experienced sufficient warmth for the melting of ice lenses and subsequent solifluction. Our inference here is further supported by the presence of rock glaciers proximal to the inner moraines.

To summarize, the scatter observed in both outer and inner moraine sets complicate interpretation and assessment of moraine ages. For the outer moraines, we are limited in our allowable inferences. Deposition pre-dates the LGM in its most classical western Cordillera sense and is 32 ± 8 ka based on the ¹⁰Be ages, regardless of the depth of exhumation, and is at least 35 ka based on ¹⁴C. Meanwhile, the tight clustering of ¹⁴C ages (11.2 ± 0.1 ka) from a single inner moraine, combined with the preponderance of young ¹⁰Be ages, suggests that the inner moraines post-date the LGM and are correlative with major retreat events in the Canadian Cordillera at the same time (Menounos et al., 2017; Darvill et al.,  2018; Lesnek et al., 2018). We are thus left with the conundrum of the lack of apparent retreat from a set of moraines at the end of LGM in the canonical sense (i.e., 26–19 ka) but later retreat concordant with widespread retreat at the Pleistocene–Holocene boundary.

6 Climatic implications

Our results suggest that the canonical LGM is absent or at least not preserved in alpine glacier deposits at Grey Hunter. We cannot, however, rule out that there was an LGM advance smaller in magnitude and therefore overridden by the advance associated with the inner moraines at the Pleistocene–Holocene transition. In contrast, the northern margin of the CIS appears to have classical LGM maxima ages and therefore there appears to be a different glacier response to past climate for Grey Hunter alpine glaciers (e.g., Margold et al., 2013a, b; Stroeven et al., 2010; Stroeven et al., 2014; Menounos et al., 2017). It is important to caveat that we cannot rule out that the outer moraines were reoccupied during the canonical LGM and thus contributing to the complex ¹⁰Be exposure age distributions but not supported by our data as well as that of a Marine Isotope Stage (MIS)$\mathrm{3}/\mathrm{4}$ advance. A possible explanation for the observed chronology is that of moisture starvation during the LGM at Grey Hunter. The presence of the LIS to the east-southeast would create anticyclonic circulation over the general Yukon region, drawing down cold–dry air as katabatic winds off the LIS and from the north largely devoid of precipitation (Manabe and Broccoli, 1985; Tulenko et al., 2020; Löfverström and Liakka, 2016). Thus, as climate deteriorated into the LGM, the expansion of the CIS and LIS hindered advance of the Grey Hunter glaciers beyond valley mouths. The expanded extent during MIS$\mathrm{3}/\mathrm{4}$ in contrast is supported by elevated elemental carbon abundance during this time in Yukon permafrost, suggestive of greater productivity and generally wetter conditions in the Yukon (Mahony, 2015). The same record also shows a marked increase in carbon content following the LGM and suggests that while temperatures likely remained colder than today during the latest Pleistocene, precipitation relative to the LGM increased because of collapse of the CIS and LIS (Fig. 9; Carlson et al., 2012; Liu et al., 2009), allowing for the expansion of glaciers at Grey Hunter.

Figure 9 (a) Median ages of the inner and outer moraines for ¹⁰Be and ¹⁴C. (b, c) TRACE 21k surface temperature difference and change in precipitation relative to present for the Grey Hunter region (Liu et al., 2009). (d) Carbon and nitrogen content for the permafrost from the Yukon Territory (Mahony, 2015).

7 Conclusions

Barring processes occurring post-deposition, most studies identified inheritance as the largest contributor to scatter in surface exposure age populations from moraines (e.g., Heyman et al., 2011; Applegate et al., 2012, 2010). Our results from Grey Hunter are in accord with these studies, at least in terms of the resulting ¹⁰Be ages. The resulting ¹⁴C concentrations (ages), however, cannot be explained by simple models of inheritance and instead require significant magnitudes of exhumation of boulders from within the moraine. Based on our exhumation model, the most parsimonious interpretation of outer moraine ages is that they pre-date the LGM. The inner moraines on the other hand have clear ¹⁰Be inheritance and young coherent ¹⁴C ages and date to the termination of the Pleistocene. Our observations at Grey Hunter discord with observations of coherent surface exposure ages from other alpine glacier locations in British Columbia and the Yukon that were formerly covered by the CIS (Darvill et al., 2018; Menounos et al., 2017). A notable exception is the Glenlyon Range that largely avoided LGM CIS inundation and similarly has evidence for advance of alpine glaciers during MIS3 or earlier and at the Pleistocene–Holocene transition (Stroeven et al., 2014; Ward and Jackson, 1992).

Our results are not only informative about the climate of the Yukon outside of the CIS limits, but also bear directly on the complications inherent in such a study. Our records from Grey Hunter massif suggest that conditions favorable for ice growth occurred during MIS$\mathrm{3}/\mathrm{4}$ and at the Pleistocene–Holocene transition when temperatures were still reduced relative to the Holocene, but when the configurations of the CIS and LIS were such that enough precipitation reached Grey Hunter to induce glacier expansion. The cold–dry climate operating at Grey Hunter also means that glaciers were likely cold based and resulted in the complicated exposure age distribution observed and that the use of ¹⁰Be and ¹⁴C yields information on processes leading to the observed exposure age distribution.