Disappearance of an ice sheet The Cordilleran Ice Sheet is thought to have covered westernmost Canada until about 13,000 years ago, even though the warming and sea level rise of the last deglaciation had begun more than a thousand years earlier. This out-of-phase behavior has puzzled glaciologists because it is not clear what mechanisms could account for it. Menounos et al. report measurements of the ages of cirque and valley glaciers that show that much of western Canada was ice-free as early as 14,000 years ago—a finding that better agrees with the record of global ice volume (see the Perspective by Marcott and Shakun). Previous reconstructions seem not to have adequately reflected the complexity of ice sheet decay. Science, this issue p. 781; see also p. 721

Abstract The Cordilleran Ice Sheet (CIS) once covered an area comparable to that of Greenland. Previous geologic evidence and numerical models indicate that the ice sheet covered much of westernmost Canada as late as 12.5 thousand years ago (ka). New data indicate that substantial areas throughout westernmost Canada were ice free prior to 12.5 ka and some as early as 14.0 ka, with implications for climate dynamics and the timing of meltwater discharge to the Pacific and Arctic oceans. Early Bølling-Allerød warmth halved the mass of the CIS in as little as 500 years, causing 2.5 to 3.0 meters of sea-level rise. Dozens of cirque and valley glaciers, along with the southern margin of the CIS, advanced into recently deglaciated regions during the Bølling-Allerød and Younger Dryas.

Beginning about 14.7 thousand years ago (ka), near the end of the Pleistocene, Earth’s climate rapidly oscillated between glacial and interglacial states (1). The first substantial warming phase triggered considerable mass loss from glaciers and ice sheets (2) and was followed by several short-lived cooling events during the Bølling-Allerød (BA) (14.6 to 12.9 ka). These short-term intervals of cooling gave way to glacial conditions during the Younger Dryas (YD) (12.9 to 11.7 ka), with an abrupt return to interglacial conditions at 11.7 ka. Climate records from Greenland (3) and the Gulf of Alaska (4) clearly record these climate oscillations, but the degree to which they affected a large proportion of the North American continent remains uncertain. This uncertainty stems in part from the presence of continental ice sheets, which limit terrestrial climate proxy records to their peripheries. The Cordilleran Ice Sheet (CIS) was one of these ice masses that, at its maximum extent, contained a sea-level equivalent comparable to that of the present-day Greenland Ice Sheet (5). Like the Greenland Ice Sheet, the CIS was underlain by rugged mountains and had marine-terminating margins. The CIS response to abrupt climate change during the latest Pleistocene provides an analog for the behavior of the Greenland Ice Sheet in the future.

Previous mapping and dating of glacial deposits from low elevations (6) and ice sheet modeling (5, 7–9) suggest that the CIS covered large portions of westernmost Canada as late as 12.5 ka. The southern margin of the CIS and alpine glaciers advanced during the YD (10–13), a time for which paleoecological records (14, 15) indicate a return to cold conditions along British Columbia’s west coast. But whether alpine glaciers and other margins of the CIS synchronously responded to abrupt climate change during the latest Pleistocene is unknown. A refined CIS chronology bears directly on changes in atmospheric circulation (16), human migration into North America (17, 18), and the timing and magnitude of meltwater transfers to the Arctic and Pacific oceans (19). Recent modeling studies (5, 7, 19, 20) differ in their simulations of the response of the CIS to climate oscillations, illustrating the need for robust geologic evidence, which can in turn be used to validate and improve numerical models that simulate ice sheet behavior.

Numerical ages from erratic boulders and bedrock surfaces are commonly used to constrain ice sheet retreat. In areas of high relief, numerical ages from cirque and valley moraines offer complementary data to constrain the time of ice sheet decay because these features can only form after regional ice sheet deglaciation. To evaluate whether the CIS responded to latest Pleistocene climate reversals, we mapped and dated moraines throughout the Canadian Cordillera [Fig. 1 and supplementary materials (SM)]. Latest Pleistocene moraines occur as both small “cirque moraines” that lie 0.5 to 1.0 km down-valley from contemporary glaciers and “valley moraines” formed by glaciers sourced from cirques but that were many times more extensive than the former (see SM).

Fig. 1 The Cordilleran Ice Sheet and sites of latest Pleistocene moraines. (A) Continental ice in North America at 16 ka (6). Yellow circles show locations of climate proxy records discussed in the text and shown in Fig. 4. (B) The Cordilleran Ice Sheet (6) at 16, 14, and 12.5 ka. Numbers refer to 10Be sites of this study (table S1), and letters denote sites of previously described glacier advances in the area (table S2). Colors denote median ages at these sites pertinent to decay of the Cordilleran Ice Sheet.

Our moraine and erratic boulder chronology includes 76 10Be surface exposure ages from 26 locations (Fig. 1). We collected samples from bedrock and stable boulders associated with moraines and from scattered erratic boulders (table S1). At some sites, we also obtained radiocarbon ages from lakes impounded by the moraines or from till, both of which provide minimum limiting ages of the moraines (Fig. 1 and table S2). The 10Be data set includes samples from eight valley and eight cirque moraines (Fig. 2). Unlike previously reported numerical ages that are generally from low elevation sites that become younger toward the center of the CIS, our 10Be ages from alpine sites reveal a considerably more complex picture of CIS deglaciation (Fig. 1).

