Research ArticleGEOLOGY

Ice-stream demise dynamically conditioned by trough shape and bed strength

See allHide authors and affiliations

Science Advances  24 Apr 2019:
Vol. 5, no. 4, eaau1380
DOI: 10.1126/sciadv.aau1380


Ice sheet mass loss is currently dominated by fast-flowing glaciers (ice streams) terminating in the ocean as ice shelves and resting on beds below sea level. The factors controlling ice-stream flow and retreat over longer time scales (>100 years), especially the role of three-dimensional bed shape and bed strength, remain major uncertainties. We focus on a former ice stream where trough shape and bed substrate are known, or can be defined, to reconstruct ice-stream retreat history and grounding-line movements over 15 millennia since the Last Glacial Maximum. We identify a major behavioral step change around 18,500 to 16,000 years ago—out of tune with external forcing factors—associated with the collapse of floating ice sectors and rapid ice-front retreat. We attribute this step change to a marked geological transition from a soft/weak bed to a hard/strong bed coincident with a change in trough geometry. Both these factors conditioned and ultimately hastened ice-stream demise.


Ice streams, fast-flowing mass conveyors within ice sheets, are important components of Earth’s interconnected ice sheet–ocean–climate system that dominate current and near-future cryospheric losses (1,2). While climatically driven mass-balance changes determine the volume of ice sheets, there is a growing body of empirical and experimental evidence showing that bed topography and bed properties exert a fundamental control on ice sheet flow dynamics and ice flux at the grounding line—the point at which ice begins to float—in ocean-terminating glaciers (35). Considerable attention has focused on ice streams with bed profiles that deepen inshore in an attempt to better understand ice sheet vulnerability (the theory of “marine ice-sheet instability”) and to better quantify future sea level rise (69). However, there is still much uncertainty regarding the determination of bed elevations beneath large parts of the Antarctic Ice Sheet (1012) and very limited information regarding subglacial bed geology or substrate properties, with access to the beds of contemporary ice streams exceedingly rare. Recent work has highlighted the importance of bed strength and force resistance at the ice-bed interface on the dynamical behavior of ice streams and large tidewater glaciers (5, 11, 13). Modeling experiments show that geometry-controlled retreat behavior is driven by the relationship between ice flux and water depth at the grounding line, with greater flux occurring in areas of deeper water (4, 6, 9). These models show that an ice stream in equilibrium will retreat or readvance periodically, adjusting to the trough width and depth upstream of the terminus (811). However, models typically use only low-resolution topographic boundary conditions (>1-km grid), and most do not include information on subglacial geology. Consequently, geometry-controlled retreat behavior and basal traction are only generally constrained or weakly parameterized in current numerical ice sheet models (12, 14), and hence, key controls on ice-stream retreat and ice-shelf stability over longer (102 to 103 years) time scales remain poorly understood.

The dynamic evolution of the Eurasian ice sheet at the end of the last Ice Age (Pleistocene era) provides an analog for the current and future behavior of the West Antarctic and Greenland ice sheets. We target a sector of the former British-Irish Ice Sheet, in northwest (NW) Europe, that encompasses both normal and reverse-sloping bed sections, as well as differing geological bed conditions, to constrain the complete record of ice-stream decay—its retreat style, geometry, and timing—since the Last Glacial Maximum (LGM) 19 to 29 thousand years (ka) before present (BP) (Fig. 1) (15). During the last glaciation, the Minch Ice Stream (MnIS) drained the NW sector of the British Ice Sheet and flowed within a 30- to 40-km-wide cross-shelf trough, terminating near the continental shelf break and feeding glacigenic sediment into a large 3500-km2 fan on the continental slope (16, 17). In places, the cross-shelf trough displays well-preserved streamlined subglacial landforms, including megascale glacial lineations (17)—the hallmarks of fast unidirectional ice flow in the past. The MnIS system is estimated to have been similar in size and discharge flux to the present-day Rutford Ice Stream in West Antarctica, with a discharge flux of 12 to 20 Gt a−1 (17, 18).

Fig. 1 Bathymetric-topographic map of palaeo ice-stream catchment and landforms offshore NW Scotland.

(A) Location of study area in NW Europe. (B) Map of grounding-zone wedges (GZWs) 1 to 17 recording grounding-line positions (thick black lines) during ice-stream/ice-shelf retreat; moraines record ice sheet margin positions (thin black lines) in nonstreaming sectors. New 10Be exposure ages (black squares) and optically stimulated luminescence (OSL) ages (black circles) in ka BP (means of several samples, without uncertainties). Previously published 10Be ages (white squares) also shown (see the Supplementary Materials). Simplified substrate (ice-stream bed) properties derived from terrestrial and marine geological mapping (see fig. S1). Hard bedrock at or near surface (dark gray); weak/soft sediments, >10 m thick (light gray); deglacial sediments (light green). New cores sites also shown (red circles). Basemap is merged grayscale elevation grid at a cell size of 250 m, hill-shaded from the NW. Bathymetric contours at water depths of 50, 100, 130, and 200 m on the shelf and at 100-m vertical intervals on the continental slope. Continental shelf break is defined by the 200-m water-depth contour.


Geology and geomorphology

Using new and existing high-resolution ship-borne bathymetry data and subseabed acoustic profiler data, we have mapped the distribution and morphology of submarine glacial landforms and characterized the seabed geology within the former ice-stream track, covering approximately 12,000 km2 of seafloor (Fig. 1). Thirty-eight sediment cores, 1 to 10 m in depth, were also taken at key locations to prove the geological composition of certain submarine landforms and basins. Existing shelf-wide datasets of seabed geology and Quaternary sediment thickness (19, 20), and newly acquired acoustic sub-bottom profiler data, were used to map and classify ice-stream bed strength (Figs. 1 and 2, figs. S1 and S2, and table S1). Rather than descriptive classifications of the subglacial geology, we express bed composition in terms of compressive strength measurements. Alongside new geochronological data, these datasets afford this former ice-stream system the potential to be one of the best characterized in the world, especially when compared to contemporary ice streams whose beds are largely inaccessible and incompletely surveyed.

Fig. 2 GZW distribution, ice-stream trough geometry, and bed geology.

(Top) Trough elevation versus distance diagram showing an along-flowline elevation profile of ice-stream bed and mapped GZW locations. Note that GZWs 1 to 9 occur on low-strength deformable-sediment bed; GZWs 10 to 17 occur on high-strength bedrock-dominated bed. (Bottom) Trough width versus distance diagram showing along-flowline variations in trough width (defined by the 100-m water-depth contour) and positions of GZWs. Note bifurcation of trough into east and west troughs at a distance of ~150 km from the shelf break. Blue bars show regions of seafloor where grounding-line features are notably absent, indicating uninterrupted ice-front retreat. GZW locations (GZW 1 to 17) shown along the distance axes. The T marks final grounding at the terrestrial transition.

