Discovery of a hypersaline subglacial lake complex beneath Devon Ice Cap, Canadian Arctic

See allHide authors and affiliations

Science Advances  11 Apr 2018:
Vol. 4, no. 4, eaar4353
DOI: 10.1126/sciadv.aar4353


Subglacial lakes are unique environments that, despite the extreme dark and cold conditions, have been shown to host microbial life. Many subglacial lakes have been discovered beneath the ice sheets of Antarctica and Greenland, but no spatially isolated water body has been documented as hypersaline. We use radio-echo sounding measurements to identify two subglacial lakes situated in bedrock troughs near the ice divide of Devon Ice Cap, Canadian Arctic. Modeled basal ice temperatures in the lake area are no higher than −10.5°C, suggesting that these lakes consist of hypersaline water. This implication of hypersalinity is in agreement with the surrounding geology, which indicates that the subglacial lakes are situated within an evaporite-rich sediment unit containing a bedded salt sequence, which likely act as the solute source for the brine. Our results reveal the first evidence for subglacial lakes in the Canadian Arctic and the first hypersaline subglacial lakes reported to date. We conclude that these previously unknown hypersaline subglacial lakes may represent significant and largely isolated microbial habitats, and are compelling analogs for potential ice-covered brine lakes and lenses on planetary bodies across the solar system.


Subglacial lakes in Antarctica mainly exist where temperatures at the glacier bed are maintained at the pressure-melting point from a combination of geothermal or frictional heating and the thermal insulation provided by the thick ice cover (13). In contrast, many of the Greenland subglacial lakes (46) receive surface meltwater input, providing heat that prevents the subglacial water bodies from freezing even at temperatures below the pressure-melting point (5, 6). Liquid subglacial water has also been found to exist at subfreezing temperatures beneath Taylor Glacier in the McMurdo Dry Valleys, Antarctica, where the composition of the brine-rich water depresses its freezing point (7, 8). Research over the past two decades has revealed increasing evidence that subglacial aquatic environments can host microbial life (911), and recent direct sampling of an Antarctic subglacial lake confirmed the presence of a viable microbial ecosystem in this dark and cold environment (12). Terrestrial subglacial water systems are therefore considered as potential analogs for microbial habitats on other planetary bodies, where the presence of liquid water beneath ice has been inferred (1315).

Devon Ice Cap (DIC) is one of the largest ice caps in the Canadian Arctic (Fig. 1A) and has been surveyed extensively by radio-echo sounding (RES), a well-established technique used to identify subglacial lakes from hydraulically flat, bright, and specular (mirror-like) reflections from the glacier bed (3, 16). The interior region of DIC consists of a rugged terrain with mountain ridges reaching 1700 m above sea level (asl) and several 100- to 200-m-deep bedrock troughs (Fig. 1B). This region has previously been inferred to be cold-based, where the ice is frozen to the bed (17, 18). Here, we present a combined interpretation of RES and geological data, which reveals evidence for the existence of a previously unknown hypersaline subglacial lake complex near the ice divide of DIC.

Fig. 1 Radar evidence for a subglacial lake complex on DIC.

(A) Landsat image overlain with the location of the subglacial lakes (blue) and ice divides (black). (B) Bedrock topography of the DIC summit area overlain with selected radar profiles. The subglacial lakes are located within the bedrock troughs T1 and T2. (C) Relative reflectivity along the radar transects atop of bedrock contours. (D) Radar transects A to F. The top panel shows the radargrams. The middle panel shows the ice surface and bedrock elevation (black) and the hydraulic head (blue). The bottom panel shows the relative basal reflectivity along the transect (black) and the recorded specularity content along cross-profile F (blue). The estimated lake extents are highlighted by blue shading.


Radar evidence for two subglacial lakes

