South-central Laurentide Ice Sheet dynamics and the formation of proglacial Lake Vita during MIS 3

Abstract

An understanding of the growth and demise of ice sheets over North America is essential to inform future climate models. One poorly studied subject is the glacial dynamics during interstadial Marine Isotope Stage (MIS) 3 (57–29 ka). To better constrain the southern margin of the Laurentide Ice Sheet during this time period, we re-examined a stratigraphic sequence in southeast Manitoba, Canada, and provide robust evidence for advance and retreat of ice. Around 46.6 ± 5.1 ka (1σ error), fluvial sands were deposited under similar precipitation and significantly cooler summer temperatures than present-day. Ice then advanced south over the area, before retreating once again and a return to boreal forest and grassland conditions. The area was then covered by proglacial Lake Vita, dammed by ice to the north. Geochronology constraints indicate Lake Vita existed from ca. 44.3 ± 3.6 to 30.4 ± 2.3 ka (1σ error), although gaps in the optical and finite radiocarbon ages suggest either a lack of data or plausible temporary ice-margin advances during this time period. Ice covering most of Manitoba during MIS 3 is in line with global δ¹⁸O records, and glacially influenced sediment deposition in the Mississippi River basin.

Introduction

An understanding of the growth and demise of ice over North America is essential to inform future climate models (e.g., Batchelor et al., 2019; Gowan et al., 2021; Lowell et al., 2021; Dalton et al., 2022a). Due to the role that ice masses play in global sea level, reconstructing the behavior of past ice sheets from geological, geomorphological, and paleoecological archives is also of critical importance for predicting the behavior of modern ice sheets. One such event that is poorly understood is the dynamics of the Laurentide Ice Sheet (LIS) during interstadial Marine Isotope Stage 3 (MIS 3, 57–29 ka; Lisiecki and Raymo, 2005). Reconstructions of LIS margins have long suggested the possibility of significant reductions in the extent of the ice sheet during MIS 3, after initial growth during MIS 4 or 5 (Dredge and Thorleifson, 1987; Batchelor et al., 2019; Dalton et al., 2019, 2022a, b). However, an MIS 3 ‘minimum’ scenario, involving complete deglaciation of Hudson Bay, is difficult to reconcile with the detailed stratigraphy, nonfinite radiocarbon ages from wood and accompanying optical ages in the western Hudson Bay Lowland (Hodder et al., 2023), the timing of Heinrich events (Miller and Andrews, 2019), and weathering inputs into the Labrador Sea (Parker et al., 2022). Moreover, recent cosmogenic and non-bulk radiocarbon ages (Figure 1, Supplementary Table 1) suggest more of a ‘maximum’ MIS 3 LIS that may have extended down into the upper Mississippi River basin after ca. 42 ka cal BP and ca. 31 ka cal BP (Kerr et al., 2021), and south of the Great Lakes around 35 ka cal BP (Arnott moraine; Carlson et al., 2018; Ceperley et al., 2019, between 39.1 and 30.4 ka cal). This conceptualized maximum extent is supported by deposition of outwash sediment in the Mississippi River basin by ca. 34.5 ka (Carson et al., 2019) and deposition of the Roxana Silt loess 60–30 ka (Forman and Pierson, 2002; Markewich et al., 2011; Muhs et al., 2018).

Figure 1. The study area (Site A) and its relationship to the possible extent of the Laurentide Ice Sheet near the end of MIS 3. The reconstructed ice margin is from around 31.1–34.5 cal ka BP, assuming a 0.1 ka ¹⁴C error for the ice-margin line (Dyke et al., 2002). The dashed line in the A–D area shows the extent of Figure 2. Index map shows North America, the area shown on Figure 1 (black square) and the MIS 2 ice sheet limits at 21.5 ka from (Dalton et al., 2023, excluding Greenland and the northern extent of all ice sheets).

It is simplest to interpret the maximum southern extent of the LIS during MIS 3 as being to the north of vetted MIS 3 radiocarbon ages from organics in nonglacial sediments (Marshall et al., 2000; Dyke et al., 2002; Batchelor et al., 2019; Gowan et al., 2021; Kerr et al., 2021). It is also possible that ice advanced and retreated before sufficient organics accumulated (e.g., Halsted et al., 2024). The relative ages of ‘flow systems’ have been used as evidence for MIS 3 build-up of ice in Quebec and Ontario (Kleman et al., 2010, Dalton et al., 2022b), but more detailed geomorphic studies show the ‘flow systems’ that cover Manitoba are actually numerous separate palimpsest flowsets of different indicator types (striae, landforms, till fabrics) that cannot be used to constrain the absolute timing of glacial events (Trommelen et al., 2012; Gauthier et al., 2019; Hodder et al., 2023, 2024). All other ice-margin reconstructions are based on the work of Dredge and Thorleifson (1987), who originally interpreted a possible MIS 3 ice margin that approximately followed the margin of the Precambrian shield. This was because the shield was presumed to have a rougher (hard bed) that provided more friction and hence slower ice. Another reason was that some areas of the Precambrian shield are covered by multiple unoxidized till units, which the authors suggested were deposited during continuous Wisconsin (MIS 2–4) ice cover. This theoretical margin was modified by Dyke et al. (2002, fig. 3) to include a lobe of ice drawing southward into the Red River valley and abutting the Manitoba escarpment (Figure 1).

An important region for constraining the southern LIS margin and improving our understanding of glacial dynamics during MIS 3 is southeastern Manitoba (Figure 1, sites A–D). This region lies near the center of former continental ice sheets and stratigraphic records, where preserved, contain information that is critical for continental-wide ice dynamics. Here, we perform a detailed re-investigation of a key stratigraphic record from this region. Site A, situated along the Roseau River, preserves a remarkable stratigraphic record of tills and organic-bearing sorted-sediment units (~12 m exposed; Fenton, 1974). The uppermost intertill sorted-sediment unit consists of yellowish-brown sand and silt that overlie dark gray clayey silts (Vita Formation; Fenton, 1974). These sediments are organic bearing; wood from the sand yielded a conventional radiocarbon age of 33.8 ± 2.3 ka cal BP (31.5–35.9 ka cal BP, BGS-625; Morlan et al., 2000), in addition to nonfinite ages (Supplementary Table 1). Shallow (1.5–3.5 m depth) organic matter at three other nearby sites also returned conventional MIS 3 radiocarbon ages, for wood and charcoal, between 33 and 46 ka cal BP (Figure 1, sites B–G, Supplementary Table 1). To determine if these ages accurately reflect a retreat of the LIS from southeast Manitoba during MIS 3, we conducted a detailed investigation of the stratigraphy at Site A, including paleobotanical characterization and optical age determination of the intertill sorted-sediment units.

Near-finite radiocarbon reliability