Fig. 2 Distribution of 10Be ages. Box-and-whisker plots of moraines ages. (Top) Age group. (Bottom) Locality. Sites are arranged from south to north (refer to Fig. 1 and table S1). Black lines and boxes denote, respectively, median and IQR, and whiskers show the 5 to 95 percentile spread. Black diamonds denote outliers (outside the 5 to 95 percentile range).

Four 10Be ages from the moraine data set exceed 22 ka, likely due to inheritance from previous exposure (table S1). We exclude these ages and use the median and interquartile range (IQR) of the other 48 sample ages to define their central tendency and dispersion (Fig. 2). The IQR provides a conservative estimate of the uncertainty of the moraine’s true age in light of environmental factors that could bias the age of the deposits, such as inheritance or incomplete exposure. The moraine ages fall into two groups, with a combined median age (±1 quartile) of the older group of 13.9 ± 1.1 ka, and a combined median age of the younger group of 11.4 ± 1.6 ka. Within each of the age groups, no significant age difference exists between the two moraine types (Fig. 2). Geomorphic factors such as moraine stabilization, snow cover, or erosion may account for apparent earliest Holocene ages for some of the samples (see SM).

The moraine ages suggest that after substantial CIS decay, alpine glaciers in western Canada advanced during the BA and YD. Our results and ages from previous work (Fig. 1) indicate that several mountainous regions in western Canada emerged from the CIS during the BA (Fig. 1 and SM). Two recent glacio-isostatic adjustment (GIA) models supported by data calibration from records of sea level, paleo lake shorelines, and present-day geodetic measurements are consistent with this view because they reveal significant thinning of the CIS between 14.5 and 14.0 ka (Fig. 3). Some discrepancies nonetheless exist between the GIA models and our deglaciation ages (Fig. 1). The Laurentide and Cordilleran ice sheets, for example, remain joined until 13.5 ka in one model (21), whereas separation is complete by 14.5 ka in the other (8). However, both GIA models show contiguous ice at 12.5 ka (9) in areas in which our data indicate that early BA deglaciation and correlative moraines postdate regional CIS decay (Fig. 1 and SM).

Fig. 3 Rapid ice sheet thinning and climate anomalies within the Bølling-Allerød and Younger Dryas. (A) GIA-modeled (8, 21) surface ice elevation change for the Cordilleran Ice Sheet and western sector of the Laurentide Ice Sheet between 14.5 and 14.0 ka. (B) (Left) TraCE-simulated (22) anomalies in surface air temperature between 14.0 to 13.8 ka and 14.4 to 14.2 ka. Dark gray shading is the extent of ice cover at 14 ka. (Right) Thicker and dashed contours denote, respectively, zero and negative anomalies in surface air temperature from 12.2 to 12.0 ka and 12.6 to 12.5 ka. Dark gray shading (8) is the extent of ice cover at 12 ka. (C) Same as (B) but for precipitation anomalies.

The age equivalence of moraines of greatly differing down-valley extent that formed after regional CIS decay complicates their use for reconstructing the latest Pleistocene climate in the Canadian Cordillera (see SM). We instead use TraCE-21ka, a transient climate simulation experiment (22), to estimate the magnitude and pattern of temperature and precipitation changes that may have caused alpine glaciers and the southwestern margin of the reduced CIS to advance during the BA and YD (Fig. 3). TraCE-21ka includes time-varying changes in greenhouse gases, insolation, and the extent of continental ice sheets (23). The simulation also includes freshwater inputs to the North Atlantic that slow thermohaline circulation and cause two notable intervals of cooling in western North America (46° to 65°N, 158° to 98°W), one within the BA and the other within the YD (Fig. 3). The modeled amplitudes of BA and YD cooling for sites with the latest Pleistocene moraines exceed, respectively, 1.5° and 0.5°C. Simulated warming in central Canada during these events can be explained by ice-albedo feedback as the Laurentide Ice Sheet retreated eastward. The simulation also reveals wetter conditions for the mid- to southern latitudes of the CIS for both climate anomalies (Fig. 3).

Abrupt warming at the onset of the BA caused significant loss of mass of the CIS and provided high-elevation space for the construction of younger, latest Pleistocene moraines reported in this study. Gulf of Alaska sea surface temperature (SST) reconstructions (24) show that early BA temperatures were nearly as warm as those of the earliest Holocene when Northern Hemisphere summer insolation peaked (Fig. 4). SST likewise rose by 4°C off the British Columbian coast (15) between 15.5 and 14 ka. Temperature increases of this magnitude provided the mechanism to rapidly melt the CIS and place its equilibrium line potentially above all but the highest peaks in the region (19, 20). Substantial mass loss during this period may also explain depletion of δ18O in planktic foraminifera (25) from the Gulf of Alaska via ocean surface freshening (Fig. 4). GIA results (8, 21) show that between 14.5 and 14.0 ka, the CIS alone contributed 2.5 to 3.0 m of sea-level rise, representing about 17.5 to 21% of meltwater pulse 1A (Fig. 4). Those models also reveal that the CIS lost nearly one-half of its Last Glacial Maximum mass in as little as 500 years (Fig. 4) and changed from a continental-scale ice sheet to a complex network of alpine glaciers, ice fields, and ice caps. Our data also confirm that the CIS had disappeared by the end of the YD (Fig. 1), likely in response to abrupt warmth that signaled the end of the Pleistocene (1).