We identify a sequence of broad, wedge-shaped deposits with markedly asymmetrical profiles (steeper down-ice slopes) crossing the former ice-stream trough (Figs. 1 and 2). On the basis of their acoustically structureless character and morphological affinities with similar features on other glaciated continental margins, these landforms are interpreted as grounding-zone wedges (GZWs) (fig. S3) (21). GZWs form by subglacial sediment deposition at the junction between grounded and floating ice typically within the reduced vertical accommodation space beneath ice shelves (21). We map 17 GZWs of differing sizes, but similar characteristic asymmetric morphologies, along the bed of the MnIS (Fig. 1). The mapped GZW distribution charts the migration of the ice-stream grounding zone over time and indicates the presence of an ice shelf (or floating ice tongue) at the time of deposition. The GZWs generally fall into two size categories (figs. S3 and S4): GZWs in the outer part of the trough are longer in the direction of ice flow, wider, and considerably greater in total sediment volume than GZWs further inshore. The mean volume of GZWs 1 to 9 is an order of magnitude larger than the mean volume of GZWs 10 to 17. There is no linear decrease in GZW size or volume with distance inshore, as might be expected as the ice stream recedes and shrinks in size; instead, we find a marked (>100%) size reduction at a distance of ~160 km along the trough centerline (Figs. 1 and 2 and fig. S4). This landform change spatially coincides with a major geological substrate transition at the ice-stream bed—from weak soft sediment to strong hard bedrock (around 58°10′N)—and a change in trough geometry from a single wide trough to two smaller feeder troughs (around 58°15′N) (Figs. 1 and 2). Typical compressive strengths of substrates range over three orders of magnitude across this geological transition: from 0.1 to 1 MPa for soft clay-rich sediments in the outer trough to 100 to 300 MPa for hard to very hard crystalline bedrock further inshore (fig. S2 and table S1). Trough widths vary considerably across this transition, from 30 to 35 km in the main trough to 10 to 15 km in each of the separate feeder troughs, with a high degree of trough-width variability moving inshore. This inner part of the ice-stream trough is also associated with a marked increase in bed relief, with numerous isolated bedrock highs and seabed deeps close to the Scottish mainland (Figs. 1 and 2).

The distribution of GZWs in the MnIS trough reflects progressive upstream migration of the grounding line during overall ice-stream retreat, with the three-dimensional (3D) size, spacing, and location of GZWs providing empirical information on past ice-stream/ice-shelf behavior. We find that GZW spacing is quasi-regular, decreasing with distance in the soft-bedded part of the trough with a high degree of predictability (GZWs 1 to 9) (fig. S4). Seven of the nine GZWs on a soft bed occur at trough narrowings (Fig. 2 and fig. S4). This colocation indicates that grounding-zone stabilizations on a soft/weak bed were probably governed by subtle but significant changes in lateral-drag forces and ice-discharge flux associated with variations in trough width (9). Over this entire distance (~100 km), water depths along the trough axis remain generally constant (±25 m), showing that changes in water depth–driven ice-discharge flux were not a major driver of grounding-line stabilization in this part of the ice-stream trough. We propose that the quasi-regular spacing of larger GZWs (19) on a low-relief soft bed highlights episodic grounding-line stability best explained by ice-shelf buttressing, peaks in lateral drag, and slowdowns at trough narrowings during overall ice-stream retreat.

The size, plan shape, and spacing of GZWs change markedly in the hard-bedded section of the ice-stream trough coinciding with a major change in trough geometry. GZWs 10 to 17 are irregular in planform shape, more closely spaced, and significantly smaller in volume than those on a soft bed, with a reversal of the predictable spacing trend seen in GZWs 1 to 9 (fig. S4). South (i.e., upstream) of GZW9, the west trough is conspicuously free of grounding-line landforms (buried or at seabed), indicating a period of uninterrupted ice-stream retreat or frontal collapse over a flowline distance of 40 km (~500-km2 area). Acoustic sub-bottom profiler data (Fig. 3) also show that GZW9 is buried in places beneath 10 to 20 m of strongly layered conformably draped sediments, in contrast to GZWs 6 to 8, which are exposed at seabed. These fine-grained draped muds and silts, recovered in seabed cores, are laminated on a millimeter-to-centimeter scale with numerous dropstones and interbedded with gravelly diamicton units on a centimeter-to-decimeter scale (fig. S9), indicating ice-proximal open-water glaciomarine conditions with abundant debris-charged icebergs. These sediment facies are characteristic of a high-sedimentation “calving bay” environment (22, 23), typically found when the floating sections of marine-terminating glaciers calve rapidly as the grounding line stabilizes in shallower water. Large calving bays are also seen adjacent to modern-day ice-shelf remnants, such as Larsen A and B in Antarctica, immediately after collapse (24, 25).

Fig. 3 Marine geophysical evidence of ice-stream bed properties and retreat behavior.

(A) Location map for data panels [red boxes enlarged in (B) to (D)]; ice-stream flow configuration at maximum stage shown as white arrows (17). (B) Seabed bathymetry of the west Minch trough, showing GZWs 6 to 9. Sub-bottom acoustic profiler line [sections shown in (E) and (F)] highlighted in purple. Long-dashed line is boundary between soft and hard bed; short-dashed line is edge of “calving-bay” deglacial sediments. (C) High-resolution bathymetry of hard-bed region, south of 58°10′N. (D) High-resolution bathymetry of GZW15 cut by abundant large iceberg scours at time of floating ice-front collapse. (E and F) Sub-bottom acoustic profiles with interpretive panels showing internal character and sediment relationships of GZWs 8 and 9 (E) and GZW6 (F). TWTT, two-way travel time.


To place the newly mapped glacial landform record in a temporal context, we have derived a robust multiproxy chronology of ice sheet deglaciation. We obtained 26 beryllium-10 (10Be) terrestrial cosmogenic-nuclide (TCN) ages and three optically stimulated luminescence (OSL) ages from onshore glacial deposits unequivocally associated with ice-stream retreat (Figs. 1 and 4, figs. S5 and S6, and tables S2 to S12 in data file S1). Where possible, we obtained three or more TCN analyses from boulders on each glacial landform. The entire suite of 43 geochronological ages (including 14 previously published) was integrated using a Bayesian sequence model to produce modeled boundary ages that document a continuous absolute chronology of ice-stream retreat from maximal extent to terrestrial transition (~30 to 15 ka BP) (fig. S7 and tables S2 to S12 in data file S1). The modeled chronology displays excellent coherence along the whole 270-km-long ice-stream corridor with a high agreement index (26). In addition to these terrestrial ages, we obtained 18 radiocarbon dates from shells (mollusks and foraminifera) in marine cores taken within the path of the former MnIS. The dated shells were taken from sediment facies immediately overlying subglacial tills or from the base of cores where no subglacial material was recovered. The radiocarbon-dated shells independently cross-check the timing of ice-stream recession offshore NW Scotland by providing minimum dates for the timing of glaciomarine sedimentation shortly after the retreat of grounded ice (Fig. 1 and tables S2 to S12 in data file S1).