Radiocarbon dating, like all dating methods, produces estimated age ranges that are calculated from organic matter. For near-finite MIS 3 radiocarbon ages (35–55 ¹⁴C ka BP), concerns about age validity include lab reproducibility (e.g., Bélanger et al., 2014; Ward and Clague, 2019) and the more problematic contamination by modern carbon (e.g., Reyes et al., 2020) or old carbon (e.g., Shotton, 1972; Nambudiri et al., 1980). This latter issue has led some researchers to conclude that near-finite ¹⁴C ages, especially on shells without specialized pre-treatment, greater than ca. 35–40 ka to be questionable (Walker, 2005; Douka et al., 2010; Miller and Andrews, 2019). A reconstructed limit for MIS 3 ice by Dyke et al. (2002) in Manitoba (Figure 1) was based on 10 ¹⁴C ages at or near the Saskatchewan border (Dalton et al., 2019), although six are ¹⁴C ages on organic-rich sediments that should be disregarded due to their high potential for contamination (Grimm et al., 2009; Bayliss and Marshall, 2019; Young et al., 2021). A lack of reproducibility and age reversals in apparently similarly aged stratigraphic sequences are also problematic (cf., Dalton et al., 2016; Miller and Andrews, 2019). Radiocarbon ages on charcoal from the Zelena site in western Manitoba (Figure 1, site F) are from a fossiliferous silt and marl zone at 3-m depth, which is underlain and overlain by different till units (Klassen, 1967). The marl contains ostracodes indicative of mesotropic or eutrophic lake conditions, which are similar to those of present-day lakes (Klassen, 1969). Two attempts to date this charcoal with conventional methods returned quite different ¹⁴C ages (41.9 ± 2.1 vs 27.9 ± 0.7 cal ka BP, Supplementary Table 1, Lowdon and Blake, 1968, 1973), but only the 41.9 ± 2.1 cal ka BP age was published. Given the contrasting ages, and the antiquity of the methodology used to obtain these ages, this site needs further work before the intertill sorted-sediment bed at the Zelena site can be assigned confidently to MIS 3. Radiocarbon ages on charcoal from paleosols in Saskatchewan (Bélanger et al., 2014) also show mixed ages and stratigraphic age reversals, indicating recycling of organics and unclear relationships (Figure 1, site H, Supplementary Table 1). Averaging the seven finite ages from paleosol 1 (0.5–1.6 m depth) provides a tentative ¹⁴C age of ca. 44.5 cal ka BP, overlying two nonfinite ages (Supplementary Table 1). Other MIS 3 charcoal ¹⁴C ages in their study ranged between 26.5 cal ka BP and 50.1 ¹⁴C ka BP (Site H, Supplementary Table 1). Hence, to confirm an interstadial designation identified by ¹⁴C ages, other dating methods supported by stratigraphy and paleobotanical datasets are needed (e.g., Palstra et al., 2021; Halsted et al., 2024).

Study area

The study area (Figure 1, site A) was overrun by Red River lobe ice stream during the end of the last glaciation (Dredge and Cowan, 1989; Patterson, 1997; Harris et al., 2020). This ice stream occupied the lowland between the Manitoba escarpment and the Precambrian shield (Elson, 1956; Clayton and Moran, 1982), and flowed south into Minnesota and Iowa, USA (Patterson, 1997). The head of the ice stream, source of ice, and duration of this ice stream are largely unknown. The ice stream remained in southern Manitoba at ca. 13.0 cal ka BP, although the timing of retreat is largely unconstrained and the study area may have been ice-free earlier (Gauthier et al., 2022). The entire region was covered by Lake Agassiz during the Emerson phase (11.57–11.29 to 10.69–10.34 cal ka BP; Young et al., 2021), and possibly during the earlier phases.

The studied portion of the Roseau River lies between ~252 m and 288 m above sea level and is situated 10.5 km east of the low-lying (240 m asl) Red River valley of central Manitoba (Figures 1 and 2). The rise in topography is due to the thickness of Quaternary sediment (70–90 m in this poorly studied area; Keller and Matile, 2021) and not bedrock topography (Figure 2).

Figure 2. West-east cross-section through southern Manitoba and the Roseau River study area highlighting the regional topography, simplified surficial materials, and drift thickness (modified from Matile and Keller, 2012).

Methods

Stratigraphy

To conduct a detailed re-investigation of the stratigraphy at Site A (Fenton, 1974), we visited six sections in the Roseau River area (Figure 3). Section exposures were cleared over a width of at least one meter and at least 30 cm into the section face to remove colluvium and expose in-situ sediment. Sections were described in detail and assessed for lateral variability. Lithostratigraphic units were defined on the basis of texture, color, sedimentary structures, the presence or absence of organic matter, and the nature of the contacts between each unit (section logs and lithostratigraphic descriptions can be found in figures S1 to S12). Given our different methodology, we chose not to fit our observations into the pre-existing correlated stratigraphy (Fenton, 1974), and instead developed a local stratigraphy.

Figure 3. Sections studied along, and next to, the Roseau River in southeast Manitoba, Canada, (a) shown on a hillshade derived from LiDAR, and (b) as stratigraphic sections. The study area is the same as ‘A’ in Figure 1.

Clast lithology

Clast-lithology counts were conducted on the 2–8 mm size fraction of clasts for 10 till samples and one gravel sample, to help identify the source area that the ice was flowing from when till was deposited (e.g., Lee, 2017). Simplified categories include Paleozoic carbonate and sandstone, intrusive igneous rock, greenstone and greywackes, and clasts of unknown source (quartzite, quartz, chert).

Clast fabric

To determine the shear-stress orientation imparted on diamict beds, we measured the strain signature of the sediment by analyzing pebble macrofabric (Benn, 1995). This was completed by measuring the trend and plunge of a-axes of clasts; data sets consisted of at least 29 clasts (rod, tabular-rectangle, or wedge-shaped) at each site. For each clast, the a-axis was at least 1.5 times longer than the b-axis and both the a-axis and b-axis dipped less than 60° (Supplementary Tables 2, 3). Clast-fabric measurement sites were chosen based on uniformity of the diamict, where no sand lenses or discontinuous bedding were present. Furthermore, clasts measured were selected from places where they were presumed to be free to rotate in the matrix at time of deposition (i.e., not in places where the diamict is clast-supported or where measured clasts are close to much larger clasts). At each site, a horizontal step was excavated at least 30 cm into the section face. Clasts were then carefully excavated and measured from within a ‘box’ consisting of three vertical faces of different orientation, over a maximum volume of 30 × 30 × 30 cm.

The strain signature of diamict interpreted as till is generally thought to be related to ice-flow orientation (Holmes, 1941). Clast measurements were first graphically displayed on equal-area, lower hemisphere projection stereonets using Rockware Sterostat v1.6.1. The principal eigenvector (V₁) and eigenvalue (S₁) of the clast fabric was calculated (Mark, 1973) and the modality described (Hickock et al., 1996). While the plunge direction of V₁ previously has been interpreted as the ‘up-ice’ direction (Mark, 1974; Kjaer and Kruger, 1998), this may only be true about 60% of the time (Andrews and Smith, 1970; Saarnisto and Peltoniemi, 1984; Larsen and Piotrowski, 2003; Gauthier et al., 2019). As such, we also used modality and eigenvalue, together with rose-diagram patterns, till-clast lithology, and surface geomorphology when assigning ice-flow direction to till fabrics.

Pollen analysis

Sediment samples were analyzed for pollen from 22 samples across four nonglacial units at section 115-21-002 along the Roseau River (Supplementary Table 4). In each sample, pollen was concentrated using standard laboratory techniques (Faegri and Iversen, 1975) along with use of a heavy-liquid solution (sodium polytungstate) to remove any remaining sediments from the residues (Zabenskie, 2006; Campbell et al., 2016). During processing, ceramic palynospheres were added to each sample to aid in estimating the pollen concentration (Kitaba and Nakagawa, 2017). Identification of pollen grains followed the key of McAndrews et al. (1973) and only samples containing more than 5000 pollen grains/cm³ and less than 10% broken and/or unidentified grains were considered appropriate for statistical analyses. For each of the samples satisfying these requirements, at least 150 grains of herb, arboreal, and shrub groups were counted. Pollen sums are based on tree, herb, and shrub groups. Following identification of pollen grains in each sample, we reconstructed paleoclimate variables using the modern analogue technique (Overpeck et al., 1985), which statistically compares the fossil pollen data to a modern-day dataset and extracts relevant climate parameters. The modern-day dataset consists of 4882 sites spanning most biomes across North America (Whitmore et al., 2005; Dalton et al., 2017). Previous analysis suggests this dataset is appropriate for reconstructing mean summer temperature (June, July, and August) and total annual precipitation (Dalton et al., 2017). Present-day climate data for the Roseau River area are interpolated from gridded climate data and estimated to be 18°C (summer temperature) and 520 mm (mean annual precipitation). Modern-day sites were considered analogues to the fossil interval if they had a squared chord distance dissimilarity of less than 0.25, as recommended for studies using a comparable approach (Williams and Shuman, 2008). Finally, paleoclimate reconstructions were generated using the R package ‘analogue’ (Simpson, 2007; Simpson and Oksanen, 2014) along with the three closest analogues 500 bootstrap iteration cross-validation.