Fig. 4 Latest Pleistocene climate proxies and moraine record from western Canada. (A) North Greenland Ice Core Project (δ18O) (28) and mean 60°N summer insolation (29). (B) Sea-level equivalent of ice mass lost from the CIS between 18 and 8 ka based on GIA models Ice-6G (8) and LW-6 (21). (C) CCSM3 TraCE surface air-temperature anomalies relative to 1960 to 1990 CE averaged over 46° to 65°N, 158° to 98°W (22). (D) Alkenone-based SST (24) and planktonic foraminifera δ18O (25) from Gulf of Alaska core EW0408-85JC. (E) Median and interquartile age ranges (horizontal bars) for (i) 10Be exposure ages (n = 48) from cirque (blue) and valley (green) moraines and (ii) two valley moraines from an earlier study (30) (thin green horizontal bars). Shaded curves are normal kernel density functions of cirque (blue) and valley (green) ages. YD, Younger Dryas; BA, Bølling-Allerød; MWP1A, meltwater pulse 1A. Also shown are calibrated radiocarbon age ranges (±2 SD) for glacier advances in western Canada (brown rectangles) (table S1).

Ironically, strong loss of mass from the CIS and the adjacent Laurentide Ice Sheet during the early BA may have helped trigger climate deterioration that slowed its demise and allowed alpine and cirque glaciers to advance (Fig. 4). As shown in the TraCE-21ka experiment (22) and in more recent experiments of a combined general circulation model and hydrologic model that simulates Laurentide and Cordilleran Ice Sheet saddle collapse at ~14.5 ka (20), introduction of meltwater into the Arctic and Atlantic oceans can weaken the Atlantic meridional overturning circulation (AMOC). This weakening, in turn, can lead to widespread hemispheric cooling that could cause alpine glaciers and margins of the CIS to advance into newly deglaciated terrain. Meltwater-induced changes in AMOC are, likewise, hypothesized to have initiated the YD (26). The greatly reduced volume of the CIS before the inception of the YD (Fig. 4) argues against CIS meltwater playing any substantial role as a trigger for the YD.

Our study reveals that the CIS and associated alpine glaciers responded to abrupt climate change (Fig. 4) and challenges the traditional view that the CIS covered large regions of westernmost Canada as late as 12.5 ka. Our results also support GIA models that place major mass loss of the CIS in line with other Northern Hemisphere ice sheets at ~14.5 ka, in response to abrupt, hemisphere-wide climate amelioration (3, 4). Alpine areas in the region had clearly emerged from the ice sheet before the BA, and alpine glaciers advanced soon after the onset of the BA and during the YD (Fig. 4).

The current understanding of the peopling of the Americas requires humans to have migrated to the south of the ice sheets after the end of the Last Glacial Maximum (~18 ka) but before 14.6 ka (18). Existing geologic evidence (17) and ice sheet models (9) rule out migration between the Cordilleran and Laurentide ice sheets before 13.4 ka. Although some intermediate and high elevation sites became ice free during the BA, lower elevations in the interior of British Columbia did not become ice free until the end of the YD (Fig. 1), making the route for human migration across the Cordillera unlikely during the latest Pleistocene.

The complexity of CIS decay can explain the age equivalence of moraines (Fig. 2) constructed by glaciers of different lengths (see SM). The thinning CIS exposed rugged mountainous areas that were transformed into a labyrinth of valley glaciers. In such a scenario, cirque glaciers might reform and construct small moraines during climate reversals such as the BA and YD. Simultaneously, valley glaciers left from the decaying ice sheet might reinvigorate and advance to positions many kilometers down-valley from cirque headwalls during climate reversals. Although marginal retreat was common during the demise of the CIS (6), especially over lower relief topography, our hypothesis for substantial ice loss at high elevations accords with both conceptual (27) and GIA-based models (8, 21). Numerical ice flow models (5, 7, 19) support widespread mass loss through marginal retreat and thinning of the CIS during the period 14 to 10 ka, and moderate resolution (5 km) models (5) simulate ice expansion during the YD. However, these simulations fail to capture the complexity of CIS decay, especially once the ice surface elevation approaches the height of the region’s many mountain ranges. Simulating the complex pattern of ice sheet decay in montane environments such as in western Canada or Greenland will require novel approaches that couple sophisticated surface mass balance and ice dynamics models at high spatial scales.

Supplementary Materials www.sciencemag.org/content/358/6364/781/suppl/DC1 Materials and Methods Supplementary Text Figs. S1 and S2 Data S1 to S4 References (31–76)