Fig. 4 Geochronology sampling sites.

A selection of the 12 terrestrial sampling sites analyzed using TCN and OSL to constrain the timing of ice-stream retreat in NW Scotland (see Fig. 1 for general locations). Field photographs relate to the location maps in the right-hand images. Location maps are hill-shaded digital terrain models, with a cell size of 5 m (contains NEXTMap GB elevation data). Note the different scales. All ages are in ka BP, with 1σ analytical uncertainties. Those in italic font are outliers. See Methods, Supplementary Materials, and supplementary tables for further information on field sampling, exposure-age calculations, and Bayesian modeling methodology.

Our geochronology shows that episodic, more predictable, grounding-line retreat between GZW1 and GZW9 was long lasting and occurred between 30.2 (±1.7) and 18.5 (±0.9) ka BP. Mean ice-stream retreat rates calculated from our Bayesian sequence model do not take into account the residence time spent at each of the GZWs, which introduces an uncertainty term. Unfortunately, the formation time of GZWs in modern and Pleistocene settings still lacks chronological control. However, previously published estimates of sediment fluxes at ice-stream grounding lines in Antarctica and Greenland (21, 23) indicate GZW formation of this size on decadal time scales rather than longer. Therefore, we assume that GZW formation was a relatively rapid process in this ice-stream system (<100 years), and given the uncertainties in the dating techniques (~0.5 ka), we calculate ice-stream retreat rates as net “reference” values without quantifying the duration of pauses during GZW formation. Mean net retreat rates of the ice-stream grounding line from 30.2 to 18.5 ka BP were between 10 and 20 m a−1. Following this, and shortly after the deposition of GZW9, the ice-stream terminus in the west trough experienced a period of uninterrupted grounding-line retreat over a flowline distance of 40 km. Submarine geomorphological, geophysical, and sedimentological evidence suggest that this ice-front collapse (equivalent to an area of ~500 km2) created a calving bay in the west trough between east Lewis and the mid-trough bedrock high. Our chronology places this first frontal collapse and associated ice-shelf loss between 18.5 ± 0.9 and 16.9 ± 0.6 ka BP.

Ice-front collapse and ice-stream evolution

Multibeam echo-sounder bathymetry data reveal well-preserved fields of large, closely spaced, cross-cutting iceberg scours in present-day water depths of 70 to 100 m (Fig. 3) in the east trough of the Minch. These iceberg scours are concentrated on the surface of GZWs 15 and 16 and are notably absent in deeper water (>150 m below present-day sea level). The iceberg furrows’ large size (length and width) but narrow water-depth range, strong linearity, preferred orientation, and predominantly u-shaped cross-profiles are the hallmarks of an abrupt ice-front breakup or mass-calving event (27, 28), rather than of cumulative iceberg scouring over a sustained period of time (fig. S9). We infer this breakup to have occurred when the retreating ice-stream calving front reached the grounding line (GZW15) during a second ice-front collapse—probably of a partially confined ice shelf or ice tongue. This sequence of events also briefly led to the formation of a deep calving bay between GZWs 15 and 17 as the grounding line temporarily stabilized and the system evolved from an ice stream with an extensive floating ice front to a large tidewater glacier. Evidence of proximal-glaciomarine calving-bay sediments are seen in the thickly draped acoustic facies and rhythmically laminated dropstone-rich muds recovered in cores adjacent to GZW15 (figs. S8 and S9). GZW17 is the innermost and smallest grounding-zone feature in the sequence that, we suggest, records transient grounding of the large tidewater-glacier ice front shortly after the loss of the ice-shelf or ice-tongue section.

The deepening trough bathymetry in the Minch, inshore of GZWs 16 to 17, is optimal to trigger and perpetuate rapid grounding-line destabilization via marine ice sheet instability (69) (Figs. 2 and 5). High-resolution bathymetry data show a pronounced overdeepening, no major bed obstacles, and an initially widening trough upstream of GZW17, indicating that once underway, grounding-line retreat would proceed in an uninterrupted fashion until the terrestrial transition. We estimate the collapsed ice-front area in the east trough, based on the absence of grounding-line landforms and the inshore bathymetry, to have been ca. 800 km2 (±10%) (Figs. 1 and 2 and fig. S10). Recent higher-order modeling experiments have simulated this accelerated retreat on the reverse-sloping section of the ice-stream bed ~180 to ~240 km from the continental shelf break, highlighting the importance of water depth and bed topography on grounding-line (in)stability and ice-stream retreat (29). These modeling experiments confirm that once the ice shelf was removed, the grounding line had no mechanism for restabilizing until it reached shallow water close to the present-day coastline (29). Our dating evidence places this second, larger, ice-front collapse event between 16.9 ± 0.6 and 15.9 ± 0.4 ka BP, with a mean net retreat rate of 70 m a−1 (Figs. 5 and 6). This ice-front retreat rate is considerably greater (>350%) than any frontal losses experienced by the ice stream up until this point. Furthermore, we suggest that this event, once underway, would have been rapid and irreversible—supported by modeling simulations (29)—perhaps resulting in significant ice-front losses (>400 km2) within a few years or decades. We see striking parallels between this unstable “threshold behavior” and the recent collapse of the Larsen A and B ice shelves fringing the Antarctic Peninsula. There, widespread ice-front collapse was immediately followed by dynamic adjustment of large tidewater glaciers in the resulting calving bays (25, 30). Satellite data of the Larsen B embayment demonstrated that increased ice discharge and accelerated mass loss in the wake of ice-shelf collapse were largely in response to a changed glaciological stress field and new boundary conditions (25, 30, 31). We envisage a similar situation in the Minch, albeit on a smaller scale, during the latter stages of ice-stream demise (~18.5 to 16 ka BP).

Fig. 5 Ice-stream retreat rates over time and on different strength beds.

Time-distance graph showing rate of ice-stream grounding-line retreat from 30 to 15 ka BP based on new multiproxy Bayesian-modeled chronology (see fig. S7). Inferred retreat behavior and ice-shelf presence/absence shown in the top panel; differing responses in west (W) and east (E) troughs are indicated. BMBs (A to J) (thin vertical lines) are Bayesian modeling boundaries. Bed elevation and simplified bed strength along trough also shown. Gray boxes indicate periods of ice-shelf collapse from geological and geomorphological evidence.

Fig. 6 Late Pleistocene climate and sea level reconstructions compared against MnIS retreat record and trough properties.