Geochronology

Chronological constraints were obtained using both radiocarbon and optical dating. Organic matter was submitted for accelerator mass spectrometry (AMS) ¹⁴C dating to the A.E. Lalonde AMS laboratory (University of Ottawa). One charcoal and two freshwater-shell samples, encountered at just one of the seven studied sites, were submitted in 2021. All radiocarbon ages herein are 2σ median ages (cal yr BP) calibrated using CALIB 8.2 (Stuiver and Reimer, 1993). Ages from terrestrial samples were converted to calendar years using the IntCal20 Northern Hemisphere database (Reimer et al., 2020). For optical dating, quartz extracts of sediment were mounted on aluminum discs, each aliquot containing 50–100 grains, and equivalent dose (D _e) values were measured using a single-aliquot regenerative-dose (SAR) protocol (Murray and Wintle, 2000, 2003) following the same laboratory protocols described in Hodder et al. (2023; see also Figure 4). Environmental dose rates were determined by drying and milling a representative portion of the bulk sediment used for dating.

Figure 4. Luminescence decay curve (main graph) showing the ‘natural’ signal (red line), and a dose-response curve (inset) for sample 21-2 which is typical of those for all samples herein. Boxes indicate individual data points that are fitted with the dose-response curve. Note that the initial part of the luminescence decay is dominated by the desired ‘fast’ signal component. Measurements were made using a Risø TL/OSL DA-20 reader. Beta dose given using a ⁹⁰Sr/⁹⁰Yr source that delivered beta radiation at ~4.7 Gy/min to the aliquots. Following irradiation, aliquots were preheated at 220°C for 10 seconds, while the treatment following test doses (~2.3 Gy) involved heating them to 160°C at 5°C/second and then turning off the heat. Aliquots were stimulated with 56 mW/cm² of blue light and ultraviolet (~350 nm) luminescence was measured using an Electron Tubes 9235QB photomultiplier tube placed behind a 7.5-mm-thick Hoya U-340 optical filter. In each case, exposure to infrared (870 ± 40 nm) light at 130 mW/cm² for 100 seconds at 50°C was made before stimulation with blue light to reduce or eliminate any luminescence from contaminating feldspar. For each measurement, the final 20 seconds of the signal was subtracted from the initial 0.4 seconds, and this value (L_i) was divided by that measured from a subsequent test dose (T_i) to produce the normalized (i.e., sensitivity corrected) signal, which is plotted on the vertical axis of the dose response curve graph. For all the samples in this study, the dose response is best fitted by an exponential + linear function. The equivalent dose (D _e) is estimated by interpolation of the natural signal onto the dose response curve, as shown.

Results

Roseau River stratigraphy

This Roseau River area study has eight lithostratigraphic units and one allostratigraphic unit. The depositional record in the Roseau River area is best preserved and exposed at three sections within just 200 m of each other (sections 115-21-002, 115-21-005, and 115-19-001; Figure 3). Observations from these three sections were used to construct a stratigraphic framework and the other sections were correlated to this framework.

Lithostratigraphic unit A: till

Massive olive-brown diamict (unit A), >2 m thick, is exposed at the base of section 115-21-002 (Figures 3, 5a). The diamict is highly consolidated and has a silty-sand matrix, minor joints that are not stained, and is more compact than the diamict units above, with 3–5% clasts that are granule- to small pebble-sized and striated or faceted (Supplementary Table 5). The clast-lithology is dominantly carbonate (75.1 count % Paleozoic carbonate clasts, 24.5 count % Precambrian shield clasts, n = 1). Unit A diamict is interpreted as basal till, based on the texture, lack of stratification, mixed-lithology clast content, consolidation, and presence of faceted and striated clasts. A spread-unimodal clast fabric (S₁ = 0.68, n = 39) was measured from unit A, and is interpreted to have been formed by southwest-flowing ice (~219°, Supplementary Tables 2, 3).

Figure 5. The stratigraphically oldest sediments in the study area are exposed at section 115-21-002. (a) Unit A diamict and unit B silt; (b) lithofacies C-1 sand and gravel, lithofacies C-2 sand with detrital organics, and lithofacies D-1 diamict; (c) lithofacies C-2 sand with detrital organics, lithofacies D-1 diamict, and lithofacies D-2 diamict; (d) unit E silt with detrital and visible organics overlain by unit F, coarse-grained sand, fine-grained sand, and silt.. Scale bar is roughly 7 cm long, divided by cm.

Lithostratigraphic unit B: nonglacial pond or floodplain

Unit B consists of massive to horizontally bedded, gray-brown to gray silt, 0.9–1.3 m thick, that is highly consolidated and has a sharp, horizontal lower contact with unit A diamict (Figures 3, 5a, 6, and Supplemental Figure 4). Unit B is interpreted to have been deposited in a nonglacial quiet-water pond or floodplain environment based on the fine-grained texture and organic content of the sediment. The general lack of visible laminations may be due to the well-sorted nature of the sediment and/or the minerology (color).

Figure 6. Detailed stratigraphy of the lower intertill sorted-sediment (nonglacial) units at site 115-21-002. See Figure 3 for the location of this site.

Lithostratigraphic unit C: nonglacial fluvial

Gravelly sand (lithofacies C-1), 0.5 m thick, has a sharp lower contact with unit B silt at section 115-21-002 (Figures 3, 6). The matrix is fine- to coarse-grained sand, very friable, and contains 5–20% granule- to small pebble-sized clasts. Fine-grained sand (lithofacies C-2), 0.25–0.35 m thick, sharply overlies lithofacies C-1 (Figure 5b). This upper lithofacies C-2 is well sorted, has laminations that conform to the lower contact, is compact, and contains detrital organics. One sample of each lithofacies was taken to analyze for pollen and an optical age was determined for lithofacies C-2 (Figure 6). Unit C is interpreted to have been deposited in a fluvial environment based on the relatively coarse texture of the sediment, changes from gravel to sand up-section, and the organic content of the sediment.

Lithostratigraphic unit D: till

Two diamict lithofacies sharply overlie unit C sand at section 115-21-002 (Figures 3, 5c). The lower diamict (lithofacies D-1, 0.4 m thick) is dark gray with a sandy silty-clay matrix and 10% clasts that are striated and faceted (Supplementary Table 5). Unit D is massive, jointed, and consolidated. The upper diamict (lithofacies D-2, 1.3 m thick) sharply overlies lithofacies D-1 (Figure 5c), is grayish brown with a sandy-silt matrix, and has 10–15% clasts that are striated and faceted (Supplementary Table 5). It is massive and consolidated, although it is strongly fissile parallel to the slope, which may be the result of mass movement (Bed L, Supplementary Figure 6). The clast lithologies of both diamicts dominantly carbonate (~74.0 count % carbonate clasts, 26.0 count % Precambrian shield clasts, n = 2). Unit D-1 diamict is interpreted as basal till based on the texture, lack of stratification, mixed-lithology clast content, consolidation, and presence of striated and faceted clasts. A spread-unimodal clast fabric (S₁ = 0.60, n = 30) was measured from lithofacies D-1 and interpreted to have formed by south-flowing ice (~187°; Supplementary Table 3). Unit D-2 diamict is significantly sandier and the fissility makes interpretation of the specific till-type more difficult. More regional data are needed to definitively identify till units.