A total of eight radar transects show relative bed reflectivity that is 10 to 15 dB higher over two troughs, T1 and T2, than in the surrounding area (Fig. 1). These reflectivity anomalies are 1.6 to 2.4 SDs above the mean of all bed reflectivities measured on DIC (fig. S1). The larger dielectric contrast of an ice-water interface compared to an ice-dry rock interface results in significantly stronger bed echoes, where the relative reflected power of fresh water and seawater with respect to unfrozen bedrock is expected to be ~10 and ~12 dB higher, respectively (19). We therefore infer that the observed reflectivity anomalies over T1 and T2 result from the presence of subglacial water.

The character of the ice bottom interface is further evaluated by computing the radar specularity content (16, 20, 21) along profile F over T1 (Fig. 1). The specularity content is an expression of the relative proportions of specularly reflected (mirror-like) and scattered radar energy and allows the discrimination of smooth, flat ice-water interfaces (high specularity) and comparatively rough ice-rock interfaces (low specularity) (16). High specularity anomalies have typically been used to identify areas of subglacial water and lakes (2, 21, 22), and very high values have been used to infer areas of basal melting beneath an ice shelf (22). We therefore propose that the observed specularity anomaly over T1 indicates subglacial water. The exceptionally high specularity content (>0.75; Fig. 1D) in the center of T1 may result from basal melting, further smoothing the ice-water interface. Although no specularity data are available over T2, the good agreement between the elevated specularity content and reflectivity anomalies over T1 (Fig. 1) supports the interpretation that the similar reflectivity anomalies over T2 result from the presence of subglacial water. On the basis of the extent of the radar signatures indicating subglacial water, and their confinement within two distinct topographic basins, we conclude that T1 and T2 each represent subglacial lakes. These subglacial lakes are located 2 to 5 km on either side of the west-east running ice divide, with mean ice thicknesses of 560 m over T1 and 740 m over T2 (fig. S2). We estimate subglacial lake areas of about 5 and 8.3 km2 for T1 and T2, respectively. These sizes are comparable to the majority of known (1, 3) and predicted (23) subglacial lakes in Antarctica and Greenland; however, the full extent of the water bodies beneath DIC may be underestimated due to limited radar coverage.

Subglacial hydraulic head over the lake area

Hydraulic heads calculated along the radar transects as well as from surface and bed digital elevation models (DEMs) are mostly flat across T2, with a local hydraulic minimum in the western part of the trough (Figs. 1D and 2A). A flat hydraulic head indicates hydrostatic equilibrium between the subglacial water and the overlying ice (24), indicating that the ice is afloat above the basal water. The criteria of hydraulically flat reflections are commonly used to identify subglacial lakes (3), further supporting the evidence for a subglacial lake in T2. Although the hydraulic head across T1 is relatively flat, no local hydraulic minimum is observed. Instead, the gradients of the hydraulic head suggest that water flows out of T1 in a northwest direction (Fig. 2A). If water outflows from T1, it is possible that the lake is replenished by active subglacial melting, which could explain the very high (>0.75) specularity content over T1. However, due to the relatively small horizontal lake extents and their confinement in bedrock troughs, bridging stresses in the overlying ice may prevent a fully developed hydraulic equilibrium. Other possible explanations for gradients in the hydraulic head include uncertainties due to undetected bedrock features, where no RES data are available, or radar anomalies that arise from shallow water or water-saturated sediments.

Fig. 2 Basal temperature, hydraulic head, and projected salt-bearing outcrops beneath DIC.

(A) Bedrock elevation overlain with the hydraulic head (blue) and ice divides (black). The blue star indicates the location of the hydraulic head minimum in T2, whereas the arrow indicates the direction of water potentially outflowing T1. (B) Bedrock elevation contours overlain with modeled basal temperatures. The black triangle is the location of the ice temperature profile, where a basal temperature of −18.5°C was measured (26). (C) The area where the Bay Fiord Formation containing the bedded salt sequence is projected to outcrop beneath the ice is marked in yellow. Radar profiles are marked in blue/red. (D) Cross-section along radar transect CC′, revealing the projected geology consisting of the Bay Fiord Formation (yellow) and underlying sedimentary rocks (brown).

Modeled basal ice temperatures indicating hypersaline water