(A) Relative abundance (%) of Neogloboquadrina pachyderma sin. (Nps) foraminifera in core MD01-2461 (51.5°N)—a proxy for ocean climate (33). (B) Alkenone-derived sea surface temperatures (SST) from the western Mediterranean Sea (36°N) (34). (C) North Greenland Ice Core Project (NGRIP) δ18O record (35); (D) N. pachyderma sin. abundance (%) from deep-water cores MD04-2829 and DAPC2—a proxy for ocean climate at the latitude of the MnIS (59°N) (36). (E) Glacio-isostatic–modeled relative sea level (RSL) predictions for NW Scotland (58.5°N, Assynt; 58°N, Coigach; 57.5°N, Applecross) (40). (F). Mean global sea level curve (magenta line); blue shading is 95% probability envelope (39). (G) Changes in MnIS trough width at the grounding line over time (black rectangles), plotted at Bayesian model boundaries, with uncertainties shown. (H) Simplified bed strength at the grounding line and inferred timings of GZW formation. (I) Net ice-front retreat rate (in m a−1) between ~30 and ~15 ka BP, plotted at Bayesian model boundaries with uncertainties shown as horizontal bars. Cumulative ice-front retreat (in km) also shown on a separate axis (orange line). In this study (G to I), all data are plotted on the same time axis to aid visual comparison. Periods of abrupt Northern Hemisphere warming (Greenland interstadials; orange shaded boxes) (35) and abrupt global sea level rises (meltwater pulses) (39) are shown as vertical long-dashed and short-dashed lines and labeled. Duration of LGM shown (long black bar) (15). Duration of “local Last Glacial Maximum” in the British Isles (BIIS LLGM) also shown (short blue bar) (36).


Ice-stream retreat rates can be controlled by many factors—some external, such as climate and sea level, and others internal, such as bed topography, flow dynamics, and calving processes. Comparing our Bayesian-modeled ice-stream chronology and mapped grounding-line history from 30 to 15 ka BP with a suite of probable internal and external drivers allows a more detailed appraisal of these drivers over the last glacial cycle. Few climate records span the entire period, with terrestrial archives from Europe >55°N being particularly rare (32). In light of this, we compare our reconstructed record with ocean-climate proxy data from the northeast (NE) Atlantic and the Greenland ice-core record (33, 34, 35, 36), forming a latitudinal transect from 40°N to 75°N and balancing winter-biased and summer-biased records (Fig. 6) (37, 38).

Following the brief and abrupt climate oscillations of Greenland Interstadials 3 and 4 (29 to 27.5 ka BP), the climate of the NE Atlantic became very cold for ~3 to 4 ka during Greenland Stadial 3, as evidenced in ocean and ice-core records (3337). We note that retreat of the MnIS is already underway with the grounding line located on the mid-shelf (GZW2) by 27.5 ka BP, somewhat earlier than the traditionally viewed timing of the local LGM in the British Isles (i.e., 22 to 27 ka BP) (36). This indicates that the ice sheet in this sector reached its maximum extent before the LGM and, by implication, may have been driven by internal dynamics and was therefore less sensitive to climatic forcing. During the sustained cold phase of Greenland Stadial 3, ice-stream retreat was relatively slow and steady, resulting in the deposition of GZW3 and a prominent moraine in North Lewis (~25 ka BP). Episodic, predictable ice-stream retreat continued, with the grounding line stabilizing around the time of Greenland Interstadial 2 (GZWs 4 to 5). We suggest that the formation of large GZWs 6 to 8 at narrowings in the Minch trough, during the first half of Greenland Stadial 2, was promoted by the existence of a topographically confined ice-shelf buttressing mass flux and a weak deformable bed.

We find no apparent agreement between the net ice-stream retreat rate, or timing of grounding-line stabilizations, and climate forcing in the latter part of the deglaciation record (after 19 ka BP) (Fig. 6). According to our chronology and geomorphological evidence, a major change in ice-stream retreat behavior occurred between ~18.5 and 16 ka BP. The start of this key time period (~19 to 18 ka) is characterized by a marked cold phase, seen in the Greenland isotope record (35), and North Atlantic sea surface temperatures 10° to 12°C lower than present (3234, 36, 37). We note that, although NE Atlantic surface ocean waters experience a slight amelioration at 50°N to 60°N around 16.5 ka BP (33, 36), it is part of overall ocean cooling from ~19 ka BP culminating at ~15 ka BP at this latitude (Fig. 6) (33, 36). All our evidence, supported by recent modeling experiments (29), indicate that it was during this period, ~18.5 to 16 ka BP, that the Minch ice-stream/ice-shelf system underwent unpredictable and widespread change, in the form of grounding-line oscillations (GZWs 10 to 15) and ice-shelf loss, followed by accelerated, irreversible retreat (after GZWs 16 to 17).

Although global sea levels were close to their minimum at this time (15, 39), ice sheet [glacio-isostatic (GIA)] loading resulted in relative sea levels at or above the present day in NW Scotland (40). Empirical constraints on former sea level in NW Scotland are not available for the whole of this time period; however, the latest GIA models of the British Isles extend as far back as 20 ka BP (40). These model predictions show stable, higher-than-present, relative sea levels in NW Scotland from 20 to 17 ka BP followed by a fall in relative sea level from 17 to 15 ka BP (Fig. 6). Collectively, the evidence indicates only modest relative sea level changes (<10 m) over the latter stages of deglaciation (20 to 16 ka BP) owing to the complex but competing GIA and eustatic effects (40). There is no evidence of a rise in sea level in NW Scotland associated with the global meltwater pulse identified elsewhere at 19 ka BP (41).

The absence of a distinct atmospheric or oceanic warming event or sea level “trigger” at the time of the step change in ice-stream behavior ~18.5 to 16 ka BP argues against external forcing as the primary driver (Fig. 6). That notwithstanding, we cannot rule out the possibility of an external trigger (rapid atmospheric or ocean warming or sea level event) occurring on a time scale below the resolution of the currently available proxy records. However, the temporally synchronous but spatially differing behavioral response of the grounding line within the MnIS sector reinforces our view that catchment-specific internal (dynamical), rather than external, factors were the major drivers of change at this time.

We suggest that, as the Minch ice-stream/ice-shelf system retreated from the continental shelf, changes in trough width, depth, and bed strength led to marked variability in ice-stream recession rates. A step-change increase in recession, accompanied by threshold behavior, occurred between ~18.5 and 16 ka BP, as the ice-stream grounding zone encountered a major substrate change, from a soft to hard bed, and a divergent trough geometry. The marked increase in bed strength and trough bifurcation, around the mid-trough bedrock high—a large “sticky spot”—would have greatly increased frictional drag at the ice-stream bed triggering a fall in upstream ice velocity. Under these conditions, the ice stream/ice shelf would have switched from a low–basal traction state on a weak deformable bed, with most of the driving stress supported at the lateral margins, to a slower-flowing high–basal traction state with drag forces exerted widely across a strong rigid bed. We propose that this slowdown restricted ice flux to the downstream region leading to dynamic thinning and extension of the ice-stream/ice-shelf terminus region as it adjusted to new boundary conditions. Extensional (dynamic) thinning of the partially floating ice front increases buoyancy to a critical point whereby ice-shelf collapse becomes possible. These dynamically driven changes, loss of ice-shelf buttressing coupled with loss of grounding-line stability in a shoreward deepening bedrock trough, represent prime conditions for marine ice sheet instability (69).