A third diamict lithofacies crops out at the bottom of section 115-21-010 (Figure 3). This diamict is dark gray with a silty-sand matrix and 3–5% clasts, many of which are striated and faceted (Supplementary Table 5). It is massive and compact but not jointed. The clast-lithology is dominantly carbonate (61.4 count % carbonate clasts, 38.3 count % Precambrian shield clasts, n = 1). Because it is the lower diamict at that section, and situated below organic-bearing, sorted sediments, we have included it here as lithofacies D-3. This diamict is also interpreted as till, based on the texture, lack of stratification, mixed-lithology clast content, consolidation, and presence of striated and faceted clasts. A spread-unimodal clast fabric (S₁ = 0.72, n = 29) was measured from lithofacies D-3 and interpreted to have formed by south-southwest-flowing ice (~193°; Supplementary Table 3). Given the similar characteristics and fabric, it could be the same as lithofacies D1.

Lithostratigraphic unit E: nonglacial pond or floodplain

Unit E consists of massive, brown (oxidized) to blue-gray silt (Figure 5d), 1.6 to more than 2.5 m thick, that sharply overlies lithofacies D-2 till or is at river level at three sections on the Roseau River (Figure 3). The silt is highly consolidated, well sorted, organic-bearing (disseminated and small gastropods), and smells like sulfur—indicative of a reduced environment. Seven samples of this silt were taken in 2021 to analyze for pollen (Figure 7). However, these pollen data could not be used for subsequent quantitative analyses since pollen grains were largely broken, which indicates significant reworking and poor preservation. Regionally, unit E silt was mapped at two additional sections, and thus extends along the modern Roseau River for ~2.5 km. The fine-grained texture and organic content indicate deposition within a nonglacial quiet-water pond or floodplain environment. The lack of visible bedding could be the result of the well-sorted nature of the sediment, its minerology (color), or that the bedding was destroyed by cryoturbation, or by pervasive shearing during subsequent glacial advance.

Figure 7. (a) Detailed stratigraphy of the upper sorted-sediment units at site 115-21-002; (b) the location of pollen samples; (c) charcoal within lithofacies F-2; (d, e) ¹⁴C dated freshwater-shell fragments. See Figure 3 for the location of this site.

Allostratigraphic unit F: proglacial to ice-marginal lake with intermittent glaciofluvial deposition

Due to the complex stratigraphy and lack of lateral continuity over as little as 200 m, unit F consists of all sediments above unit E silts and below unit G till or unit H gravel (Figure 3b). These sediments infill a topographic low near the top of section 115-21-002 that deepens on the northwest side of the section reaching a maximum thickness of 5.6 m (Supplemental Figure 6). The lowermost sediments consist of intercalated coarse-grained sand with 20–30% granules to small pebbles, silty sand, and laminated fine-grained sand with sparse disseminated organics. The contact between unit F and unit E is sharp (Figure 5d). At section 115-19-001, the lowermost sediments consist of a fine- to medium-grained sand that is well sorted, massive to horizontal or ripple-cross bedded, with sparse disseminated organics and rarer coarser beds that contain up to 5% clasts (Figure 8a). This is the dominant lithofacies throughout unit F and is present at all studied sections (Figures 3b, 8a–d). Eight samples of this fine-grained sand, sampled near the top of section 115-21-002, were analyzed for pollen (Figure 7a, b). These sandy lithofacies are interpreted to have been deposited in a nearshore lacustrine environment, as evinced by horizontal to ripple cross bedding, well-sorted sediment, detrital organic material, and regional setting at height of land.

Figure 8. Unit F contains (a–d) well-sorted, massive to horizontal or ripple cross-bedded sand, with sparse disseminated organics; (e, f) a few gravel and coarse-grained sand beds occur within the finer-grained sand at section 115-21-002.

Gravel and coarse-grained sand beds are also present in unit F (Figure 3, Supplemental Figures 5, 6). Massive to horizontally bedded sandy gravels contain 20–50% clasts (granule to medium pebble-sized, rounded to angular) and are matrix- to clast-supported and poorly sorted (Figure 8e, f). The clast-lithology composition of one gravel sample (8.9 m depth) is dominantly carbonate (75 count % carbonate clasts, 25 count % Precambrian shield clasts). Some gravel and sand beds contain scattered organic lenses, as well as rare charcoal and freshwater-shell fragments (Figure 7c–e).

Unit F at section 115-19-001 is particularly interesting, because it contains chaotically bedded diamict and gravel (Figure 9a–c). These sediments are variably injected into, ripped up within, and dropped into unit F sands. The diamict is brown and contains 10–40% granule to small pebble-sized clasts, and rip-up clasts of gravel. A very strong spread-unimodal clast fabric (S₁ = 0.85, n = 30, Supplementary Table 3) was measured from this diamict, which can form in a variety of genetic environments (Bennett et al., 1999). This diamict could be an immature till, based on its mixed-lithology clast content (70.4 count % carbonate clasts, 29.6 count % Precambrian shield clasts, Supplementary Table 5), presence of some faceted clasts, and strong clast fabric. It also could have been formed elsewhere, and then was frozen and ice-rafted to this site, which would account for both its strong clast fabric and its contorted bedding. This diamict was not observed at other sections, nor in previous stratigraphic investigations (Fenton, 1974), which suggests it is a localized lithofacies occurrence. The gravel is very poorly sorted, massive to weakly stratified to chaotic, and includes syngenetic cryoturbation structures, such as diapirs of clast-supported pea gravel (Figure 9b). The gravel-clast lithology is dominantly carbonate (2.6 m depth, 77 count % carbonate clasts, 23 count % Precambrian shield clasts, n = 1). A pod of the gravel, 4.5 × 1.5 m, appears dropped into the surrounding laminated sand (Figure 9c), and includes several granitoid boulders with draping and curved lower beds that are interpreted as dropstones (Figure 9d). Deposition of the gravels and diamict within the unit F sedimentary succession requires a change in sediment source (larger clasts) and higher energy environment(s) to transport the clasts. Both gravel and diamict contain Precambrian clasts in similar proportions, which indicates a partial northeastern provenance for the sediment. The contorted bedding, involutions, and diapirs at section 115-19-001 resemble cryoturbation types 4 and 6, which indicate cold conditions but not necessarily permafrost (Vandenberghe, 2006).

Figure 9. (a) Chaotically bedded diamict, (b, c) gravel, and (d) dropstones variably injected into, ripped up within, and dropped into unit F sands.

Lithostratigraphic unit G: till