Using a one-dimensional (1D) steady-state advection-diffusion model (25), we calculate basal temperatures of around −15.5° ± 3.5°C in T1 and −14.3° ± 3.75°C in T2 (Fig. 2B). These modeled temperatures are in good agreement with a nearby measured ice temperature profile in which the basal temperature was −18.5°C at a depth of 300 m (26), but are well below the pressure-melting point for fresh water derived from the overlying ice (about −0.5°C). Surficial water (and heat) input that would prevent the subglacial lakes from freezing is improbable because only limited surface meltwater is produced at elevations near the DIC summit (27). We therefore conclude that the subglacial lakes beneath DIC must consist of hypersaline water with a significantly depressed freezing point temperature.

Underlying geology and source of salinity

Geological mapping shows outcrops of the Archean shield (igneous rocks) surrounding the eastern part of DIC and Cambrian-Ordovician sediments, including an evaporite unit within the Bay Fiord Formation, to the west of DIC (fig. S3) (28). Drill cores from Bathurst Island, located west of Devon Island, show that the bottom of the Bay Fiord Formation contains a bedded salt sequence, with the salt consisting almost entirely of halite (98%) (29). We use published outcrop data from Devon Island (28) to construct a 3D geological model (fig. S3), which allows us to constrain the bedrock geology beneath the ice. The modeled geology is consistent within the uncertainty for depth to magnetic basement (DMB) solutions derived from airborne magnetics data (fig. S4). Our geological model shows that the salt-bearing Bay Fiord Formation outcrops near the DIC summit area, where it encloses the subglacial lake in T1, and outcrops about 100 m above the bottom of T2 (Fig. 2). However, given the estimated uncertainty of 200 m in the projection of the geological model, it is possible that the subglacial lake in T2 is also situated within the Bay Fiord Formation. The presence of the salt-bearing unit underlying the DIC summit area supports the existence of a substantial hypersaline subglacial lake complex, where the halite contained in the Bay Fiord Formation likely provides the main source of salinity. The highest modeled ice temperatures within the uncertainty ranges for T1 and T2 are −12° and −10.5°C, respectively. Empirical studies show that to depress the freezing point of distilled water by 10.5° to 12°C with the addition of NaCl, a salinity of 140 to 160 practical salinity units (psu) (2.4 to 2.7 M NaCl) is required (30). These temperature and salinity values are similar to those of the brine-rich water body found beneath Taylor Glacier (−7.8°C, 125 psu) (7) and the ice-covered Lake Vida (−13°C, 200 psu) (31), both located in the McMurdo Dry Valleys, Antarctica.

The observation of radar signatures indicative of subglacial water in an area, where modeled basal ice temperatures are well below the pressure-melting point, and the likely outcrop of a bedded salt sequence beneath the ice in the same area provide compelling evidence for the existence of a hypersaline subglacial lake complex beneath DIC. These subglacial lakes are likely confined to bedrock troughs where the basal temperatures are higher due to a thicker, more insulating ice cover (fig. S2). The geological model indicates that other troughs near the DIC summit are also underlain by the Bay Fiord Formation (Fig. 2C); however, more radar data are required to investigate their basal hydrological conditions. In addition, outcrops of the Bay Fiord Formation are common around Ellesmere Island (28) and may be found beneath other Canadian Arctic ice caps.


Although numerous subglacial lakes have been discovered beneath the large ice sheets, the subglacial lake complex identified beneath DIC is unprecedented. These subglacial lakes and their surrounding environments are very different from those reported in Antarctica and Greenland: The DIC lakes are situated within bedrock troughs in mountainous terrain, exist at temperatures well below the pressure-melting point, do not receive surface meltwater input, and likely consist of hypersaline water derived from dissolution of a surrounding salt-bearing geological formation. The origins of the hypersaline subglacial lakes, the processes by which they formed, and their potential interactions with shallow or deep groundwater flow remain unclear. However, possible mechanisms for the formation of the hypersaline subglacial water include melting of basal ice upon direct contact with the salt in the surrounding rocks, and modulated by double-diffusion processes within the water bodies. The subglacial lakes might also be relict water bodies that formed subaerially during an interglacial period, with the water becoming increasingly briny through interactions with underlying saline rocks, and as a result of cryoconcentration after they became permanently ice covered.