We have reconstructed the retreat history and grounding-line movements of the Minch ice-stream/ice-shelf system, NW Scotland, over 15 millennia, in a well-constrained trough with known physiographic and geological bed properties. We identify a major behavioral step change around 18.5 to 16 ka ago—out of tune with external forcing factors—associated with the collapse of floating ice sectors and rapid ice-front retreat. We conclude that the step change in ice-stream recession was in response to new boundary conditions established during the latter stages of deglaciation in the Minch. A marked increase in bed strength, from soft deformable sediment to hard rigid bedrock, coupled with changes in trough width and a shoreward deepening bed profile, would have been instrumental in promoting instability of the ice stream’s floating terminus, which until then had provided force buttressing and a check on mass flux at the grounding line. Following the subsequent breakup of supporting ice shelves, our data show that ice-stream losses after ~17 ka BP were much more rapid under the new stress conditions. Driven by twin factors—(i) a major change in bed geology increasing subglacial traction and (ii) greater variability in trough width decreasing lateral support—ice-stream demise would have progressed at a faster, more unpredictable rate. Final collapse or “runaway” retreat ~17 to 15 ka BP would have been hastened by a reverse sloping bed (1° to 5°), enhanced calving, and, in the absence of an ice shelf, enhanced grounding-line sensitivity to changes in ice thickness and flotation threshold. The nonlinear response of the ice-stream grounding line to bed properties elucidated here underlines the future requirement for numerical ice sheet models to include sufficiently high-resolution 3D basal topography (<1 km) and bed strength (geological) properties to further reduce uncertainties in predicted ice sheet dynamical behavior.


Bathymetry data

Almost the entire seabed within the study area, focusing on the trough of the former MnIS, is covered by good-quality digital bathymetry data (fig. S1). The datasets have been compiled in this study to form a seamless elevation model of the submarine landscape around NW Scotland. The four main bathymetric dataset sources were as follows: (i) European Marine Observation and Data Network (EMODnet) Digital Terrain Model (DTM) data, (ii) UK Hydrographic Office (UKHO)/Maritime and Coastguard Agency (MCA) (INSPIRE) data, (iii) Olex data, and (iv) Natural Environment Research Council (NERC) cruise JC123 data.