Unit G, documented at three sections along the Roseau River, sharply overlies unit F sand (Figure 3). Unit G is a 2.2–2.7-m-thick diamict, consolidated, but not jointed, and classified as massive (G-1) or laminated (G-2) lithofacies that are separated by an erosional, horizontal to undulatory lower contact. This diamict is missing from the top of sections 115-21-002 and 115-19-001 but is present within a gully just 90 m to the south of those sections (Figure 3). The diamict is also absent from most surface sites in the surrounding region (Fenton, 1974; Matile and Conley, 1979; Mihychuk, 1997). The lower contact of lithofacies G-1 is undulatory over 0.25 m at section 115-21-005 and is injected into the underlying sand at section 115-21-010 (Figure 10a). Sheared beds, thrusted beds, and drag folds at the lower contact were also noted by Fenton (1974), who interpreted that the folds were formed by west-trending shear stress. Lithofacies G-1 diamict is variably gray-brown, brown, or light olive-brown, compact, and has a clayey sandy-silt or silty-sand matrix and 10–15% clasts that are striated or faceted (Figure 10a–c, Supplementary Table 5). The clast lithology is dominantly carbonate (63–77 count % carbonate clasts, 23–37 count % Precambrian shield clasts, n = 3). Lithofacies G-2 diamict is variably light olive brown or light yellowish brown, compact, has a clayey sandy-silt or silty sand matrix, and 10–15% clasts that are striated and faceted (Figure 10d, e, Supplementary Table 5). The clast lithology is dominantly carbonate (59–70 count % carbonate clasts, 30–41 count % Precambrian shield clasts, n = 2). The horizontal laminations are 0.002–0.005 m thick, with variable lateral thickness. The lower lithofacies G-1 includes sand and gravelly sand beds of different orientations at section 115-21-003, where a spread-bimodal clast fabric (S₁ = 0.56, n = 30) was measured and interpreted to have been formed by shear stress to the west-northwest or east-southeast (287–107°, Figure 3). Contrastingly, strong spread-unimodal clast fabrics at section 115-21-010 and 115-21-005 are interpreted to have been formed by shear stress to the south (181° [S₁ = 0.74, n = 31] and 169° [S₁ = 0.67, n = 30], Supplementary Table 3). Unit G diamict is interpreted as till, based on the texture, degree of consolidation, inclusion of striated or faceted clasts from multiple source areas, and moderate to strong clast fabrics regardless of structure.

Figure 10. Upper stratigraphy. (a–c) Unit G-1 massive diamict overlain by (d, e) unit G-2 laminated diamict, (f–h) which sometimes is overlain by unit H gravel and unit I sand.

Regional work has suggested that deposition of a till sourced from the northeast (Senkiw till) was followed by deposition of stratified sands (Bedford Formation) and then a till sourced from the northwest (Roseau till; Fenton, 1974; Teller and Fenton, 1980). Our limited fieldwork has not yet verified this pattern, and more data are needed. The concentration of Precambrian shield clasts in the regional-surface diamict is 7.9–37.5% (n = 39); both the 1st and 95th percentile values occur in the surface tills within 15 km of the Roseau River sections (Thorleifson and Matile, 1993; Gauthier and Hodder, 2023). This spatial pattern seems to support the idea of two or more surface diamicts that are present as a fragmented patchy mosaic of uncertain genesis.

Lithostratigraphic unit H: glaciofluvial gravels

Unit H, documented at most sections, sharply overlies unit F sands or unit G-2 laminated diamict (Figure 3). Unit H is a massive matrix- to clast-supported, poorly sorted gravel, 0.1–0.8 m thick, with a very fine-grained sand to granule matrix, and 10–50% clasts. The clast are subangular to rounded, commonly faceted or striated, granule- to boulder-sized, and of both Precambrian shield and carbonate lithologies (Figure 10f–h). The underlying lithofacies F-8 sands are folded at section 115-21-001 (Figure 10h). The nature of this contact, together with a lack of topography available to initiate a debris flow, suggests that unit H is glacially sourced. Unit H gravel, which occurs at or near surface through the study area, has been interpreted as a lag deposit due to thorough wave-washing of till by Lake Agassiz (Fenton, 1974; Mihychuk, 1997). The gravels sharply (Figure 10f–h) overlie sand at numerous sites, which argues against formation from the wave-washing of a till. Instead, we interpret the gravels as glaciofluvial in origin.

Lithostratigraphic unit I: glaciolacustrine sands

Unit I is a massive well-sorted fine-grained sand, 0–1.5 m thick, which is documented at two sections (Figure 3) and at the surface throughout most of the study area. Regionally, sand is discontinuous and generally found in beach ridges or within inter-ridge lows below 270 m asl (Mihychuk, 1997). As such, unit I sands were deposited during and after deglaciation, within glacial Lake Agassiz (Fenton, 1974). Regionally, wood collected from organics overlying thin gravel and till at ~300 m asl provided radiocarbon ages of 11.56 ± 0.30 ka (GSC-5296) and 10,100 ± 90 ¹⁴C BP (Morlan et al., 2000). It is unknown how the unit I sands, ~20 m lower in elevation, relate to those ages.

Pollen analysis

Sediment samples analyzed for pollen included silts (8 from unit B, 7 from unit E) and sands (2 from unit C, 8 from unit F) from two separate horizons (Figures 6, 7) at section 115-21-002 (Figure 3).

Unit B: nonglacial pond or floodplain

Fenton (1974) documented ~0.05 m of peat or organic-bearing soil within the silt bed of unit B, with low pollen and spore content, indicative of a closed coniferous forest. The low pollen content suggests a soil context with inadequate waterlogging to preserve pollen well, with at least partial drainage and aeration, within a closed-canopy coniferous forest. Peat was not encountered during our field work in 2022, but pollen grains were adequately abundant and well-enough preserved in two silt samples for further analyses (Supplementary Table 4). Our analysis confirms the presence of a coniferous forest and a nearby grassland (Figure 11). The pollen assemblage is comprised primarily of Picea (~15%), Poaceae (~25%), and Salix (~15%), with smaller amounts of Betula and Pinus (~10% each). A small presence of fern spores (~5% Polypodiaceae) suggests occasional openings of the forest canopy, and all samples contained Pediastrum, a green alga usually found in sediments of freshwater lakes. Paleoclimate reconstructions from these two pollen samples suggest cooler average summer temperatures of 16.1 ± 1.6°C, compared to present-day (18°C), along with similar total annual precipitation (455 ± 125 mm), as compared to present day (520 mm). Modern-day analogue sites were found in the grassland-boreal ecotone. The remaining six pollen samples processed from this unit contained more than 10% broken grains and enumeration was not pursued.

Figure 11. Pollen groups reaching at least 1% abundance from organic-bearing units B, C-2, and E. Unit F samples all contained >10% broken grains. Pollen-derived paleoclimate reconstructions are also shown for unit B. Vertical lines represent present-day climate in the region, which is 18°C and 520 mm/year, respectively. Organic content for each interval is also shown. Black dots to the left of the stratigraphic column represent intervals that were sampled for pollen (details within Figures 6 and 7).

Unit C: nonglacial fluvial

Only the uppermost fine-grained sand sample contained sufficiently preserved pollen for interpretation (Figures 6, 11, Supplementary Table 4). The assemblage is dominated by Picea (~15%), Betula (~25%), and Salix (~15%). This sample also contained an abundance of fern spores (~70% Polypodiaceae), which suggests a more open forest canopy than the earlier interval (unit B), along with wetland indicators Cyperaceae (~150%) and Sphagnum (~70%). The lack of herbaceous taxa suggests a paleoenvironment dominated by boreal forest, more homogeneous than the one represented by the samples of unit B. Paleoclimate reconstruction was not possible at this site owing to a lack of modern-day analogues.

Unit E: nonglacial pond or floodplain

Unit E silts are the middle ‘pond and floodplain deposits’ of the so-called Vita Formation (Fenton, 1974). Based on four grab samples examined at the Canada Center for Inland Waters, macrofossils include ostracodes, terrestrial and aquatic mollusks, plant material, and insects (Fenton, 1974). Those macrofossils suggested deposition within an oxbow lake or muddy stream bank, in a grassland or tundra environment, that was potentially cooler and drier than today (Fenton, 1974). Pollen concentration was below 5000 g/cm³, for six of the seven pollen samples, and there were too many broken or unidentified grains for appropriate quantitative analysis (Figures 7, 11, Supplementary Table 4). The sole sample with barely enough preserved pollen was dominated by Picea (~60%), Poaceae (~20%), and Abies (~15%). While this indicates the presence of a boreal forest proximal to the depositional environment, it also supports the presence of a nearby grassland, supporting at least partially the local environmental interpretation from prior macrofossil analysis. Paleoclimate reconstruction was not possible for unit E, owing to a lack of modern-day analogues, however five of the unit E samples (PM9-12, 15) contained Pediastrum, a green alga usually found in sediments of freshwater lakes or wetlands. Overall, pollen data suggest a boreal environment with nearby grassland.

