INTRODUCTION
Ice wedges and associated tundra polygons are a ubiquitous feature of permafrost terrain (French, Reference French2013). The density, size, and overall activity of ice wedges vary as a function of climate, vegetation, and polygon morphometry (Mackay, Reference Mackay1974, Reference Mackay1992, Reference Mackay1993; Fortier and Allard, Reference Fortier and Allard2005; Abolt et al., Reference Abolt, Young, Atchley and Harp2018). Early research on ice wedges focused on their growth and the mechanics of cracking in relation to winter conditions (Lachenbruch, Reference Lachenbruch1962). Ice wedges tend to develop when winter ground surface temperatures are < −10°C and when the tensile strength of the ice is exceeded during thermal contraction following rapid cooling (Mackay, Reference Mackay1993; Fortier and Allard, Reference Fortier and Allard2005). Most ice wedges crack between mid-January to late March near the center, with a crack width in the 2–4 cm range and a crack depth that approximates the depth of the wedge (Mackay, Reference Mackay1974). The frequency of ice wedge cracking is variable, ranging from annual to decadal, and depends primarily on winter air and ground temperatures, the latter being dependent on snow depth, vegetation, and polygon morphology (Liljedahl et al., Reference Liljedahl, Boike, Daanen, Fedorov, Frost, Grosse and Hinzman2016). Mackay (Reference Mackay1974) observed that narrow (<1 m) and wide (>2 m) ice wedges have lower cracking frequencies due to competing cracking and a reduction of thermal stress caused by the insolating effect of a greater snow cover in a wider trough, respectively. Once a crack is formed, it can be infilled by a variety of moisture sources, the most common being snow meltwater, but vapor condensation (hoarfrost) and snow have been observed as well (Lauriol et al., Reference Lauriol, Duchesne and Clark1995; French and Guglielmin, Reference French and Guglielmin2000; St-Jean et al., Reference St-Jean, Lauriol, Clark, Lacelle and Zdanowicz2011). The annual growth increment of an ice wedge can vary from <1 mm to a few millimeters, much smaller than the width of a winter crack due to thermal expansion (Mackay, Reference Mackay1974; Lewkowicz, Reference Lewkowicz1994).
With recent advances in radiocarbon dating of dissolved organic carbon (14CDOC), Lachniet et al. (Reference Lachniet, Lawson and Sloat2012) demonstrated that measurements of 14CDOC greatly reduced the potential contamination of older particulate organic carbon in the ice wedge. To date, 14CDOC chronologies of ice wedge formation are all from western Arctic Canada and they provide insights into periods when winter temperatures and snow depth were favorable to both ice wedge cracking and growth (see Grinter et al., Reference Grinter, Lacelle, Baranova, Murseli and Clark2019; Holland et al., Reference Holland, Porter, Froese, Kokelj and Buchanan2020). For example, the 14CDOC time series of ice wedges on the Blackstone Plateau in central Yukon showed that the wedges were developing when winter temperatures were cold and with sufficient snow to allow their growth and were not active during the warm early Holocene interval (Grinter et al., Reference Grinter, Lacelle, Baranova, Murseli and Clark2019). Despite the chronology of ice wedge formation in western Arctic Canada, the timing of the formation of ice wedges over the late Quaternary in other permafrost regions is largely unknown. In high Arctic Canada, ice wedge polygons occupy ~50% of the Fosheim Peninsula (Nunavut, Canada), with most situated below the Holocene marine limit in the Eureka Sound Lowlands (ESL) (Couture and Pollard, Reference Couture and Pollard1998; Bernard-Grand'Maison and Pollard, Reference Bernard-Grand'Maison and Pollard2018), Despite studies documenting the recent degradation of ice wedge polygons in the region (Becker et al., Reference Becker, Davies and Pollard2016; Ward Jones et al., Reference Ward Jones, Pollard and Amyot2020), the timing of the formation of ice wedges in this polar desert region is unknown, especially in relation to the marine regression, permafrost aggradation, and Holocene climate.
This study investigates the timing Holocene ice wedge formation in the ESL based on 32 14CDOC measurements from 15 wedges situated along an elevation gradient from near sea level to the Holocene marine limit. Considering that air temperatures in the high Arctic were always favorable to the growth of ice wedges (mean annual air temperature < −12°C; mean winter air temperature < −30°C; Lecavalier et al., Reference Lecavalier, Fisher, Milne, Vinther, Tarasov, Huybrechts, Lacelle, Main, Zheng, Bourgeois and Dyke2017; Buizert et al., Reference Buizert, Keisling, Box, He, Carlson, Sinclair and DeConto2018), we hypothesize that the formation ice wedge reflects changes in snow, which would affect the winter ground surface temperature. The objective is achieved by: (1) modeling cracking of ice wedges based on simulations of stress buildup in ice wedges following rapid cooling for a range of ground surface temperatures; (2) determining the source of moisture and DOC in the ice wedges by comparing δD-δ18O, DOC/Cl, and δ13CDOC to those in residual snow and leached active layer soils; (3) assessing the cracking frequency and period of growth of the ice wedges based on the width and 14CDOC measurements of the wedges and the Mackay (Reference Mackay1974) cracking probability model; and (4) comparing the 14C age distribution of ice wedges with regional Holocene air temperature and precipitation reconstructions and periods of peat accumulation.
Study area
The ESL is an intermontane basin on the Fosheim Peninsula situated in the west-central section of Ellesmere Island and the southeastern section of Axel Heiberg Island (Fig. 1). The region consists of flat to gently rolling terrain situated at elevations <200 m above sea level (m asl) and is underlain by the Sverdrup Basin, a folded and faulted Carboniferous to Paleogene sedimentary bedrock composed of sandstone, siltstone, and shale strata (Bell, Reference Bell1996; Hodgson and Nixon, Reference Hodgson and Nixon1998).
(color online) Location of ice wedges sampled in the Eureka Sound Lowlands, high Arctic Canada. (A–C) Maps showing sampling sites around Eureka (Gemini, Nunavut, Blacktop, and Dump slumps) and Mokka Fjord. (D and E) Field photographs of ice wedges sampled in headwall of thaw slumps at Blacktop slump (D) and Mokka Fjord (E).
The late Quaternary glacial history of the Queen Elizabeth Islands has long been debated with contrasting views. The most recent glacial history reconstruction suggests that the Innuitian Ice Sheet covered most of the Queen Elizabeth Islands and reached its maximum extent at about 15,800 14C yr BP, with the ice being 1.6 km thick in the ESL (Dyke et al., Reference Dyke, Andrews, Clark, England, Miller, Shaw and Veillette2002; England et al., Reference England, Atkinson, Bednarski, Dyke, Hodgson and Ó Cofaigh2006). Deglaciation started between 10,300–8800 14C yr BP (England et al., Reference England, Atkinson, Bednarski, Dyke, Hodgson and Ó Cofaigh2006), and remnants of the ice sheet are still found at the nearby Agassiz Ice Cap. Marine transgression followed until 8800–7800 14C yr BP, with the Holocene marine limit found between 139 and 150 m asl (Bell, Reference Bell1996). The surficial sediments on the Fosheim Peninsula reflect the regional glacial history, with weathered bedrock, till, and glaciofluvial sediments found above the marine limit and early Holocene–age marine deposits found along the lowlands below the marine limit (Hodgson and Nixon, Reference Hodgson and Nixon1998; Bell and Hodgson, Reference Bell, Hodgson, Garneau and Alt2000). Solifluction, rill washing, and slope failures reworked these deposits into colluvium deposits along hillslopes (Bell and Hodgson, Reference Bell, Hodgson, Garneau and Alt2000).
The ESL is characterized by a polar desert climate, with cold winter and cooler summer air temperatures than inland locations. The 1980–2016 mean annual air temperature recorded at the Eureka weather station, situated within 50 km of all study sites, was −18.5 ± 1.4°C; mean winter air temperature was −38.5 ± 1.7°C; mean summer air temperature was 4.3 ± 1.3°C (Environment Canada, 2019). Total precipitation reaches 77.6 ± 25.8 mm/yr, with 80% falling as snow. Late winter snow depth ranges from 5 to 35 cm, although it is likely thicker in well-developed polygonal troughs. The low precipitation compared with other sites in the high Arctic is attributed to the surrounding topography, which blocks cold ocean air masses and creates a rain shadow (Edlund and Alt, Reference Edlund and Alt1989). Vegetation is composed of patchy graminoid and prostrate dwarf-shrub (Edlund et al., Reference Edlund, Alt, Garneau, Garneau and Alt2000). Holocene air temperatures were reconstructed from δ18O measurements from the Agassiz Ice Cap. The 25-yr mean annual air temperature record shows a rapid early Holocene warming with temperatures being 6°C warmer than today at 10 ka followed by a gradual cooling to AD 1700 (Lecavalier et al., Reference Lecavalier, Fisher, Milne, Vinther, Tarasov, Huybrechts, Lacelle, Main, Zheng, Bourgeois and Dyke2017). Presently, air temperatures are rapidly warming and are now at their warmest in the past 6800–7800 yr (Lecavalier et al., Reference Lecavalier, Fisher, Milne, Vinther, Tarasov, Huybrechts, Lacelle, Main, Zheng, Bourgeois and Dyke2017; Copland et al., Reference Copland, Lacelle, Fisher, Delaney, Thomson, Main and Burgess2020). The 20-yr winter air temperature reconstruction follows a similar pattern to the mean annual one, with winter temperatures reaching about −30°C in the early Holocene, although this reconstruction does not include the correction for change in surface elevation (Buizert et al., Reference Buizert, Keisling, Box, He, Carlson, Sinclair and DeConto2018).
The climate and vegetation conditions ensure that permafrost is cold, continuous, and ~500 m thick (Pollard, Reference Pollard, Garneau and Alt2000). The active layer thickness averages 60 cm (Hodgson and Nixon, Reference Hodgson and Nixon1998; Couture and Pollard, Reference Couture and Pollard2007) and the mean annual ground temperature (measured at 15.4 m, the depth of zero annual amplitude) is −16.5°C (Pollard et al., Reference Pollard, Ward and Becker2015). Ice wedge polygons are the main periglacial landform and are mostly found within the Holocene marine limit (Couture and Pollard, Reference Couture and Pollard1998; Bernard-Grand'Maison and Pollard, Reference Bernard-Grand'Maison and Pollard2018). The ice wedges are mostly epigenetic (Pollard, Reference Pollard1991), but syngenetic wedges can be found in places of recent colluviation, with ice wedge growth estimated at 1.9–4.8 mm/yr from tritium analyses (Lewkowicz, Reference Lewkowicz1994). Based on a survey of 150 ice wedges, their average width and depth are 1.46 ± 0.56 m and 3.23 ± 0.8 m, respectively (Couture and Pollard, Reference Couture and Pollard1998).
METHODS
Modeling cracking of ice wedges
Numerical simulations were conducted using the REGO model to test for the effects of past warmer ground surface temperatures (either warmer winters or thicker snow cover) on ice wedge cracking in the ESL. REGO is a soil environmental model that uses as input a series of parameters related to the soil temperature, porosity, and permeability and soil water chemistry (Fisher, Reference Fisher2005; Fisher and Lacelle, Reference Fisher and Lacelle2014; Fisher et al., Reference Fisher, Lacelle and Pollard2020). The model solves for coupled frozen soil stresses using the thermophysical characteristics of the soils (i.e., dry porosity, density, tortuosity, thermal diffusivity), the mean ground surface temperature and its amplitude, the period of the temperature cycle, and the geothermal gradient. The theory behind temperature waves of a given period to generate various stress fields in icy soils has been explained in detail in other studies (e.g., Lachenbruch, Reference Lachenbruch1962).
The cracking depth in ice wedges developed in icy soils was modeled for initial ground surface temperatures of −38°C, −28°C, −18°C, and −8°C. REGO uses as input a cyclic temperature; therefore, to simulate a rapid decrease in temperature, the ground surface was cooled by 8°C or 17°C with a periodicity of 1/16 of a year (or 23 days, equivalent to cooler temperature below the mean lasting ca. 6 days). This periodicity represents the typical time scale for weather systems to move through a given region (ca. 10 days; e.g., Hoskins et al., Reference Hoskins, James and White1983). Simulating a faster or slower thermal change does not greatly affect the stress fields. The simulations were performed with a soil thermal diffusivity of 62 m2/yr and treat the icy soil and ice as a Maxwellian viscoelastic solid subject to elastic and thermal strains. The [1] decaying temperature wave that propagates into the icy permafrost and ice wedge, [2] the associated thermal expansion, and [3] the resulting elastic and creep strains were used to determine the nonhydrostatic stress field. The ice wedge (pure ice) will typically crack when tension exceeds ~0.4–0.6 MPa; icy soils will typically crack when tension exceeds 1.5–2 MPa (Williams and Smith, Reference Williams and Smith1989).
Field sampling of ice wedges
Epigenetic ice wedges were sampled in July 2019 from the headwall of five thaw slumps situated between the Holocene marine limit and the coastline: Gemini (120 m), Nunavut (73–100 m), Black Top (60–66 m), Mokka Fjord (55–58 m), and Dump (34–36 m) (Fig. 1). At the Nunavut, Blacktop, and Mokka Fjord slumps, the ice wedges developed in ice-rich, silty-clay marine sediments, whereas at Gemini and Dump slumps, the ice wedges developed in tabular bodies of massive ice.
Depending on accessibility, two to five ice wedges were sampled from each site. The ice wedges were sampled by first removing the overlying active layer soil and exposing the full width of the wedge. Each ice wedge was then cored from its center and near the edges using a mini-CRREL coring kit. The ice wedges were cored from their surfaces to 30 cm depths, and all samples were wrapped in aluminum foil and kept frozen in a chest freezer in Eureka. The active layer soil above each ice wedge was sampled at 5 cm intervals in Ziploc bags. Residual snow patches were also collected in well-sealed Ziploc bags, melted, and transferred to 20 ml HDPE bottles. All samples were shipped to the University of Ottawa for analyses.
Laboratory analyses
In the laboratory, the ice wedge samples were first melted in prebaked 500 mL glass beakers covered with plastic wrap and sealed to ensure no evaporation. Samples were then filtered in prerinsed 0.45 μm cellulose nitrate filters, with subsamples transferred in 20 mL plastic vials for major ions and δD-δ18O analyses, in 40 ml amber glass vials for DOC and δ13CDOC analyses, and in prebaked 1 L glass amber bottles for 14CDOC analyses. The total water-soluble ions in the active layer were determined in the laboratory from sequential extractions of dried soils (<2 mm fraction) using a 1:10 soil:water ratio (e.g., Conklin, Reference Conklin2005). The soil–water mixtures were shaken for 1 h at room temperature and subsequently centrifuged and filtered with 0.45 μm prerinsed cellulose nitrate filters and stored in 40 ml amber glass vials until analysis for major ions, DOC, and δ13CDOC.
Major cations and anions were analyzed in the Geochemical Laboratories at the University of Ottawa. The cations (Ca2+, Mg2+, Na+, and K+) were acidified with ultrapure nitric acid using an Agilent 4100 Microwave Plasma Emission Spectrometer. The anions (Cl−, SO42−, and NO3−) were analyzed unacidified by ion chromatography using a Dionex ICS-2100. Analytical precision is ±5%. To facilitate a comparison between ice wedges and leached active layer samples, the DOC results were converted to molar ratios with respect to Cl−, a conservative tracer.
The source of moisture forming the ice wedge was inferred from δD-δ18O analysis. The 18O/16O and D/H ratios of the melted ice wedges and residual snow were determined using a Los Gatos Research liquid water analyzer coupled with a CTC LC-PAL autosampler and verified for spectral interference contamination (Berman et al., Reference Berman, Gupta, Gabrielli, Garland and McDonnell2009). The results are presented using the δ-notation, where δ represents the parts per thousand differences for 18O/16O or D/H in a sample with respect to Vienna Standard Mean Ocean Water. Analytical precision for δ18O and δD was ±0.3‰ and ±1‰, respectively.
The source of DOC in the ice wedges was inferred from DOC (DOC/Cl) and δ13CDOC of melted ice wedge, residual snow, and leached active layer water. The [DOC] and δ13CDOC were measured by a wet total organic carbon (TOC) analyzer interfaced to a Thermo DeltaPlus XP isotope-ratio mass spectrometer at the Ján Veizer Stable Isotope Laboratory, University of Ottawa (St-Jean, Reference St-Jean2003). The 13C/12C results are presented using the δ-notation relative to Vienna Pee Dee Belemnite. The analytical precision is ±0.5 ppm for DOC and ±0.2‰ for δ13CTOC.
Radiocarbon analysis of the DOC of the ice wedges was performed at the A.E. Lalonde Accelerator Mass Spectrometry Laboratory, University of Ottawa. Sample preparation, including the extraction of DOC from waters and graphitization for 14C analysis, is described in Murseli et al. (Reference Murseli, Middlestead, St-Jean, Zhao, Jean, Crann, Kieser and Clark2019). Graphitized samples were analyzed on a 3MV tandem mass spectrometer (Kieser et al. Reference Kieser, Zhao, Clark, Cornett, Litherland, Klein, Mous and Alary2015). The 14C/12C ratios are expressed as fraction of Modern Carbon (F14C) and corrected for spectrometer and preparation fractionation using the 13C/12C ratio measured using accelerator mass spectrometry (Crann et al. Reference Crann, Murseli, St-Jean, Zhao, Clark and Kieser2017). Radiocarbon ages are calculated as −8033ln(F14C) and reported in 14C yr BP (BP = AD 1950) (Stuiver and Polach, Reference Stuiver and Polach1977). The measurements were then converted to calendar years (cal yr BP) using OxCal v. 4.2.4 and the IntCal13 calibration curve (Bronk Ramsey, Reference Bronk Ramsey2009; Reimer et al., Reference Reimer, Bard, Bayliss, Beck, Blackwell, Ramsey and Buck2013). In the text, radiocarbon ages with a standard error between 50 and 1000 yr are rounded to the nearest 10 yr.
RESULTS
Winter temperatures and cracking of ice wedges
At Eureka, the 1980–2016 winter air temperature was −38.5 ± 1.7°C, with late winter snow depths of 14.9 ± 6.3 cm. The freezing n-factor (ratio of winter temperature at the ground surface to that in the air) for these conditions would be about 0.85 (see Riseborough and Smith, Reference Riseborough and Smith1998), which would give a winter ground surface temperature of about −33°C. The resulting stress field following rapid cooling of a ground surface with a temperature of −33°C can generate cracks to a depth of 3.1–3.6 m, which is the depth where tension exceeds 0.4–0.5 MPa (Fig. 2). In the early Holocene, the winter air temperature was ca. 6−8°C warmer than today, and assuming a thicker snowpack (30–50 cm, corresponding to freezing n-factor of 0.65 to 0.5; Riseborough and Smith, Reference Riseborough and Smith1998), winter ground surface temperatures could reach about −18°C. For this ground surface temperature, cracking depths would be in the 2.7 to 3.2 m range. These simulated cracking depths in ice wedges appear reasonable, as they are within the range of measured ice wedge depths in the region (average of 3.23 ± 0.8 m; Couture and Pollard, Reference Couture and Pollard1998). The simulations also predict little to no cracking at winter ground surface temperatures > −10°C, as reported in the literature, because the icy permafrost is less brittle and can respond by creep and elastic strain, so tension is poorly developed (see Mackay, Reference Mackay1993; Fortier and Allard, Reference Fortier and Allard2005). In the ESL, winter ground temperatures > −10°C could be reached if snow depth exceeds 1 m (corresponding to freezing n-factors in the 0.4 to 0.3 range; Riseborough and Smith, Reference Riseborough and Smith1998).
(color online) Stress field and cracking depth in icy permafrost as a function of mean winter ground surface temperatures. (A) An example of simulated winter ground surface temperature with average of −18°C subjected to a decrease of 17°C around a period of 1/16 of a year. (B) The stress field that results from the rapid cooling of winter ground surface temperature averaging −18°C. Tensions are indicated as positive numbers, and pressures are negative numbers. The stresses are generated in the icy permafrost (thermal diffusivity of 62 m2/yr), and these are assumed to be the same in the ice wedge. The icy permafrost itself will not crack until the tension stresses reach about 1.5 MPa, but the ice wedge would crack when tension reaches 0.4–0.6 MPa. (C) Maximum depth at which cracking may occur in ice wedges for two possible tension thresholds (0.5 and 0.4 MPa) as a function of average winter ground surface temperatures of −38°C, −28°C, −18°C, and −8°C.
Isotope chemistry of ice wedges and active layer
At Dump slump (33–36 m asl), three ice wedges were sampled. The δ18O values in the six cores ranged from −35.1‰ to −29.2‰ (Table 1). The value of the regression slope between δD-δ18O was 6.5 (δD = 6.5 δ18O – 42.3; r 2 = 0.98), which is lower than the Eureka local meteoric water line (LMWL: δD = 7.4 δ18O – 9.1; IAEA/WMO, 2015; Fig. 3). The three ice wedges had similar geochemical composition. The concentration of major cations averaged 15.0 ± 5.2 mg/L, 6.2 ± 3.1 mg/L, and 42.7 ± 26.7 mg/L for Ca2+, Mg2+, and Na+, respectively; major anions averaged 36.6 ± 22.7 mg/L and 60.4 ± 39.3 mg/L for Cl− and SO42−, respectively (Table 1). The [DOC] averaged 1.0 ± 0.8 mg/L, with an average δ13CDOC of −25.7 ± 1.2‰. The DOC/Cl molar ratios in the ice wedges (0.16 ± 0.14) are similar to that of snow (1.5); the DOC/Cl in the active layer is much higher (5.55 ± 2.37) (Fig. 4). The 14CDOC ages of the ice wedges range from 28,950 to 12,130 cal yr BP, and the age in the center was older relative to the edges (Table 2). A 14CDOC age of 15,170 cal yr BP was obtained from the surrounding body of massive ice.
(color online) δD-δ18O composition of ice wedges sampled at five sites in the Eureka Sound Lowlands, high Arctic Canada. Also shown is the composition of two snow samples and the Eureka local meteoric water line (LMWL: δD = 7.4 δ18O – 9.1; IAEA/WMO, 2015).
(color online) Dissolved organic carbon (DOC)/Cl (molar ratio) and δ13CDOC composition in the ice wedges, snow, and leached active layer at four sites in the Eureka Sound Lowlands, high Arctic Canada. (A) Dump slump, (B) Blacktop slump, (C) Nunavut slump, (D) Mokka slump.
δ18O, D-excess (d), and major ion (Ca2+, Mg+, Na+, Cl−, SO42−) results for ice wedge samples at five sites in the Eureka Sound Lowlands, high Arctic Canada.a
a Distances from the centers of sampled ice wedges are provided when measurements were made; otherwise, they are indicated relative to the center (C) or edge (E) of the wedge. n/d, no data.
Dissolved organic carbon (DOC), δ13CDOC, and 14CDOC results for ice wedges sampled at five sites in Eureka Sound Lowlands, high Arctic Canada.a
a Distances from center of ice wedges are provided when measurements were made, otherwise, they are indicated relative to center (C) or edge (E) of the wedge. n/d, no data.
At Mokka Fjord slump (55–58 m asl), five ice wedges were sampled. The δ18O values in the 10 cores ranged from −31.9‰ to −28.9‰ (Table 1). The value of the regression slope between δD-δ18O was 7.5 (δD = 7.5 δ18O – 9.6; r 2 = 0.87), which is close to the Eureka LMWL (Fig. 3). All three wedges had similar geochemical composition. The concentration of major cations averaged 14.3 ± 9.6 mg/L, 5.7 ± 5.8 mg/L, and 22.8 ± 23.7 mg/L for Ca2+, Mg2+, and Na+, respectively; major anions averaged 27.6 ± 21.8 mg/L and 47.6 ± 65.8 mg/L for Cl− and SO42−, respectively (Table 1). The [DOC] averaged 8.4 ± 15.1 mg/L, with δ13CDOC averaging −25.9 ± 0.5‰. The DOC/Cl ratios in the ice wedges (0.31 ± 0.26) are similar to that of snow; the DOC/Cl in the active layer is higher (13.44 ± 9.58) (Fig. 4). The five ice wedges were analyzed for 14CDOC; the ages ranged from 8100 to 2586 cal yr BP, and the 14CDOC age was younger in the center relative to the edges (Table 2).
At Blacktop slump (60–66 m asl), three ice wedges were sampled. The δ18O values in the seven cores ranged from −31.7‰ to −27.6‰ (Table 1). The value of the regression slope between δD-δ18O was 4.1 (δD = 4.1 δ18O – 108.5; r 2 = 0.95), which is lower than the Eureka LMWL (Fig. 3). The concentration of major cations averaged 18.8 ± 8.3 mg/L, 6.6 ± 2.7 mg/L, and 69.6 ± 105.9 mg/L for Ca2+, Mg2+, and Na+, respectively; major anions averaged 102.8 ± 145.4 mg/L and 33.1 ± 19.6 mg/L for Cl− and SO42−, respectively (Table 1). The [DOC] averaged 5.5 ± 2.9 mg/L, with the δ13CDOC averaging −24.5 ± 0.1‰. The ice wedges have DOC/Cl (0.16 ± 0.13) ratios similar to that of snow; the DOC/Cl in the active layer is higher (5.76 ± 3.91) (Fig. 4). The 14CDOC of the ice wedges ranged from 33,130 to 5998 cal yr BP and the ice was younger in the center relative to the edges (Table 2).
At Nunavut slump (73–100 m asl), three ice wedges were sampled. The δ18O values in the nine cores ranged from −32.4‰ to −28.3‰ (Table 1). The value of the regression slope between δD-δ18O was 6.2 (δD = 6.2 δ18O – 50.7; r 2 = 0.91), which is lower than the Eureka LMWL (Fig. 3). All three wedges had similar geochemical composition (Fig. 4). The concentration of major cations averaged 11.9 ± 6.4 mg/L, 4.7 ± 3.2 mg/L, and 21.9 ± 14.7 mg/L for Ca2+, Mg2+, and Na+, respectively; major anions averaged 20.9 ± 14.7 mg/L and 34.4 ± 21.1 mg/L for Cl− and SO42−, respectively (Table 1). The [DOC] averaged 3.6 ± 1.5 mg/L, with the δ13CDOC averaging −25.1 ± 0.5‰. The ice wedges have DOC/Cl (0.72 ± 0.60) ratios similar to that of snow; the DOC/Cl in the active layer is higher (12.61 ± 11.90) (Fig. 4). The 14CDOC of the three ice wedges ranged from 8100 to 6885 cal yr BP, and with the exception of EU19-W5, a younger age was obtained in the center relative to the edges (Table 2).
At the Gemini site (120 m asl), two ice wedges were sampled. The δ18O was measured in one sample only and yielded a value of −28.7‰ (Fig. 3). The concentration of major cations in the two wedges averaged 12.5 ± 1.9 mg/L, 5.4 ± 0.9 mg/L, and 11.5 ± 0.6 mg/L for Ca2+, Mg2+, and Na+, respectively; major anions averaged 13.5 ± 0.9 mg/L and 5.9 ± 1.3 mg/L for Cl− and SO42−, respectively (Table 1). The [DOC] averaged 4.9 ± 1.5 mg/L. The two ice wedges have DOC/Cl (1.07 ± 0.26) ratios similar to that of snow; the active layer was not sampled at this site (Fig. 4). The 14CDOC of the two ice wedges ranged from 5412 to 2753 cal yr BP, and the younger age was obtained in the center relative to the edges (Table 2).
DISCUSSION
Source of moisture and DOC in the ice wedges
On the Fosheim Peninsula, the glacial and marine sediments contain reworked Tertiary material, including carbon, that could yield older 14C ages if it leached into DOC. The source of DOC in the ice wedges was first assessed from δD-δ18O and DOC/Cl before the 14CDOC measurements were used to infer the timing of ice wedge formation. Once a crack forms in winter, it can be infilled by one or more sources of moisture: (1) snow, (2) hoarfrost accretion, and (3) snow meltwater (Mackay, Reference Mackay1992; Lauriol et al., Reference Lauriol, Duchesne and Clark1995; St-Jean et al., Reference St-Jean, Lauriol, Clark, Lacelle and Zdanowicz2011; Boereboom et al., Reference Boereboom, Samyn, Meyer and Tison2013). The δ18O composition of the ice wedges ranges from −35.1‰ to −27.6‰, which is in the range of the two sampled residual snow patches (Fig. 3). The δD-δ18O composition of the ice wedges is distributed along a regression slope value of 6.6 (δD = 6.6 δ18O – 37.8; r 2 = 0.91), which is lower than the LWML (7.4) and suggests ice wedges developed from the freezing of snow meltwater that infiltrated the cracks (see Lacelle, Reference Lacelle2011).
For ice wedges forming from the freezing of snowmelt, the DOC in the wedges originates from the snowmelt, but additional DOC may be leached from the active layer as the meltwater infiltrates the cracks. The [DOC] of the ice wedges in the ESL averaged 3.3 ± 1.7 mg/L (0.88 to 8.7 mg/L), which is within the range of the snow sample (2.1 mg/L) and concentrations measured in ice wedges on Hershel Island (Yukon Territory, Canada) (2.4 to 8.8 mg/L; Tanski et al., Reference Tanski, Couture, Lantuit, Eulenburg and Fritz2016), in central Yukon Territory (8.1 to 16.3 mg/L; Grinter et al., Reference Grinter, Lacelle, Baranova, Murseli and Clark2019), and in Alaska (USA) and Siberia (Russia) (1.6 to 28.6 mg/L; Fritz et al., Reference Fritz, Opel, Tanski, Herzschuh, Meyer, Eulenburg and Lantuit2015). Considering that the concentration of solutes can increase during freezing, the DOC/Cl molar ratio was used to assess the source of DOC in the ice wedges. At the five sites, the DOC/Cl ratios of the ice wedges (0.36 ± 0.29) are closer to that of the snow (1.5) but are significantly lower than in the leached active layer samples (8.87 ± 8.59; two-tailed t-test, P < 0.05) (Fig. 4). The δ13CDOC cannot distinguish between the source of the DOC in the ice wedges (−25.4 ± 0.8‰), as it is within the range of both the snow (−26.8‰) and the leached active layer (−24.9 ± 0.6‰). Therefore, the DOC/Cl suggests that the DOC in the ice wedges is sourced from the snowmelt and not leached from the active layer, which can release Tertiary-age organics contained in the frozen marine sediments. Although the 14CDOC of snow meltwater was not measured, it is assumed to have a modern DOC signature at the time of infiltration. For example, Grinter et al. (Reference Grinter, Lacelle, Baranova, Murseli and Clark2019) showed that snow meltwater and spring freshets of the Ogilvie River in central Yukon had similar [DOC], both with modern 14CDOC ages.
Cracking probability, width, and period of growth of ice wedges
Epigenetic ice wedges, such as those sampled in the ESL, develop in pre-existing permafrost and grow progressively wider from the center of the wedge (Mackay, Reference Mackay1990). The ice wedges in this study yielded younger 14CDOC ages near the center with older ages near the edges, and positive relationships were observed between their widths and ranges of 14CDOC ages (Table 2). However, while one can count layers to establish a robust chronology for glaciers, the same cannot be done with individual veins forming the ice wedges. The frequency of ice wedge cracking is variable and depends primarily on winter air and ground temperatures, the latter being dependent on snow depth, vegetation cover, and polygon morphology (Mackay, Reference Mackay1992, Reference Mackay1993; Fortier and Allard, Reference Fortier and Allard2005; Abolt et al., Reference Abolt, Young, Atchley and Harp2018).
Based on measurements from Garry Island (Mackenzie Delta region, Northwest Territories, Canada), Mackay (Reference Mackay1974) suggested that the probability of cracking (P) of ice wedges follows a Gaussian distribution and can be estimated from the average width (μ) and standard deviation (σ) of a population of ice wedges:
This equation can then be combined with the annual ice wedge growth increment (ΔX, typically 1–5 mm) to estimate the number of years (Tji) required to increase the width from Xi to Xj:
where, Pji is the mean probability of cracking for that growth interval.
The time period of ice wedge growth of width Aj, can then be estimated from Eq. 3 and compared with our measurements of 14CDOC and width for the ice wedges:
A cracking probability curve was established for ice wedges in the ESL using the average widths of 1.46 ± 0.56 m (n = 150; Couture and Pollard, Reference Couture and Pollard1998), and the estimated age curve was then established using three possible growth increments (2, 2.5, and 3 mm) that are within the range of that reported by Lewkowicz (Reference Lewkowicz1994). The widths and 14CDOC age ranges of the sampled ice wedges from the different sites show a good fit with the theoretical age distribution curve as a function of ice wedge width using the 2–2.5 mm growth increment, especially considering the sensitivity of the age curve to the cracking probability (an improved fit is obtained using an average width of 1.36 m instead of 1.46 m; Fig. 5). The good fit between the widths of ice wedges and 14CDOC age range appears to validate the concept of cracking probability, width, and growth period proposed by Mackay (Reference Mackay1974) and suggests that the age of ice wedges cannot simply be determined by interpolating between measurements, as the cracking frequency is modified as the wedge grows.
(color online) Comparison of radiocarbon dating of dissolved organic carbon (14CDOC) of ice wedges sampled at different elevation with the regional isostatic curves. Marine regression curves are from Hodgson and Nixon (Reference Hodgson and Nixon1998) and Bell and Hodgson (Reference Bell, Hodgson, Garneau and Alt2000).
Paleo-geography and development of ice wedges in the ESL
The latest glacial reconstruction suggests that the ESL was glaciated by the Innuitian Ice Sheet and deglaciation started between 10,300 and 8800 14C yr BP (England et al., Reference England, Atkinson, Bednarski, Dyke, Hodgson and Ó Cofaigh2006). Marine transgression occurred until 8800–7800 14C yr BP with the Holocene marine limit near Eureka reaching 146 m asl (Bell, Reference Bell1996). The 14CDOC ages obtained from ice wedges across elevations of 120 to 33 m asl are compared with the local marine emergence curves (Fig. 6). The majority of the 14CDOC ages of the ice wedges sampled at Gemini, Nunavut, and Mokka Fjord slumps fall along or above the marine regression curves. This suggests that, at these sites, conditions became suitable for the aggradation of permafrost and subsequent ice wedge development in the marine sediments occurred immediately or shortly after permafrost aggradation. This is analogous to contemporary studies that observed the development of ice wedges in recently drained lake beds, for example, at the Illisarvik drained-lake basin in the Mackenzie Delta (Mackay, Reference Mackay2000).
Probability of ice wedge cracking and age–width relation in ice wedges in the Eureka Sound Lowlands, high Arctic Canada. Dashed line represents the probability of ice wedge cracking, which is assumed to follow a Gaussian distribution for (A) a mean of 1.46 and standard deviation of 0.56 m and (B) a mean of 1.36 and standard deviation of 0.56 m (see Eq. 1). Envelope represents the theoretical widths and ages of ice wedges for growth increments of 2 to 3 mm (see Eqs. 2 and 3). Dots represents the measured width and radiocarbon dating of dissolved organic carbon (14CDOC) age range in ice wedges in the Eureka Sound Lowlands.
The 14CDOC ages of ice wedges at Blacktop slump (28,950 to 5251 14C yr BP) and Dump slump (24,870 to 10,310 14C yr BP), both situated closest to the coast, are distributed below the marine regression curves. This was an unexpected finding, either at these sites: (1) ice wedges were developing in an ice-free region during the last glaciation (see England, Reference England1987; Lemmen, Reference Lemmen1989); (2) the 14CDOC ages were contaminated with Tertiary organic material; or (3) these ice wedges formed from a mixture of late Pleistocene–age glacial meltwater and snow meltwater. The development of ice wedges in an ice-free region can be ruled out because the wedges at Dump slump developed in massive ice with DOC dated to 12,730 14C yr BP (similar to the age obtained in the ice wedges; Table 1), whereas those at Blacktop slump developed in Holocene-age ice-rich marine sediments. A contamination from Tertiary-age organic material can also be ruled out as the DOC/Cl ratios of the ice wedges do not support the leaching of DOC from the frozen active layer. The most likely explanation is that these ice wedges formed from a mixture of late Pleistocene–age glacial meltwater and snow meltwater, with the glacial meltwater providing the majority of the source water. Support for this hypothesis comes from the lower δ18O composition of these ice wedges (−33.6‰ to −29.1‰; Fig. 3), and their [DOC] and major ions (DOC: 2.7 ± 3.0 mg/L; Ca2+: 19.3 ± 5.4 mg/L; Mg2+: 8.1 ± 2.3 mg/L; Na+: 97.2 ± 104.4 mg/L; Cl−: 108.3 ± 157.2 mg/L; SO42−: 75.2 ± 28.9 mg/L), which are slightly higher relative to the other ice wedges but in the range of late Pleistocene glacial ice at Penny Ice Cap (Fisher et al., Reference Fisher, Koerner, Bourgeois, Zielinski, Wake, Hammer and Clausen1998). The landscape following deglaciation was vastly different from today. The ice wedges at Blacktop and Dump slumps could have developed shortly after marine regression (ca. 6000 14C yr BP for the 30–70 m elevation range), with the water being sourced from the melting of residual ice caps remaining on top of Blacktop Ridge or other nearby locations (see Bell, Reference Bell1992).
Holocene ice wedge activity in the ESL
The 14CDOC age distribution of ice wedges can be used as proxies for paleoenvironmental conditions, because their development requires: (1) the presence of icy permafrost; (2) mean winter air temperatures < −10°C; (3) surface conditions that allow cracking of the ice wedge when air temperatures rapidly decrease; and (4) a moisture source to fill the crack (Washburn Reference Washburn1980; Harry and Gozdzik Reference Harry and Gozdzik1988; Christiansen et al. Reference Christiansen, Matsuoka and Watanabe2016). Figure 7 compares the 14CDOC age distribution of ice wedges in the ESL with Holocene paleotemperature reconstructions from the nearby Agassiz Ice Cap (25-yr mean annual air temperature record [Lecavalier et al., Reference Lecavalier, Fisher, Milne, Vinther, Tarasov, Huybrechts, Lacelle, Main, Zheng, Bourgeois and Dyke2017] and 20-yr mean winter air temperature record [Buizert et al., Reference Buizert, Keisling, Box, He, Carlson, Sinclair and DeConto2018]), Holocene precipitation reconstructions on Victoria Island that included modern analog sites from the Fosheim Peninsula (Peros and Gajewski, Reference Peros and Gajewski2008), and the period of peat accumulation in the Fosheim Peninsula (Garneau, Reference Garneau, Garneau and Alt2000). The 14CDOC ages, beyond having a degree of analytical and calibration uncertainties, represent a range in ages, as samples were collected across multiple ice veins using a core barrel with diameter of 8.8 cm and 14CDOC ages were grouped in 500 yr bins.
Frequency distribution of ice wedge activity in the Eureka Sound Lowlands, high Arctic Canada, compared with various paleoclimate proxies. (A) 25-yr mean annual temperature anomaly derived from Agassiz Ice Cap δ18O record (Lecavalier et al., Reference Lecavalier, Fisher, Milne, Vinther, Tarasov, Huybrechts, Lacelle, Main, Zheng, Bourgeois and Dyke2017); (B) 20-yr mean winter air temperature derived from Agassiz Ice Cap δ18O record (Buizert et al., Reference Buizert, Keisling, Box, He, Carlson, Sinclair and DeConto2018); (C) annual precipitation reconstruction from a pollen record on Victoria Island (Peros and Gajewski, Reference Peros and Gajewski2008). The two curves are derived from distinct climate transfer functions. MAT, modern analogue technique; WAPLS, weighted averaging partial least squares. Distribution of peat accumulation age on the Fosheim Peninsula (Garneau, Reference Garneau, Garneau and Alt2000); (E) radiocarbon dating of dissolved organic carbon (14CDOC) age distribution from ice wedges in the Eureka Sound Lowlands.
The distribution of the 32 14CDOC ages of ice wedges in the ESL cluster in three periods of activity: (1) immediately to shortly after marine regression and subsequent permafrost aggradation in the marine sediments (9000–6500 cal yr BP), (2) between 6000 and 4500 cal yr BP, and (3) between 4000 and 2500 cal yr BP. However, given that the ice wedges were sampled at discrete locations and reported ages are not continuous across their width, it is likely the ice wedges were active over the 9000 to 2500 cal yr BP period (Fig. 7). The 14CDOC age distribution of ice wedges is very similar to the period of peat accumulation on the Fosheim Peninsula, where all sites with peat accumulation over 50 cm thick yielded 14C ages between 8500 and 2500 cal yr BP (Garneau, Reference Garneau, Garneau and Alt2000). The temperatures were likely always favorable to cracking and the growth of ice wedges (even if snow reached 1 m depth), considering that Holocene winter temperatures in the ESL reached a maximum of about −30°C (Fig. 2). However, in this polar desert environment, both the ice wedges and peatlands require sufficient moisture to develop. Reconstructed Holocene precipitation from Victoria Island suggests that precipitation during the early and mid-Holocene was 10%–15% higher than in the present day. The increased precipitation was attributed to the marine transgression and reduced sea ice cover, which imparted a stronger maritime influence on the climate (Koerner, Reference Koerner1977). A transition to a near present-day level of precipitation occurred near the onset of the late Holocene (4200 cal yr BP) with the cooler and drier conditions lasting through the Little Ice Age (Peros and Gajewski, Reference Peros and Gajewski2008). Therefore, the transition to the drier climate conditions near the start of the late Holocene corresponds to the reduced ice wedge formation and peat accumulation, as an insufficient amount of moisture was present to allow for their development. An estimate of the volume of snow meltwater required to fill a crack in an ice wedge can be made with knowledge on crack geometry and snow depth. Approximately 4.5–5 L of snow meltwater is required to fill a 3- to 3.5-m-deep crack over a 1 m section with an ~3-mm-wide vein. For the average width of polygonal troughs in the ESL (1.46 ± 0.56 m), sufficient meltwater is provided if the snow depth is about 30 cm. This calculation assumes that the snow meltwater is sourced directly from the trough above the crack and ignores the fact that snowmelt may also originate from the center of polygons or from nearby hillslopes. However, it does suggest that snow depth <30 cm may not provide sufficient meltwater to fill a 3- to 3.5-m-deep crack. The 1980–2016 snow depths at Eureka are near that threshold, with late winter snow depth ranging between 5 and 35 cm (Environment Canada, 2019). The observed rejuvenation of certain ice wedges in the ESL (see Lewkowicz, Reference Lewkowicz1994) could be associated with the increase in precipitation over the past few years (e.g., Copland et al., Reference Copland, Lacelle, Fisher, Delaney, Thomson, Main and Burgess2020). Therefore, the relationship between current ice wedge formation and snow conditions is likely symbiotic, with the distribution of snow affecting the ground thermal regime and the formation of ice wedge troughs affecting the snow distribution and moisture sources.
CONCLUSION
This study investigated the timing of formation of ice wedges in the ESL, a high Arctic polar desert environment. Based on the results, five conclusions can be made.
First, simulations of stress induced in ice wedges following rapid cooling suggest that the cracking of ice wedges may occur when ground surface temperatures are < −10°C, with cracks to depths of 2.5 to 3.5 m. These results correspond well with the observed depths of epigenetic ice wedges in Eureka. Ice wedges observed at greater depths (>4 m) are likely of syngenetic origin and developed in colluviated sediment along hillslopes.
Second, the δ18O composition, [DOC], and DOC/Cl molar ratios of the ice wedges match those of residual snow patches. This suggests that ice wedges developed from the freezing of snow meltwater, with the DOC sourced from snowmelt and not from the leaching of organics in the active layer.
Third, the majority of ice wedges had 14CDOC ages younger than the regional isostatic curves, suggesting that conditions became suitable for the aggradation of permafrost and subsequent ice wedge development in the marine sediments after marine regression.
Fourth, Mackay (Reference Mackay1974) proposed a cracking probability model of ice wedges from which their growth periods could be estimated from their widths. The model was tested against the average widths of ice wedges in the ESL and three possible growth increments (2, 2.5, and 3 mm). A good fit was observed between the widths of ice wedges and 14CDOC age ranges, which supports the concept of cracking probability, width, and growth period proposed by Mackay (Reference Mackay1974) and suggests that the age of ice wedges cannot simply be determined by interpolating between measurements.
Finally, the distribution of 14CDOC of ice wedges in the ESL follows closely that of regional peat accumulation, with ice wedge active mainly between 9000 and 2500 cal yr BP. Given that temperatures were always favorable to ice wedge growth, the corresponding periods of ice wedge activity and peat accumulation reflect increased precipitation.
Acknowledgments
We would like to thank C. Roy and B. Faucher for valuable field assistance and the staff at the Eureka Weather Station for accommodations. We thank T. Porter, an anonymous reviewer, and the editorial office (D. Booth and T. Lowell) for their constructive comments.
Financial support
This project was supported by a Natural Sciences and Engineering Research Council of Canada (NSERC) Discovery Grant, with logistical support provided by NSERC Northern Supplement, Polar Continental Shelf Project (project no. 653-19), and the Northern Scientific Training Program.