The EMODnet DTM ( is a freely available harmonized bathymetric dataset integrating all best-available survey data for maximum areal coverage. Around NW Scotland, the EMODNet DTM data are gridded at a resolution of 0.125 arc min (approximately 200 m) and comprise single-beam echo-sounder and multibeam echo-sounder data with small areas of General Bathymetric Chart of the Oceans (GEBCO) 2014 data, on the outermost continental shelf, to fill any gaps (< 1000 km2). The raw single-beam data typically have a vertical resolution of 0.1 m and a variable data density (spatial resolution) of 100 to 200 m (equivalent to one point per 0.02 km2). These best-available moderate-resolution data have been regridded at a harmonized spatial resolution. They cover around 4400 km2, predominantly the outer and mid sections of the ice-stream trough (fig. S1), north of 58.3°N and west of 5.3°W. The single-beam data were collected by various vessels with various specification echo-sounder equipment.

The central and inner parts of the ice-stream trough are covered by high-resolution (hydrographic quality) bathymetric data. These datasets were collected using various specification multibeam echo-sounder systems between 2000 and 2017 as part of the UK Civil Hydrography Programme coordinated by the MCA. Geotiffs of the data can be freely accessed from the UKHO data portal ( Typically, these data were collected using hull-mounted multibeam systems operating at ~200 kHz, using 80 beams across a 120° swath and providing data coverage of approximately three times the water depth. Sound velocities were collected from a hull-mounted sound velocity probe. Vessel attitude, heading, and positional data were collected using sea-going motion sensors, a gyrocompass, and a dynamic global positional system. To remove data artifacts and to prepare basic mean bathymetric surfaces and combined-uncertainty evaluation surfaces, post-acquisition primary data processing was conducted by the UKHO. Subsequent data processing and regridding were carried out at the British Geological Survey (BGS) using CARIS HIPS and SIPS and Fledermaus software. The processed bathymetric datasets were viewed and analyzed in ESRI ArcMap 10.4 GIS software.

Olex AS single-beam echo-sounder bathymetric datasets were used in conjunction with high-resolution multibeam data to map the glacial geomorphology in areas where no EMODnet data or UKHO (INSPIRE) data were available (fig. S1). The Olex single-beam echo-sounder data are part of a global dataset, managed, compiled, and licensed for scientific research. The echo-sounder data have a positional accuracy of 10 m or less; vertical resolution ranges from 0.1 to 1 m. The dataset is spatially discontinuous but is more robust where multiple soundings have been conducted in the same area. Bathymetric surfaces were viewed and exported using the Linux-platform software (via the BGS-NERC license).

In addition to these three regional bathymetric datasets, multibeam echo-sounder data were collected onboard the NERC research vessel RRS James Cook in 2015. The extent of these data is shown in fig. S1. The multibeam echo-sounder data were collected using a Kongsberg EM710 70- to 100-kHz system capable of high-resolution imaging in shallower waters (5 to 1500 m). Data were imported into CARIS in blocks for processing and cleaning before being exported to Feldermaus for 3D visualization purposes.

Sub-bottom geophysical data

Two main sub-bottom dataset sources were used in this study: (i) BGS-NERC geophysical data, collected on various cruises between 1969 and 2007 by the BGS UK Continental Shelf Mapping Programme, and (ii) NERC cruise JC123 data. BGS-NERC data included more than 100 single-channel, analog acoustic (seismic) geophysical lines from the Minch and continental shelf around NW Scotland (fig. S1). The entire dataset has now been digitized, and screen-resolution imagery of the original paper records is available via a web portal ( The best information regarding marine geological sediment character and thickness was gained from higher-frequency (1 to 5 KHz) acoustic sub-bottom profiles taken with boomer and pinger systems, although much use was made of lower-frequency (100 to 1000 Hz) sparker profiles, especially in the North Minch and on the mid-to-outer shelf where the Pleistocene sediment sequence is considerably thicker (>50 m). Interpreted linework was exported into ArcGIS for assimilation with other geological and bathymetric data.

NERC cruise JC123 aboard the RRS James Cook in 2015 collected high-resolution acoustic sub-bottom data using a Kongsberg SBP120 sub-bottom profiler (fig. S1). The system has a sweep frequency of 2.5 to 6.5 KHz and a ping interval of 500 ms. Data were continuously logged in seg and raw form and visualized on a digital echogram. Differential GPS (Appanix POS-MV) was used as the primary positioning and motion sensor, while Seapath200 was used as the secondary system. Data were imported into Kingdom Suite software for processing and visualization. Interpreted linework was exported into ArcGIS for assimilation with other geological and bathymetric data.

Marine geological mapping and determination of bed substrate strength

Marine geological and geophysical investigations (42) coupled with numerous 2D seismostratigraphical interpretations have allowed the general thickness of Quaternary deposits, or the depth to bedrock, to be mapped on the continental shelf around NW Scotland (19, 20). These offshore Quaternary thickness maps show both Pleistocene (i.e., predominantly glacial, glaciomarine, and nonglacial marine) and Holocene sediments (exclusively nonglacial, marine) sediments as an amalgamated unit. In this study, we determined the ice-stream bed substrates at the time of the last glaciation (Late Weichselian; MIS2-3). Using existing and newly acquired geophysical sub-bottom profiles (figs. S1 and S2), sediment cores (figs. S2, S8, and S9), and archived marine geological data (fig. S1), we separated the deglacial stratigraphic component (deposited mainly as glaciomarine sediments during ice-stream retreat), and the Holocene component (deposited during nonglacial marine conditions), from the Late Weichselian subglacial and pre-Weichselian sediments (older than MIS3) (fig. S2). Note, for simplicity and brevity in the main manuscript, that these younger sediments are sometimes referred to as simply “deglacial sediments” to distinguish them from sediments deposited at the bed of the ice stream during Late Weichselian times. This classification has enabled a simplified geological substrate map to be generated for the ice-stream bed at the time of maximal Late Weichselian ice sheet glaciation (MIS2-3) (i.e., 22 to 28 ka BP) (36, 43).

We converted the substrate map into a first-order “bed strength map” based on the uniaxial compressive strength of the dominant substrate types present at the Late Weichselian ice-stream bed (table S1). The reference compressive strengths were derived by other workers through geotechnical testing of the same rock type, or equivalent lithostratigraphic units, sampled in NW Scotland and/or elsewhere in the UK (4446). Most unlithified silt and clay/mud-rich Pleistocene sediments on the continental shelf around Scotland are classified as “very weak” or “extremely weak”: with compressive strengths typically between 50 and 500 kPa (<<1 MPa) (19, 20, 46, 47, 48). Certain marine muds and glaciomarine sediment facies within the Pleistocene sequence in the outer Minch may have compressive strengths as low as 10 kPa (0.01 MPa). However, waterlain tills, subaqueous subglacial facies, and other clay-rich matrix-supported diamictons can have strengths of up to 1 MPa, depending on water content, clast density, and compaction. Coarse granular (sand-rich and/or gravelly) facies are difficult to classify in terms of compressive strength as they have very low cohesive and intergranular strength. By contrast, lithified rocks within the ice-stream catchment are predominantly classified as “strong,” “very strong,” or “extremely strong” with compressive strengths typically >100 MPa, with several lithologies [e.g., Lewisian gneiss and Eriboll sandstone (quartzite)] in the >250 MPa range (4548). In general, the weakest lithified rocks present at the ice-stream bed would have been the Jurassic sandstones in the inner Minch (between Skye and the East Shiant Bank), typically ranging between 50 and 100 MPa (4548). However, even these less strong (not weak) rocks have been further hardened in many places by contact metamorphism during the emplacement of numerous Palaeogene (Tertiary) igneous bodies (19), hence our general classification of the bedrock bed as strong to very strong (>100 MPa). It is probable that moderately strong to moderately weak (5 to 50 MPa) bedrock, within the oldest and youngest parts of the Jurassic sequence (Hettangian and Oxfordian/Kimmeridgian Formations), may have locally occurred at the bed of the MnIS during earlier glaciations. However, no seabed or subseabed outcrops of these rocks have been identified in the Minch; with only localized onshore exposures found on the Isle of Skye and southern Raasay (19). It is very probable that any such soft Jurassic rocks (partially lithified clays, weak shales, etc.), if originally present in the path of the MnIS, would have been preferentially removed by earlier Pleistocene glaciations before the Late Weichselian (MIS2-3). Our resulting maps and cross sections highlight the major geological boundary at the former ice-stream bed: with very strong/hard substrata to the south and weak/soft substrata to the north of ~58°10′N (Figs. 1 and 3 and fig. S2).

Geomorphological mapping

We mapped the glacial geomorphology preserved on the MnIS bed using the best available digital bathymetry data (fig. S1). This work supersedes all previous compilations (16, 17, 43, 49), as much of the bathymetric data used in this study was not available at that time.

To perform geomorphological analyses, we produced the following derived datasets in ArcGIS 10.4: (i) grayscale shaded relief maps (1× and 3× vertical exaggeration) with NW and SW illuminations, (ii) bathymetric surface contour maps (at 1-, 2-, and 5-m vertical intervals), (iii) slope (gradient) maps (0° to 90°), (iv) slope aspect maps (000° to 359°), and (v) profile curvature maps (in degrees m−1). To avoid mapping errors caused by azimuth bias (50, 51) and to ensure faithful landform representation, 2D cross-sectional profiling techniques were used (in Fledermaus and ArcGIS), alongside the first- and second-order derivatives of the bathymetric grid. Color and grayscale shaded relief maps (geotiffs) of digital elevation models were primarily used for 3D visualization purposes and cartographic presentation. All maps and datasets (onshore and offshore) were projected in Geographic Coordinate System GCS-WGS-1984 to minimize geolocation errors. Submarine glacial landforms within the former flow path of the MnIS were digitized in ArcGIS 10.4. Because of the large size of the study area (~10,000 km2) and differing resolution bathymetric datasets (see fig. S1), it was not possible to map all of the study area at the same level of detail, although mapping was typically conducted at 1:10,000 scale. In this study, we paid particular attention to the following: large transverse ridges within the ice-stream trough, other ice-marginal features (moraines, etc.), other evidence of ice-stream/ice-shelf grounding, and areas of hard substrate (bedrock) at seabed.

For computational convenience, large transverse ridges (mapped as GZWs) are referred to by their distance along the ice-stream trough relative to the trough centerline, starting at the continental shelf break datum (0 km; see Fig. 2 for examples). For the purposes of this work, the MnIS trough is defined by the 100-m present-day water-depth contour (except for in the outermost shelf west of 6°30′W where the 120- and 130-m contours are used, as the 100-m bathymetric contour does not come close enough to the shelf break to generate a closed polygon). The trough is defined at the shelf break by the 200-m water-depth contour. A generalized “trough outline” polygon was generated on the basis of a 1-km regridded surface (contour) model of best-available bathymetry. From this polygon, a centerline (equidistant from both trough margins) was generated using a script in ArcGIS.


We used TCN surface-exposure dating and OSL analyses to constrain the timing of deglaciation in the former MnIS flow path. Samples for TCN analysis were taken from glacially transported boulders within a relatively small area (i.e., 102 to 104 m2) and were assumed to be contemporaneous (i.e., with no evidence to suggest otherwise, adjacent boulders are assumed to have been deposited at broadly the same time during ice sheet wastage) (52). We report a total of 26 new cosmogenic 10Be surface exposure ages from seven carefully chosen sites situated at various key points along the flowline of the former MnIS. An additional 13 10Be ages, from two other sites, can be used to constrain the MnIS retreat chronology; these have been reported previously (53, 54) but are presented here with recalculated exposure ages. The combined dataset of 39 10Be samples is presented in table S2. Other published 10Be ages from three other sites are shown in Fig. 1 as supporting evidence of ice sheet retreat in NW Scotland. These were recalculated, as above, from previous studies (55, 56). These ages are not presented in tables S2 to S12 in data file S1 or used in the Bayesian modeling analysis as they lie outside the main flow path of the ice stream. OSL dating of proglacial and ice-marginal sediments has been used successfully elsewhere to constrain ice sheet margin behavior and the timing of deglaciation around the British Isles (57, 58). We report three new OSL ages from two carefully chosen sites to augment our 10Be TCN dataset (table S5). OSL samples were taken from glacigenic sediments directly relating to the margin of the former MnIS.

Field descriptions and site-specific details of each of the TCN and OSL dating sites are not presented here for reasons of brevity but can be accessed on the project data depository or obtained from the authors upon request. Samples for TCN analysis were taken from the uppermost surfaces of glacially transported boulders using a battery-operated angle-grinder and/or hammer and chisel by experienced field scientists wearing the appropriate protective clothing. We chose boulders resting on level surfaces, where possible, >100 m away from slopes (>5°); preference was given to large boulders with b axis > 1 m to minimize the potential for disturbance and snow cover. We sampled flat or gently rounded surfaces and avoided fractured or obviously weathered/degraded surfaces to minimize the potential for variable “erosion rates.” Care was taken to avoid any unnecessary surface damage to the boulders during sampling; sharp and prominent edges were removed after sampling to minimize the visual impact of sampling. Photographs were taken before and after sampling in most cases. Boulder dimensions were measured in the field to the nearest 0.1 m. Sample locations and elevations were recorded on a hand-held GNSS (Global Navigation Satellite System) unit (with GPS and GLONASS satellite acquisition). Positional data were checked against UK Ordnance Survey 1:25,000 maps and in a GIS database for accuracy. Any natural shielding from surrounding topography was measured in the field using a compass and clinometer and corrected for using the CRONUS-Earth online calculator (59) accessed on 01 October 2017 (

Sample thickness was measured using digital calipers. Individual rock fragments were measured from which a mass-weighted average thickness was calculated for each sample. All new samples reported here were crushed and washed at the University of Glasgow (School of Geographical and Earth Sciences). Following sample crushing and washing, quartz was separated from the 250- to 500-μm fraction using standard mineral separation techniques (60) and purified by ultrasonicating in 2% HF/HNO3 to remove remaining contaminants (mainly feldspars) and meteoric 10Be. The purity of the leached samples was assessed by measuring the aluminum content using flame atomic absorption spectrometry with bulk Al content considered a proxy for presence of feldspars. Be extraction was carried out at two independent laboratories housed at the Scottish Universities Environmental Research Centre (SUERC): the NERC Cosmogenic Isotope Analysis Facility and the SUERC Cosmogenic Nuclide Laboratory, using established procedures (61). The 10Be/9Be ratios of all samples (new and previously reported) were measured on the 5 MV accelerator mass spectrometer (AMS) at SUERC (62). Locational and analytical details for all samples are given in table S2. NIST27900 [10Be/9Be = 2.79 × 10−11] standardization was used in the AMS measurements. Measured ratios were converted to concentrations of 10Be in quartz (atoms g−1). Blank corrections ranged from 1.2 to 10.3% and were propagated in quadrature with attendant AMS analytical uncertainties. Analytical, chemistry, and blank data are included in table S2. See the Supplementary Materials for more information on age calculations, choice of production rates, and other external factors taken into consideration.

We sampled for OSL using opaque 30-mm-diameter tubes hammered into the chosen sand-grade sediment facies to prevent exposure to sunlight during field sampling. Samples were capped and labeled and kept in dark conditions at all times during analysis. External beta dose rates were determined for OSL dating using inductively coupled plasma mass spectrometry (ICP-MS) and ICP atomic emission spectroscopy (ICP-AES), while the external gamma dose rates were determined using in situ gamma spectrometry in the field (table S5). Water contents of 23 ± 5% were estimated considering the field and saturated water contents and environmental history for each sample; these are expressed as a percentage of the mass of dry sediment. Samples were taken from depths of 0.7 m (T8SKIG01), 9.0 m (T8SKIG02), and 2.0 m (T8GABB01); cosmic dose rates were determined accordingly (63).

To isolate coarse-grained quartz for OSL analysis of equivalent doses (De), each sample was first treated with a 10% (v/v) dilution of 37% HCl and with 20 (v/v) of H2O2 to remove carbonates and organics, respectively. Dry sieving then isolated the 212- to 250-μm-diameter grains, and density separation using sodium polytungstate provided the (quartz-dominated) fractions (2.62 to 2.70 g cm−3). The quartz grains were etched for 1 hour in 40% hydrofluoric (HF) acid to remove the outer portion of the quartz grains that was affected by alpha irradiation and also remove any contaminating grains of feldspar. After the HF etching, grains were washed in a 10% solution of HCl to remove any fluorides that may have been produced. Grains were finally mounted into 10 by 10 grids of 300-μm-diameter holes in a 9.8-mm-diameter aluminum single-grain disc for analysis.

All luminescence measurements were performed using a Risø TL/OSL DA-15 automated single-grain system equipped with a 90Sr/90Y beta source (64). A green laser was used to stimulate the grains for 1 s, and the OSL signal was detected through a 2.5-mm-thick U-340 filter and convex quartz lens placed in front of the photomultiplier tube. A preheat plateau test performed on 5-mm aliquots of sample T8SKIG02 was used to determine the preheat temperature (220°C for 10 s) used with a cut-heat of 160°C for the single-aliquot regenerative dose (SAR) protocol (65). OSL stimulation was performed at 125°C for 1 s, and the initial and background OSL signals were summed over the first 0.1 s and final 0.2 s, respectively. Typical decay curves and a dose-response curve measured for a single grain of quartz from sample T8SKIG02 are shown in fig. S6. Dose-recovery experiments performed on samples T8SKIG02 (ratio of 0.96 ± 0.03 and overdispersion of 15 ± 1%) and T8GABB01 (ratio of 1.01 ± 0.03 and overdispersion of 0 ± 0%) suggest that the SAR protocol was appropriate for OSL dating. The overdispersion values determined from dose-recovery experiments provided estimates of the amount of scatter caused by intrinsic sources of uncertainty that is beyond measurement uncertainties incorporated into a natural De distribution (66).

Bayesian sequence modeling

In general, our TCN exposure ages and luminescence ages from any given site show good agreement and form approximately normal distributions (fig. S7). We calculated χR2 statistics for all sites (except Tolsta Head; see below). In several cases, the experimental χR2 statistic exceeded the value expected (at 95% confidence). We tested for outliers using a uniform-phase Bayesian sequence model in OxCal 4.2 (26) run in an outlier mode (62, 67) to assess for outliers in time (t). Sample sites are arranged in a pseudostratigraphical and spatial sequence of ice-marginal retreat (68), i.e., they are arranged in the order by which they would have been deglaciated (Fig. 1). This prior model is informed by the reconstruction of the flow signature of the MnIS at maximum conditions (17, 49) and mapping of the offshore moraines (Fig. 1). We include all of the data presented in this model with the exception of the site at Stoer as we were unable to confidently determine its relative position within the sequence, vis a vis the sites at Garrabost and Rainish, without reference to its age determinations. In addition to the 10Be and OSL chronology, the Bayesian model is constrained by a suite of radiocarbon dates from an organic horizon underlying till at Tolsta Head, east Lewis (see table S4) (69). These provide terminus post quem for advance of the MnIS towards its LGM limit. All data were assigned a prior probability of 0.05 (i.e., 1 in 20) for being an outlier. The model was set up to assess for outliers in time (t) and uses a Student’s t distribution to define how the outliers are distributed and a scale of 100 to 104 years.

The geochronological data produced a sequence (fig. S7) with overall agreement indices of ~35%, below the >60% threshold commonly applied (26, 68). The modeled posterior outlier probabilities are shown in table S4. To improve the coherence of the sequence model, we iteratively increased the individual prior probabilities of samples where the model gave a posterior probability that exceeded our assigned initial value of 0.05. The refined model produced a conformable sequence with an overall model agreement index of 136 (fig. S7). The results, including modeled boundary ages, are summarized in table S11. These modeled boundary ages document the timing of retreat of the MnIS. Notably, the samples that initially gave posterior outlier probabilities exceeding the initial assigned prior probability of 0.05 came from sample groups (North Rona, North Lewis, and North Raasay), where the experimental χR2 statistic exceeded the value expected at 95% confidence. Excluding these samples would yield an acceptable χR2 for these sites. In other sites where the χR2 statistic exceeded the acceptable value, the Bayesian outlier analysis did not identify any samples as outliers (e.g., Cape Wrath). We attribute this to the fact that χR2 calculations use the analytical uncertainties of the surface exposure ages. Because our full retreat sequence also includes radiocarbon and OSL age determinations, we modeled the geochronology using the full, external uncertainties of the ages. As a test, we used the quality control criteria (52) to assign prior outlier probabilities to all geochronological ages with our dataset. In these cases, ages were manually removed until an acceptable χR2was obtained. This sequence shows good coherence (A = 227) but, more importantly, is in good agreement with our refined-age model. Note that, for discussion purposes, we use our initial refined-age model as it retains a larger number of geochronological dates and thus, we believe, yields a more conservative assessment of the associated uncertainties.


Supplementary material for this article is available at

Supplementary Text

Fig. S1. Overview of offshore datasets, covering the MnIS bed, used in the mapping compilations.

Fig. S2. Summary Quaternary geological cross-section and sediment cores from the MnIS bed (west trough).

Fig. S3. Bathymetric cross-sectional profiles of GZWs identified within the MnIS trough.

Fig. S4. GZW morphological metrics derived in this study.

Fig. S5. Probability density function plots of 10Be TCN exposure ages used in Bayesian sequence model.

Fig. S6. Single-grain De distributions determined for OSL dating of glaciofluvial ice-marginal sediments.

Fig. S7. Bayesian sequence modeling of chronology derived in this study.

Fig. S8. Submarine geomorphological and geological evidence of ice-shelf breakup.

Fig. S9. Diagnostic sediment facies within two key cores (JC123-031PC and JC123-021PC) (see figs. S2 and S8 for location within the west and east troughs, respectively).

Fig. S10. Conceptual model of ice-stream/ice-shelf threshold behavior conditioned by bed strength and bed topography.

Table S1. Rock strength classifications (according to BS 5930) for geological units at MnIS bed.

Data file S1contains the following supplementary tables:

Table S2. Sample, chemistry, AMS, and blank data for all 10Be samples presented and used in this study.

Table S3. TCN exposure-age production-rate comparison tables.

Table S4. Results of TCN exposure-age dating presented and used in this study and reduced χ2 analysis.

Table S5. Results from OSL dating, including the environmental dose rates to grains of quartz.

Table S6. Radiocarbon sample data and results of radiocarbon analyses presented in this study.

Table S7. Data for exposure age determination (CRONUS).

Table S8. Input for CRONUS online calculator.

Table S9. Input for CREp (Cosmic Ray Exposure program) for calculating exposure ages.

Table S10. Age model outlier-detection analysis.

Table S11. Bayesian sequence model (OxCal output).

Table S12. Summary of Bayesian model boundaries, ages, retreat rates and uncertainties.

References (7088)

This is an open-access article distributed under the terms of the Creative Commons Attribution license, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.


Acknowledgments: We thank the captain and crew of RRS James Cook (Science Cruise 123), the BGS Marine Operations team on JC123, and S. Morgan (University of Leicester) for operating the Geotek multisensor core logger on the cruise. Bathymetry data in Fig. 3B are from EMODnet. Figure 3 (C and D) contains UKHO bathymetry data provided by U.K. Maritime and Coastguard Agency who are gratefully acknowledged. Funding: This research was supported by a NERC UK consortium grant (BRITICE-CHRONO NE/J008672/1). Author contributions: T.B., D.F., C.D.C., R.C.C., and C.Ó.C. conceived the study, with C.D.C. as overall project Principal Investigator (PI) T.B., D.S., D.F., R.K.S., and R.C.C. conducted onshore fieldwork for TCN and OSL samples. D.S. and D.F. analyzed rock samples, measured 10Be concentrations, and calculated exposure ages. R.C.C. and D.S. undertook Bayesian modeling. R.K.S. and G.A.T.D. analyzed sand samples, measured dosages, and compiled OSL chronological data. T.B., C.D.C., M.H.S., S.L.C., D.D., R.C.C., D.H.R., and C.Ó.C. conducted offshore fieldwork collecting sediment cores and geophysical data. M.S. and S.L.C. conducted sediment core analysis and compiled sediment logs. T.B. and M.S. conducted x-ray CT (computer tomography) scanning of sediment cores. T.B. and D.D. compiled bathymetric data and conducted multibeam image analysis. R.C.C. and D.D. compiled geophysical acoustic (sub-bottom) data. T.B. interpreted geophysical data (sub-bottom and multibeam) and undertook mapping compilations and morphometric analyses. S.G.M. conducted radiocarbon analyses. T.B. led the manuscript writing and drafted all figures. T.B. and D.S. compiled the Supplementary Materials. All authors contributed to manuscript writing and editing. Competing interests: The authors declare that they have no competing interests. Data and materials availability: All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. Additional data related to this paper may be requested from the authors.

Stay Connected to Science Advances

Navigate This Article