Unit F: proglacial lake

In all sand samples (8/8), pollen concentration was below 5000 g/cm³ and there were too many broken or unidentified grains for an analysis of paleoclimate (Figures 7, 11, Supplementary Table 4). As a result, further enumeration of pollen from these samples was not pursued. Nevertheless, identified grains include typical boreal trees (Picea, Pinus, Betula, Salix), herbaceous plants (Poaceae, Ericaceae, Corylus), and several grains of Juglandaceae. From a qualitative standpoint, the pollen assemblage of this interval resembles that of modern to recent samples from the Hudson Bay lowland (e.g., Dalton et al., 2017), hinting at a similar vegetation community to the one observed today in that area.

Geochronology

Radiocarbon dating

Three samples of organic matter were collected for radiocarbon dating from lithostratigraphic unit F at section 115-21-002 (Figures 3, 7). A fragment of charcoal returned a non-finite radiocarbon age (>50 ¹⁴C ka BP, UOC-16071), while radiocarbon ages obtained on freshwater-shell fragments are 43.7 ± 1.7 cal ka BP (UOC-16077) from the same bed and 42.4 ± 0.5 cal ka BP (UOC-16076) from 0.56 m higher in section (Table 1, Figure 7). The shell ages may be too old, because a freshwater reservoir effect is expected to increase their uptake of old carbon (Gauthier, 2022; Rech et al., 2023). The dated fragments were too small to identify taxa. Both the charcoal and shell fragments are detrital and could have been transported from different source beds.

Table 1. Radiocarbon ages obtained from organic matter collected in this study.

^a Calibrated using Calib 8.2 and the IntCal20 Northern Hemisphere curve

^b The uncertainty value represents the 68% confidence interval, which is reported as the difference between the midpoint and the upper or lower limit, using the greater value.

^c cf., Gauthier, 2022

Optical dating

Four samples of water-laid sediment were collected for optical dating: one from unit C and three from unit F above the radiocarbon samples (Figure 3). Optical dating results, including environmental dose rates, radioisotope contents, water contents, sample location and depth, D _e, and age estimates, are found in Tables 2 and 3. Between 18 and 24 aliquots out of 24 aliquots measured, and 26 aliquots out of 46 measured, passed all SAR quality-control criteria for each sample (Table 3). For sample 22-2, 11 aliquots were rejected due to the recuperation exceeding 5%. Other aliquots were rejected due to failing recycling ratio tests, or the natural signal interpolated beyond the highest given dose on the dose-response curve, or because the shine-down curve showed signs of a significant slow component (in total 7 aliquots). In addition, one aliquot D _e value from each of two of the samples was rejected as an outlier (Figure 12), reducing the overdispersion of these samples from 29 ± 5% to 17 ± 3% for sample 21-1 and from 41 ± 6% to 37 ± 6% for sample 22-2. In this study, outliers were identified using the InterQuartile Range method (cf., Medialdea et al., 2014). The results of dose-recovery tests (Table 3) suggest that the laboratory protocol (Figure 4) used to estimate D _e values was appropriate. Ages were calculated using both a central age model (CAM) and minimum age model (MAM), using accepted aliquot D _e values (Galbraith et al., 1999; Galbraith and Roberts, 2012). All ages use AD2022 as the datum year. Ages calculated with and without outliers rejected are statistically consistent within 2σ. Ages were calculated using the MAM and the CAM, and in most cases ages found using both models are consistent at 2σ. However, we consider the MAM to give better estimates of the true depositional age for these samples due to our interpretation that they were deposited within a fluvial or nearshore lacustrine environments where the potential for partial bleaching of the sediment grains is relatively high. The use of MAM for sample 22-2 is especially appropriate because it has a relatively high overdispersion (OD) value of 37 ± 6% (Table 3).

Figure 12. Radial plots showing the distribution of the equivalent dose values (D _e) from individual aliquots and estimations of representative D _e values used for age determination calculated using either the central age model (CAM; solid line, blue shading) or the minimum age model (MAM; dashed line, gray shading). Values plotted within the shaded regions fall within 2σ of the weighted mean values. The data points shaded black were identified as outliers and were thus not included in calculation of the CAM and MAM D _e values.

Table 2. Optical dating sample water content, radioisotope concentrations, sample depths and calculated dose rates.

¹ D^w = Water content (mass water/mass dry sediment). For the dose rate calculations, each D^w value

² Radioisotope concentrations found using neutron activation analysis (NAA).included an uncertainty of ±10% (1s) to account for fluctuations in water content over time. As-collected water contents were used because the samples were collected from well-drained sand.

² d = Sample depth beneath the ground surface, measured from the center of the collection tube (5 cm diameter).

³ Ḋ _C = Cosmic ray dose rate, found using present burial depths and the procedure of Prescott and Hutton (1994); Ḋ _T  = total dose rate (that due to cosmic rays plus that due to β and γ radiation).

Table 3. Optical dating results for quartz samples in the Roseau River area

¹ N = numbers of aliquots accepted/number of aliquots measured

² C = Determined using the central age model (CAM)

³ M = Determined using the three-parameter minimum age model (MAM)

⁴ DR = Dose recovery ratio (given dose/recovered dose)

Note: equivalent dose and OD values were calculated using the software of Liang and Forman (2019)

Using the D _e values derived using the MAM (Figure 12, Table 3), sample 22-2 from unit C-2 yielded an optical age of 46.6 ± 5.1 ka, while samples 22-1, 21-1, and 21-4 from unit F (Figures 8b, 8c, and 8d) gave optical ages of 44.3 ± 3.6 ka, 36.2 ± 2.7 ka, and 30.4 ± 2.3 ka, respectively (Table 3). Uncertainties are analytical only and are quoted as ± 1σ.

Discussion

The Roseau River exposes sediments that were deposited during at least two different ice-free phases and three different glaciated phases. There is repeated deposition of till (units A and D), organic-bearing quiet-water pond or low-energy floodplain sediments (units B and E), and organic-bearing fluvial or lacustrine sediment (unit C and unit F), followed by subsequent glaciation and retreat (Figure 13). Throughout the examined record, only a few pollen samples contained sufficient grains for paleoclimatic reconstruction, and in most places modern analogs are lacking. Nonetheless, the stratigraphy indicates a repeated change from till to sorted sediments containing pollen that suggests coniferous forest with nearby grassland, and wetland to boreal forest (Supplementary Table 5). The lower quiet-water pond or low-energy floodplain sediments (unit B) were deposited during significantly cooler summer temperatures than present-day and with similar total annual precipitation. Cooler temperatures support deposition during an interstadial, which is consistent with the MIS 3 optical age of 46.6 ± 5.1 ka (1σ, sample 22-2) for the overlying organic-bearing fluvial sands (lithofacies C-2, Table 3, Figure 13). Deposition of till by south-trending ice-flow (lithofacies D-1, D-2) directly above the ca. 47-ka organic-bearing fluvial sands, indicates that ice was active along the southern margin of the LIS at least once during MIS 3 (Figure 3). Glaciation of the Roseau River area during MIS 3 was short-lived, however (ca. 2300 years, discounting error ranges), as evinced by the 44.3 ± 3.6 ka optical age (1σ, sample 21-2) obtained on the overlying organic-bearing sands at the same section (unit F, Table 3, Figures 3, 13).