The only known subglacial fluid with temperature and salinity values comparable to this hypothesized hypersaline lake system beneath DIC is the brine beneath Taylor Glacier, Antarctica (7). However, this brine body is not constrained as a spatially isolated subglacial lake, but is connected to a salt-rich groundwater system that is sourced by ancient marine water (8, 32, 33). Brine outflows from Taylor Glacier have been shown to contain active microbial communities (9), revealing that life is possible in such hypersaline, sub-ice aquatic environments. We therefore suggest that the DIC subglacial lakes may have the potential to support microbial life and could harbor a unique ecosystem. If life exists in these lakes, it could have evolved in isolation since the area was last overridden by glacier ice, which was at least 120,000 years ago (34). We conclude that the DIC’s hypersaline subglacial lakes represent good analogs for the brine bodies inferred to exist beneath and potentially within Europa’s ice shell (13) or the Martian polar ice caps (14), and are therefore compelling targets for future exploration.


Radar data and interpretation

RES data were acquired with the University of Kansas Multichannel Coherent Radar Depth Sounder (MCoRDS) (35) during spring in the years 2011, 2012, 2014, and 2015, and with the High Capability Airborne Radar Sounder (HiCARS) (36) operated by the University of Texas Institute for Geophysics during spring 2014. Ice thicknesses and basal reflectivity were extracted using a common semiautomatic method with rough localization from manual picking. Basal reflectivities were extracted from the HiCARS low gain data and from the MCoRDS CSARP_standard combined-gain data product. Basal reflectivity values were corrected for geometrical spreading losses, englacial attenuation, and variations in radar systems (37). Englacial attenuation rates were calculated and corrections were applied along each individual radar profile via an adaptive fitting method (37). We set the minimal requirement parameters to ensure a meaningful fitting to the same values as those by Schroeder et al. (37). If these requirements were not met along a profile, we applied the overall mean englacial attenuation rate derived from all radar data on DIC (26.8 dB/km, with an SD of 7.3 dB/km) to the profile. To correct for regional changes in attenuation rates and variations in radar system parameters, a long-wavelength signal (moving average with a window length of 30 km) from the bed echo along each radar profile was removed. The specularity content from the HiCARS data was extracted by comparing the bed echo response from two different focusing aperture lengths (16, 21).


An existing bedrock DEM (38) was updated using the original bedrock depth data set [derived from radar data collected over DIC in spring 2000 (38)] in conjunction with bedrock depths from MCoRDS and HiCARS data collected between 2011 and 2015. The radar-derived bedrock elevations were interpolated via a triangular linear interpolation algorithm over a 1-km grid mesh. A crossover analysis yielded a mean error of ice thickness observation of 14 m (before interpolation).

Hydraulic head

Assuming that the water pressure is equal to the ice overburden pressure, the basal hydraulic head (h) was calculated from the gridded bed elevations (B), and surface elevations (S) obtained from the ArcticDEM (Polar Geospatial Center from DigitalGlobe Inc. imagery) viaEmbedded Imagewhere ρi is the ice density (917 kg/m3) and ρb is the brine density (1150 kg/m3, corresponding to 15 weight % NaCl) (24).

Ice temperature modeling

To model basal ice temperatures, we used a 1D steady-state advection-diffusion model (25) with the following parameters: a long-term accumulation rate of 0.19 ± 0.05 m water equivalent per year (39, 40), an estimated geothermal heat flux of 65 ± 5 mW/m2 (41), and an average annual surface temperature of −23° ± 2°C at an elevation of 1825 m (42) with a lapse rate of 4.1°C/km. The lapse rate was derived from surface air temperature records along a transect on DIC, sampling air temperatures at numerous elevations (43).

Geological model