Figure 13. Summary of stratigraphic relationships in the Roseau River area, based upon the few sections studied and the chronological data (wood and quartz grains from sand).

Near-surface fine-grained sands, likely belonging to unit F, are documented between ~258 and 288 m asl (Figure 14), over at least 850 km² (Fenton, 1974; Keatinge, 1975; Mihychuk, 1997). Similar organic-bearing silts and sands, ~22.4 m thick, were noted below till in a borehole near Pansy (Figure 14; Matile et al., 2023). While the paleotopography of this paleo-lake is unknown, modern drainage is towards the northwest where waterways drain into the Red River, then ultimately into Hudson Bay towards the northeast. If the drainage was similar, for an MIS 3 lake to form at the height of land in the study area, it must have been impounded by ice to the north. This happened in the late Pleistocene–Early Holocene too, when retreating ice blocked this natural drainage pattern forming Lake Agassiz (Thorleifson, 1996). This MIS 3 proglacial lake, we name Lake Vita in deference to earlier stratigraphic work (Fenton, 1974). Paleoclimatic inferences are scant for this lake, owing to exceedingly low pollen concentration and a high percentage of broken or unidentified pollen grains. However, the few well-preserved grains suggest a vegetation assemblage somewhat different from that observed today, with more open grassland, more boreal elements, and the possibility of limited local peat formation, all of which is consistent with a dynamic peri-glacial environment. It is possible that the cold, periglacial environment, and high rates of sediment influx during this interval contributed to low pollen production on the landscape, low concentration in the sediments, reworking and breakage of pollen grains, and poor preservation. Periodic deposition of gravel and diamict within the unit F sands requires both a sediment source for the mixed-lithology clasts and fluctuating energy environments. We suggest the gravels and diamict were likely deposited in Lake Vita when the oscillating LIS margin was situated closer to the Roseau River area and provided increased meltwater discharge and/or icebergs.

Figure 14. Radiocarbon ages on wood, charcoal, or shell (Table 1, labels on Figure 1 and in Supplementary Table 1) and optical ages on sand from unit F (Table 3) in the Roseau River area of southeastern Manitoba. Review of data collected during previous fieldwork, where available, has identified multiple sites with confirmed (yellow dots) and likely (pink dots) Lake Vita sediments at or near the surface. The background image is a hillshade generated from LiDAR digital elevation models (Manitoba Government, 2020). Streamlined landforms in the area were formed during retreat of the Red River Ice Stream at the start of deglaciation (cf., Gauthier et al., 2022).

Timing of Lake Vita

Initiation and duration of Lake Vita can be constrained by the three radiocarbon ages and four optical ages (Figure 15). The new optical age determinations on nearshore lacustrine sands suggest that Lake Vita may have existed for ca. 14,100 ± 3000 years at 1σ error (n = 3, Table 3, Figure 15). This suggests a surprisingly long duration of Lake Vita, although it is tentatively supported by the range of conventional radiocarbon ages in the same region from unit F-equivalent sediments (31.5–42.1 cal ka BP at 2σ error, Figures 14, 15; n = 3, Supplementary Table 1). Small gaps in the regional ages could be due to either the low number of data points, or additional short-lived glacial advance(s). It is unlikely that a proglacial lake could remain stable in one place for ca. 14,100 years, given the inclination of ice to surge at the margins of proglacial lakes (e.g., Quiquet et al., 2021; Hinck et al., 2022). Nonetheless, our available chronological data currently indicate that Lake Vita did exist, likely with variable spatio-temporal configurations and water levels, between roughly 44.3 ± 3.6 ka (1σ, sample 21-2) to 30.4 ± 2.3 ka (1σ, sample 21-4).

Figure 15. Distribution of cal BP radiocarbon ages on wood/charcoal (Table 1, Supplementary Table 1) and optical ages on sand (Table 3), from sediments presumably correlative to unit F. Bars shows 2σ error on radiocarbon ages and 1σ and 2σ error on optical ages. Nonfinite (‘greater than’) radiocarbon ages are shown to demarcate their possible age and depth, acknowledging that the written age is a minimum equivalent to the limitations of the method at the time it was obtained. A short-lived glaciation (green line) occurred sometime between ca. 46.6 and 44.3 ka, according to the stratigraphy at section 115-21-002 (Figure 3).

Extent of the LIS during MIS 3?

Whether ice advanced over southeast Manitoba at 44 ka and 30 ka is uncertain, but radiocarbon and optical ages together suggest that this was a long-lasting ice-free period. The global marine δ¹⁸O record indicates that MIS 3 was a time of reduced global ice volume (Lisiecki and Raymo, 2005) and modeling suggests that fluctuations in global sea level during this time were almost exclusively controlled by the LIS (Gowan et al., 2021). Our work on the Roseau River sites, therefore, offers constraints on a restricted LIS during this interval. More broadly, our work comes at a time of renewed interest in empirically constraining continental ice masses during MIS 3. For example, evidence from Scandinavia and Britain suggest ice-free conditions in those regions during MIS 3 (Finlayson et al., 2010; Kleman et al., 2021) and that the Eurasian ice sheet contributed ~20 m of global eustatic sea-level rise during MIS 3 (Mangerud et al., 2023). As noted by Mangerud et al. (2023), this quantity of melt represents almost the entire amount of sea-level rise in some reconstructions (Spratt and Lisiecki, 2016; Pico et al., 2017).

The issue of North American ice extent during MIS 3 is, however, complex. Notably, there are also data that support an LIS ice maximum during MIS 3: the Mississippi River drainage basin contains glacially influenced outwash deposition dated to ca. 34.5 ka (Carson et al., 2019) and loess accumulation between 45 and 33 cal ka (Roxanna silt; Grimley and Leigh, 2022). Sediment contribution from the Des Moines Lobe to the Roxanna loess was interpreted to be minimal (Dendy et al., 2021), suggesting that the Mississippi River basin was mainly influenced by ice sourced from the Quebec–Labrador sector of the LIS. However, there is still the outlying problem of MIS 3 loess in the Missouri River valley, which seems to have a more northern source.

Of considerable relevance to our work along the Roseau River are recent findings of a potential MIS 3 ice advance into the Midcontinent United States. Along the southern margin of the LIS, Kerr et al. (2021) proposed that ice streamed over Minnesota and into Iowa twice during MIS 3. Using stratigraphy and radiocarbon probability density curves on wood ages in tills, they suggested that ice was present in Iowa between ca. 40–49 ka and 29–34 ka. Our stratigraphy and ages now confirm ice crossed over southern Manitoba sometime between 47 ± 5 ka and 44 ± 4 ka (1σ, 22-2 and 21-2, respectively), although there is no evidence for a second advance between 44 ± 4 ka and 30 ± 2 ka (1σ, 21-2 and 21-4, respectively). If both datasets are correct, it may be possible that Lake Vita existed on the periphery of an ice stream that was active in the topographic low between our study site and northern Iowa (Figure 1). For now, our data confirm that the MIS 3 margin drawn by Dyke et al. (2002) is reasonable for ca. 37–30 ka in southeastern Manitoba (Figure 1).

A major question is whether MIS 3 ice flowing south near the Manitoba–US border was sourced from the Keewatin or the Quebec–Labrador sector of the LIS. Using two clast fabrics, we interpret that the till deposited during MIS 3 (unit D) was sourced from the north or north-northeast (187°, 193°; Figure 3). Depending on how close ice from the Quebec–Labrador sector was to Manitoba at that time, ice could have easily crossed the Precambrian shield and traveled south within what is now a topographic trough (Figure 1). Kerr et al. (2021) suggested that the Iowa MIS 3 advance was from the Keewatin sector, because the till contains abundant shale, which presumably was sourced from Mesozoic rocks of the Manitoba escarpment (Figure 2; Johnson et al., 2016). Our new work shows that tills in the Roseau River area do not contain shale, but the lateral extent of the presumed ice stream that deposited this till has not been mapped. If this ice stream extended to the Manitoba Escarpment (a natural topographic barrier; Figure 1), then this same ice-flow event could be responsible for the shale content observed within till in Iowa. As such, compositional differences could be due to the migration of an ice-stream catchment area over time—and are not an indicator of the source sector.