A geological map from the Canadian Arctic (28) was combined with surface elevations derived from the ArcticDEM to generate a 3D geological model (fig. S3): First, surfaces through outcropping geological formation boundaries were interpolated, and formation thicknesses were derived, assuming that all geological units lie parallel to each other. Finally, the geological units were stacked atop of the interpolated formations that are most adjacent to DIC, using the derived formation thicknesses and the assumption that the units lie parallel to each other. Combining all formation thickness uncertainties, we estimated a total uncertainty of 200 m for the projected elevation of the Bay Fiord Formation (salt-bearing unit) beneath DIC.

Depth to magnetic basement

We used DMB solutions to infer the thickness of sediment beneath DIC. To do so, a 2D Werner deconvolution (44) was applied to available scalar airborne magnetics data with a single pass targeting shallow sources (0.5 to 5 km). Uncertainties in DMB estimates are commonly between 20 and 40% of the distance between the source and the sensor (44). The results indicate sediment thicknesses consistent with our geological model within the range of uncertainties in the model and the DMB solutions (fig. S4).


Supplementary material for this article is available at

fig. S1. Histogram of all recorded basal reflectivity values on DIC.

fig. S2. Interpolated ice thickness near the DIC summit area.

fig. S3. 3D geology model.

fig. S4. Comparison of modeled geology and DMB solutions.

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is not for commercial advantage and provided the original work is properly cited.


Acknowledgments: We thank the NASA Operation IceBridge project team and the National Snow and Ice Data Center and acknowledge the use of MCoRDS airborne radar data from the Center for Remote Sensing of Ice Sheets (CReSIS) generated with support from the University of Kansas NSF grant ANT-0424589 and NASA Operation IceBridge grant NNX16AH54G. We thank the Polar Geospatial Center who provided the DEM under NSF Office of Polar Programs awards 1043681, 1559691, and 1542736. We thank F. Habbal, G. Ng, S. Palmer, and T. Richter for support during fieldwork; S. Grasby for providing us with an estimate for the geothermal heat flux; S. Kempf for assistance with radar data processing; and T. Benham for providing bedrock elevations from the spring 2000 radar data set. We also thank N. Wolfenbarger, N. Ross, an anonymous reviewer, and the editors M. Kelly and P. Bierman for providing feedback on the manuscript. This is University of Texas Institute for Geophysics (UTIG) contribution 3255. Funding: This study was funded by grants from UK Natural Environment Research Council (NE/K004999), NASA (13-ICEE13-00018), Natural Sciences and Engineering Research Council of Canada (Discovery Grant/Northern Research Supplement), Alberta Innovates Technology Futures, and the CRYSYS Program (Environment Canada). M.L.S. was partially supported by NSF 1543537 and NASA NNX16AJ64G. D.D.B. and J.S.G. were partially supported by the G. Unger Vetlesen Foundation. D.D.B. (assisted by A.R.) was supported as a Fulbright Scholar at the University of Gothenburg, Sweden. Author contributions: A.R. participated in the fieldwork, performed the data analysis and interpretation, and wrote the manuscript. D.D.B. and M.S. assisted with data analysis and interpretation. M.L.S. provided guidance for the salinity estimates and the geological interpretation. J.S.G. performed the magnetic depth to basement calculations and assisted with the reflectivity and specularity interpretation. C.G. provided guidance for the radar reflectivity and specularity interpretation. D.M.S. provided guidance for the radar attenuation correction and interpretation of radar reflectivity and specularity. D.A.Y. performed the radar data processing and derived the specularity content. D.D.B., D.A.Y., J.S.G., C.G., D.M.S., A.R., M.S., and J.A.D. designed the radar surveys. D.A.Y., C.G., A.R., and J.A.D. collected the data. All authors read and contributed comments for the preparation of the final manuscript. Competing interests: The authors declare that they have no competing interests. Data and materials availability: MCoRDS radar data were accessed through the CReSIS. UTIG radar data and the bedrock DEM will be made available upon request. ArcticDEM surface elevation data were accessed through the Polar Geospatial Center. All other data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. Additional data related to this paper may be requested from the authors.

Stay Connected to Science Advances

Navigate This Article