MIS 2 ice margin?

The last glaciation in south-central Canada occurred during MIS 2 (29–14 ka; Lisiecki and Raymo, 2005), and was the only time in the Quaternary when the LIS completely covered the Canadian prairies (Jackson et al., 1999; Trommelen and Levson, 2008; Paulen et al., 2021) and coalesced with the Cordilleran Ice Sheet (cf., Young et al., 1994; Burns, 2010). The timing of this advance is largely unknown, but it is constrained by a ¹⁴C age of 28.6 ± 0.45 cal ka BP from wood in gravels under till in northeast British Columbia close to the Rocky Mountains (Beta 183598, 24.4 ± 1.5 ¹⁴C ka; Levson et al., 2004) and radiocarbon ages on bones from the preglacial valleys of central Alberta (Edmonton) between > 41 ¹⁴C ka BP and 25.6 ± 0.8 cal ka BP (Young et al., 1994; latter is AECV:1664c, 21.33 ± 3.4 ¹⁴C ka). There are no conclusive ages on sediment or organics in Saskatchewan or Manitoba between 30.4 ± 2.3 ka (1σ, 22-2, this paper) and 13.0 ka (UCIAMS-101445, 11.09 ± 0.4 ¹⁴C ka; Teller et al., 2020). Two inconclusive ¹⁴C ages include a 27.9 ± 0.7 cal ka BP charcoal age near the Saskatchewan border that also dated to 41.9 cal ka BP (site G, Figure 1, Supplementary Table 1; Lowdon and Blake, 1968, 1973), and a 26.5 ± 2.8 cal ka BP age based on charcoal from a dune paleosol that is situated just 0.2 m above a charcoal age of 42.3 ± 0.9 cal ka BP (site I, Figure 1, Supplementary Table 1; Bélanger et al., 2014; Dalton et al., 2019). As such, the position to where the formerly streaming ice margin retreated between a potential ca. 29–34 ka advance into Iowa (Kerr et al., 2021) and the final advance of ice during MIS 2 is chronologically unconstrained. The second advance of ice denoted by Kerr et al. (2021) possibly did not occur during MIS 3, but rather at the start of MIS 2 during Heinrich Stadial 2 (26.6–23.6 ka; Heath et al., 2018). Indeed, the Green Bay Lobe (Quebec–Labrador sector) in eastern Wisconsin was within several kilometers of its MIS 2 maximum as early as 24.7 ka (Carson et al., 2020). A third (and final) advance of the Des Moines Lobe to its maximum position occurred by 16.2 ± 0.3 ka, near the end of Henrich Stadial 1 (19.3–15.3 ka; Heath et al., 2018).

Future Work

A key avenue for future work will be better understanding the dynamics of Lake Vita. We have identified part of this extensive proglacial to ice-marginal lake, which existed at the southern margin of the LIS during MIS 3, possibly for a duration of 14,100 years. Ice-sheet models are now beginning to include ice-marginal lakes (Quiquet et al., 2021; Austermann et al., 2022; Hinck et al., 2022), and it would be interesting to study whether Lake Vita could have existed at the periphery of an ice stream, or whether ice would have had to completely retreat back north of the Manitoba–US border. As such, we foresee our work on Lake Vita to be of interest to numerical modeling and glacio-isostatic adjustment groups.

A second avenue for future work is potentially correlating Lake Vita sediments laterally and understanding the extent of the lake basin, while also investigating the occurrence of other dynamic proglacial lakes that likely formed during the LIS retreat in MIS 3. For example, rare glaciolacustrine sediments were also encountered below a till, presumed to be Late Wisconsinan, in a borehole in southwestern-most Ontario (Bacj, 1991). Although just 0.1 m thick, understanding the occurrence of additional proglacial lakes that likely formed during MIS 3 will help inform future reconstructions of the LIS during this critical time period. Related, the effects of large, highly dynamic proglacial lakes on establishment of local vegetation must be investigated further using pollen and other paleoecological proxies.

A third line of future work is understanding the mechanism and effects of the multiple ice advances during MIS 3, as first suggested by the work of Kerr et al. (2021) and now supported by our work in Manitoba. What dynamics or instabilities of the ice sheet resulted in this ice stream? Were these ice advances sourced from the Keewatin or the Quebec–Labrador sector of the LIS? Results from this work could also have implications for the formation of the Roxana Silt, which is a loess deposit in the Midwestern US that is thought to have resulted from ice advances into the Mississippi River watershed during MIS 3, but remains enigmatic because there was no evidence of an ice advance during that interval (Winters et al., 1988; Johnson and Follmer, 1989; Leigh, 1994; Forman and Pierson, 2002). While a recent detrital-zircon study suggested little contribution from the Keewatin-derived Des Moines Lobe to Roxanna silt (Dendy et al., 2021), our new work may require revisions to the larger-scale regional conclusions about the effects of ice sheets and ice dynamics over time.

Fourth, optical dating has provided valuable information about the timing of events discussed in this paper. However, the efficacy of linear + exponential dose responses, as was observed in the samples reported here, to dating quartz in general is still being tested (see the discussion in Hodder et al., 2023). Nevertheless, consistency among our age information, radiocarbon ages, and paleoecology suggests that our optical ages are valid.

Conclusions

New stratigraphy, optical ages, and paleoenvironment data provide robust evidence for advance and retreat of ice along the southern margin of the Laurentide Ice Sheet during MIS 3. An optical age on sands underlying till indicates that ice had retreated north of the study area before 46.6 ± 5.1 ka. Those sands were likely deposited within a fluvial environment, because they overlie silts deposited in a pond or floodplain. Pollen indicates that the sands were surrounded by boreal forests and wetlands, while the underlying silts were surrounded by more open forest with grassland when the area had significantly cooler summer temperatures than present (16.1 ± 1.6°C) and similar annual precipitation (455 ± 125 mm). Temporary glaciation followed, as evinced by till interpreted to have been deposited by south-trending ice flow. Ice then retreated, and ~1.5–2.5 m of sediment was deposited in a quiet-water pond or low-energy floodplain environment and overlain by ~5.0–7.0 m of sediment deposited in a shallow ice-marginal lake environment. Formation of one or more shallow lakes, at the height of land, requires impoundment of drainage by an ice margin at least 40 km north of the Manitoba–Minnesota border; similar to Lake Agassiz during the Holocene. The quiet-water sands were deposited in a large proglacial lake, herein named Lake Vita. Organic-bearing gravel (detrital charcoal >50 ¹⁴C ka BP) and diamict are interbedded with the quiet-water sands, which requires both a sediment source for the clasts and higher energy environments. The gravels, which are of mixed composition, were likely deposited in a proglacial environment when the Laurentide Ice Sheet advanced closer to the lake and provided increased meltwater discharge and/or icebergs. Current geochronology constraints indicate that Lake Vita existed from ca. 44.3 ± 4 to 30.4 ± 2 ka (optical ages from quartz, 1σ error), although time gaps suggest either a lack of data, quiescent period of sediment deposition, or possible temporary ice-margin advances during this time period. The final ice advance across the study area occurred after 30.4 ± 2.3 ka. Our new data confirm the Roseau River area was deglaciated during MIS 3 yet support the idea of ‘maximum’ nearby Laurentide ice sheet at this time.