Open access

Identifying groundwater discharge zones in the Central Mackenzie Valley using remotely sensed optical and thermal imagery

Publication: Canadian Journal of Earth Sciences
23 June 2020


Landsat 4–5 Thematic Mapper, Landsat 8 Operational Land Imager, and RapidEye-3 data sets were used to identify potential groundwater discharge zones, via icings, in the Central Mackenzie Valley (CMV) of the Northwest Territories. Given that this area is undergoing active shale oil exploration and climatic changes, identification of groundwater discharge zones is of great importance both for pinpointing potential contaminant transport pathways and for characterizing the hydrologic system. Following the work of Morse and Wolfe (2015), a series of image algorithms were applied to imagery for the entire CMV and for the Bogg Creek watershed (a sub watershed of the CMV) for selected years between 2004 and 2017. Icings were statistically examined for all of the selected years to determine whether a significant difference in their spatial occurrence existed. It was concluded that there was a significant difference in the spatial distribution of icings from year to year (α = 0.05), but that there were several places where icings were recurring. During the summer of 2018, these recurrent icings, which are expected to be spring sourced, were verified using a thermal camera aboard a helicopter, as well as in situ measurements of hydraulic gradient, groundwater geochemistry, and electroconductivity. Strong agreement was found between the mapped icings and summer field data, making them ideal field monitoring locations. Furthermore, identifying these discharge points remotely is expected to have drastically reduced the field efforts that would have been required to find them in situ. This work demonstrates the value of remote sensing methods for hydrogeological applications, particularly in remote northern locations.


Des ensembles de données obtenues par capteur Landsat 4–5 Thematic Mapper, imageur Landsat 8 Operational Land Imager et satellite RapidEye 3 ont été utilisés pour cerner d’éventuelles zones d’émergence d’eau souterraine, sous forme d’aufeis, dans la vallée centrale du fleuve Mackenzie (VCFM) des Territoires du Nord-Ouest. Comme cette région est fait l’objet de travaux d’exploration pour le gaz de schiste et subi des changements climatiques, la localisation de zones d’émergence d’eau souterraine revêt une grande importance tant pour l’établissement de possibles voies de transport de contaminants que pour la caractérisation du réseau hydrologique. Comme suite aux travaux de Morse et Wolfe (2015), nous avons appliqué une série d’algorithmes de traitement d’images à l’imagerie pour toute la VCFM et pour le bassin versant du ruisseau Bogg (un sous-bassin de la VCFM) pour des années choisies entre 2004 et 2017. Les aufeis ont fait l’objet d’un examen statistique pour toutes les années sélectionnées afin de vérifier la présence d’une variation significative de leur répartition spatiale. Nous concluons qu’il y avait une variation significative de la répartition spatiale des aufeis d’une année à l’autre (α = 0,05), mais que des aufeis étaient récurrents en plusieurs endroits. Durant l’été 2018, la présence de ces aufeis récurrents, qui sont présumés être alimentés par des sources, a été validée à l’aide d’une caméra thermique héliportée, ainsi que par des mesures in situ du gradient hydraulique, de la géochimie de l’eau souterraine et de la conductivité électrique. Une bonne concordance a été relevée entre les aufeis cartographiés et les données de terrain estivales, ce qui en fait des sites idéaux pour la surveillance sur place. La détermination par télédétection de ces zones d’émergence aura en outre réduit considérablement les efforts de terrain qui auraient été requis pour les trouver sur place. Ces travaux démontrent l’utilité de méthodes de télédétection pour des applications hydrogéologiques, particulièrement en régions nordiques éloignées. [Traduit par la Rédaction]


The Central Mackenzie Valley (CMV) of the Northwest Territories (NWT), Canada, is a region of active shale oil and gas exploration. In advance of shale oil and (or) gas production, it is of vital importance that baseline hydrologic monitoring be established for the region so that any negative environmental effects induced by the process of oil and gas extraction can be quantified. The work presented in this article forms part of a long-term hydrologic monitoring project that aims to gain a better understanding of the interaction between groundwater and surface water resources in the CMV. In addition to being a region of shale oil and gas exploration, the CMV is undergoing rapid climate warming (Bush and Lemmen 2019). Therefore, hydrologic monitoring data, gathered through time, will provide valuable insight into the effects of climate warming that include increased fall and winter precipitation and thawing degree days (Environment Canada 2018), and permafrost thaw (Rudolph et al. 2016).
The specific objective of the work described here is to use remote sensing techniques to locate regions of expected groundwater discharge that will be utilized as priority field monitoring locations, and further, to determine whether those discharge zones are more likely to be permanent or intermittent. It was hypothesized, firstly, that icings (winter groundwater discharge features, described in detail in the next section) could be detected remotely and used as a proxy for groundwater discharge; secondly, that annually recurring icings would be more likely to represent permanent discharge features (such as springs) than would icings that occur only in 1 year; and thirdly, that in late winter, permanent icings would exhibit a warmer thermal signature than intermittent icings due to continuous water discharge. These hypotheses were evaluated in the current study by constructing icing and thermal anomaly distributions using Landsat and RapidEye-3 satellite imagery for selected years between 2004 and 2017, and by comparing the results of this remote analysis with in situ groundwater measurements for field verification. This study was performed in the greater CMV and the Bogg Creek Watershed (a smaller catchment of the CMV), as shown in Fig. 1.
Fig. 1.
Fig. 1. Central Mackenzie Valley scene used for analysis, shown as a Landsat 8 Operational Land Imager natural colour composite from 2016 with the Bogg Creek watershed as a subset, shown as a RapidEye-3 natural colour composite from 2017. Map generated using ArcGIS 10. Coordinates provided with reference to outset map. [Colour online.]
Groundwater discharge locations may represent potential pathways for contaminant transport associated with hydraulic fracturing. Although the precise relationship between fracturing operations and groundwater contamination is still contested, the U.S. Environmental Protection Agency (EPA) notes that there is greater potential for water resources nearby fracturing operations to be affected (EPA 2016). Given the remoteness of the study area and the logistical difficulty in accessing it, identifying groundwater discharge zones remotely drastically reduces the cost and time that would be required to identify them in situ. Following the work of Morse and Wolfe (2015), the methodology presented here uses remotely sensed optical and thermal imagery from Landsat 4–5, Landsat 8, and RapidEye-3 satellites to identify expected groundwater discharge zones via icings (also called aufeis).
Icings are sheet-like masses of ice that can form on pre-existing river ice or on the land surface. Three main types of icings are identified in existing literature: ground type, spring type, and river type (Carey 1973; Yoshikawa et al. 2007). River icings, which form when water discharges through a body of river ice and laps onto a frozen river surface, are not considered in this study. Only land-fast (spring and ground type) icings are considered. These two types of icings, conceptualized in Figs. 2 and 3, are groundwater sourced and can be differentiated by whether they occur from a permanent spring or from a temporary seep (Carey 1973; Pollard and van Everdingen 1992). Spring icings are formed from groundwater springs, where water tends to be sourced from sub-permafrost, deep groundwater reserves; these are often annually recurring features (Carey 1973; Yoshikawa et al. 2007). Ground icings, which are temporary or intermittent features, are formed when water in the active layer becomes “trapped” between a downward propagating freezing front and the top of the permafrost table (Carey 1973). This is expected to occur when differences in surface topography result in an uneven freezing front, when soils of differing frost susceptibility interact, or when supra-permafrost taliks are formed. Shallow groundwater flow in permafrost regions often occurs through supra-permafrost taliks that exist in higher permeability, lower ice content soils (Woo and Xia 1995; Jepsen et al. 2015). Although the physical mechanics of ground icings have not, to the knowledge of this work, been examined in any detail, it is understood that supra-permafrost taliks (distinct from lake taliks) transmit groundwater to surface water bodies (Jepsen et al. 2015; You et al. 2017) and it is therefore expected that such taliks could also transmit water directly to the ground surface, forming ground icings given a strong enough hydraulic gradient. Connon et al. (2018) report that supra-permafrost taliks thicken in response to a subsiding permafrost table and a thinning active layer; thus, it is also expected that continued climate warming will result in more supra-permafrost groundwater flow and greater availability of water for discharge.
Fig. 2.
Fig. 2. Conceptualization of a spring icing, which is sourced from sub-permafrost groundwater (GW) flowing under confined Artesian conditions. (Artesian wells included as a conceptualization for the purpose of demonstrating Artesian conditions.) A through talik remains open throughout the year; therefore, the spring icing has a continuous source of water to supply it. The downward progression of the freezing front may increase pressure head in the Artesian well later in the winter. [Colour online.]
Fig. 3.
Fig. 3. Conceptualization of a ground icing where Sw is the saturation of water and Sθ is the residual saturation of water (maximum saturation of ice). In early winter, an upward vertical hydraulic gradient is imposed in sandy sediment (low frost susceptibility) by the progressing freezing front and relatively less permeable sediment (e.g., silt) surrounding it. The freezing front has progressed further in the silt due to its higher frost susceptibility. By mid-winter, the icing has reached peak formation and the freezing front has progressed to such a depth that vertical flow no longer exists. GW, groundwater. [Colour online.]
Remote geophysical instruments have been effective at monitoring and mapping icing distribution in regions where field instrumentation is difficult. The vast majority of satellite remote sensing research focused on icing distribution in North America has examined river icings in the Brooks Range of Alaska and northern Yukon. Li et al. (1997) used single-look synthetic aperture radar from the ERS-1 satellite to monitor the seasonal growth of icings in the Inishak and Echooka river valleys. That study concluded that icings could be separated from other land cover types by observing low coherence values, noisy phase patterns, and large changes in radar backscatter. Yoshikawa et al. (2007) used a combination of near-infrared (NIR), shortwave infrared (SWIR), and thermal infrared data from Landsat and aircraft flights, and synthetic aperture radar data obtained from RADARSAT and ERS-1/2, to examine icing dynamics in another region of the Brooks Range. One of the most extensive remote sensing studies pertaining not specifically to river type icings was carried out by Morse and Wolfe (2015) with the goal of determining icing recurrence over a 24-year period. That work utilized an entire Landsat scene in the Great Slave Lake Region of the NWT — collected for 24 consecutive years — and analyzed the annual recurrence and sizes of icings.

Study area


The CMV is bound by the Mackenzie River in the east and the Mackenzie Mountains in the western portion of the study area. The valley is characterized by undulating topography that slopes gently upwards from the base of the Mackenzie River to the Mackenzie Mountains. The CMV is composed of three distinct ecosystems that, although similar overall, do possess some unique topographic characteristics. These regions are the North Mackenzie Plains, Carcajou Plains, and Mackenzie Foothills regions (AMEC 2014). The North Mackenzie Plains comprise the majority of the study area, which is the undulating valley parallel to the Mackenzie River. The Carcajou Plains represent a transition zone between the flat-wetland-dominated plains in the east and the higher elevated regions in the west (AMEC 2014). Finally, the Mackenzie Foothills represent the region parallel to the Mackenzie Mountains, which has the greatest variability in elevation throughout the entire study area, spanning from 200 m a.s.l. at the eastern side to 900 m a.s.l at the western side (AMEC 2014).


The works of Rouse et al. (1997) characterize the CMV as a subarctic climate, owing to its below freezing mean annual air temperatures. Mean annual air temperature, rainfall, and snowfall are given by Environment Canada (2018) for 1943–2012 at climate stations in Norman Wells and Tulita: For Norman Wells and Tulita, respectively, the mean air temperature is −5.7 °C and −5.3 °C, the mean rainfall is 181.0 and 183.5 mm, and the mean snowfall is 151.8 and 114.9 mm. For the Mackenzie Valley, the 2005 Arctic Climate Impact Assessment reports annual mean temperature increases of 2–3 °C during the past 50 years, and up to 4 °C during the winter for a period from 1958 to 2012 (NWT Government 2012). The same study reports an increase in precipitation from 1958 to 2012 in the fall and winter and a decrease in precipitation during the spring and summer. These trends are slight in all seasons but summer, when the decrease in mean annual precipitation is significant (∼25 mm over the course of the 50-year measured period). Zhang et al. (2019) document similar findings: a mean annual temperature increase of 2.3 °C in northern Canada occurred between 1948 and 2016 with winter having the highest mean temperature increase (4.3 °C); the authors expect that more than half of this temperature increase can be attributed to anthropogenic sources. Zhang et al. (2019) projected mean annual temperature increases in northern Canada (under worst-case emission scenarios) to be 1.8 and 5.7 °C for years 2031–2050 and 2081–2100, respectively. The same study also predicts that precipitation will occur in increasing proportions as rainfall, owing to warmer temperatures in the winter causing inhibition of snowfall. Such climatic changes may have a significant impact on the groundwater flow regime and occurrence of icings.

Permafrost and hydrogeology

Permafrost in the CMV is a combination of extensively discontinuous and intermediate discontinuous (AMEC 2014). It marks the transition zone between these two permafrost classifications. In the Norman Wells region, some areas are not underlain by permafrost at all, hence the classification of “discontinuous”, while some areas have permafrost thicknesses reportedly varying from 50 to 143 m. Permafrost is generally within 0.5–1 m of the surface, but this depth quickly increases when approaching surface water bodies. Measurements of permafrost depth and subsurface temperature monitoring indicates that taliks reside beneath many water bodies. These taliks may extend for several metres under small water bodies while permafrost may be entirely nonexistent beneath larger bodies such as the Mackenzie River or large lakes (Burgess and Smith 2000). Recent work suggests that permafrost in both the northern and southern Mackenzie Valley is thawing. Derksen et al. (2018) report that permafrost temperature in the Norman Wells area has been increasing in temperature since the mid-1980s. Active layer thickening also lends evidence to suggest that permafrost in this region is degrading; a thaw tube network throughout the CMV indicated that the active layer has been steadily thickening since 2008 (Derksen et al. 2018).
Several regional aquifers have been identified in a report by AMEC (2014), but the groundwater flow system and recharge/discharge dynamics are not well understood. It is expected that regions that have lower ice content permafrost (coarse-grained materials) experience more precipitation-sourced recharge than do regions that have higher ice content permafrost. The study area resides on the Northern Foreland Belt, which consists of the Mackenzie Plain structural domain, sandwiched between the Mackenzie Mountains to the southwest and Franklin Mountains to the northeast (Hayes and Dunn 2012). Bedrock geology consists of Cambrian to Cretaceous aged sedimentary rocks, which has undergone folding and minor faulting in a general northwest–southeast trend. Glaciation during the Quaternary led to deposition of discontinuous deposits of till, glaciolacustrine, and glaciofluvial sediments that vary in thickness.
Major aquifers in the region include the Early Cretaceous Martin House Formation, a thin, coarse-grained sandstone and conglomerate unit and the Late Cretaceous fine-grained sandstone of the Little Bear Formation (Hayes and Dunn 2012). Other potential aquifers are the sandstone and conglomerates of the Late Cretaceous Summit Creek Formation and localized surficial deposits of coarser-grained Quaternary sediments.
AMEC (2014) reports that groundwater isotopic signatures are highly variable: The shallow Summit Creek aquifer contains relatively modern groundwater, whereas the deeper Little Bear Formation contains groundwater carbon dated to around 20 000 years (AMEC 2014). Further details pertaining to the hydrogeology and geology of this region can be found in the literature (AMEC 2014; Rudolph et al. 2016).
Of particular interest for this work is the Bogg Creek Watershed, a region of proposed shale oil development and a focus area for hydrologic monitoring within the CMV. The geology of the Bogg Creek Watershed consists mostly of shales and sandstones mantled by 2–40 m of glacial sediment. Structural features in the watershed include a regional syncline to the southwest and minor anticline to the north (Fallas and MacNaughton 2013). The catchment itself resides east of the center and partially on the eastern limb of this syncline. Oil shales of the Canol Formation (located 1650 m below ground surface) were the target of oil and gas exploration activities prior to 2014 by Husky Energy (Rudolph et al. 2016). The oldest subcropping formation in the watershed is the Devonian-aged Imperial Formation, consisting mostly of shale and mudstone. Overlying this are the early Cretaceous sandstone and conglomerate of the Martin House and shale of the Arctic Red Formations. The thickest unit, the Slater River Formation, overlies the Arctic Red Formation. This is characterized by shale and mudstone with minor sandstone. The Little Bear Formation sandstone and mudstone is the youngest rock unit in the watershed, and forms the center of the syncline (Fallas and MacNaughton 2013). Till and glaciolacustrine sediments blanket the bedrock through most of the watershed, a remnant of the last glaciation. Deposits of modern fluvial and glaciofluvial sediments are also found closer to streams, while peatlands and organic deposits are common throughout the region (Côté et al. 2013). Both the bedrock and surface geology, as mapped by the Canadian Geological Survey, are illustrated in Figs. 4 and 5. The groundwater flow system within the Bogg Creek Watershed is expected to consist of regional to local flow regimes controlled largely by permafrost distribution. Regional flow systems are hypothesized to originate along the southwestern limb of the geologic syncline and in the foothills of the Mackenzie Mountains. Groundwater is hypothesized to flow through taliks, and discharge as springs and into waterbodies in the watershed. Regional flow may also originate in the northeast from the Franklin Mountains, which form the eastern limb of the regional syncline. Local flow is hypothesized to occur in supra-permafrost zones and shallow taliks, through mineral and organic soils. This system is likely seasonally dependent.
Fig. 4.
Fig. 4. Bedrock geology formations and structural features in the Central Mackenzie Valley and Bogg Creek watershed (shown as a red outline). The composition of geologic formations as they relate to groundwater discharge are discussed in the Introduction. Geologic information obtained from Fallas et al. (2013). Map generated using ArcGIS 10. [Colour online.]
Fig. 5.
Fig. 5. Primary surficial geology in the Central Mackenzie Valley and Bogg Creek watershed (show as a yellow outline). Geologic information obtained from Côté et al. (2013). Map generated using ArcGIS 10. [Colour online.]


Remote sensing methods

In this study, the identification of icings was first performed for the entire CMV using 30 m Landsat 4–5 Thematic Mapper (TM) and Landsat 8 Operational Land Imager (OLI) optical imagery, and Landsat 4–5 and Landsat 8 120 m thermal imagery for years 2004, 2009, and 2016. The acquisition dates of the imagery utilized in this work are summarized in Table 1. These three years were selected as imagery was available for the desired season and had few atmospheric interferences (i.e., low presence of clouds and haze). The study scene and geographic situation is shown in Fig. 1. Following the algorithm process developed by Morse and Wolfe (2015), late spring imagery from each of these years was used to discriminate areal icing coverage as shown in Fig. 6. Spring imagery was used as icings are most easily viewed in late spring when the majority of the snow pack has deteriorated. To obtain the icings, the Normalized Difference Snow Index (NDSI) was first applied to the spring image:
Table 1.
Table 1. Image acquisition dates for all images used in remote analysis.

Note: TM, Thematic Mapper; OLI, Operational Land Imager.

Winter thermal anomalies were not calculated for 2017.
As given for Norman Wells Climate Station A by Environment Canada (2018) on the acquisition date.
Fig. 6.
Fig. 6. Detailed process summary used to extract icings and thermal anomalies. MDSII, Maximum Difference Snow and Ice Index; NDSI, Normalized Difference Snow Index; NDVI, Normalized Difference Vegetation Indices; NDWI, Normalized Difference Water Index; TOA, top of atmosphere. [Colour online.]
Using a threshold value of 0.4, ice, snow, water, and marl water were separated from the rest of the land cover types. This algorithm is founded on the principle that snow and ice will absorb most of the incoming short-wave radiation and reflect in the visible portion of the electromagnetic spectrum. Multiple studies have confirmed the effectiveness of the NDSI for this purpose (e.g., Hall et al. 1995; Salomonson and Appel 2004; Cao and Liu, 2006). All values exceeding 0.4 were isolated from the image and used to mask the original multi-spectral top of atmosphere conversion. Then, ice was separated from snow, water, and marl water using the Maximum Difference Snow and Ice Index (MDSII):
In a similar manner to the NDSI, threshold values of the MDSII can be used to discriminate ice from other land cover types. For a study area near Great Slave Lake, NWT, Morse and Wolfe (2015) determined a value of 0.144 to be used for discrimination of ice when some snow is still present in the scene, and a value of 0.040 to be used when no snow is present. These threshold values were used for the current study as well. All late-spring scenes used to extract icings were snow free; therefore, the MDSII threshold value separated ice from water and marl water (ice–water mixtures). Values in excess of the threshold were considered to be ice and were extracted from the image.
The MDSII extracts all ice in the study area, including ice on water bodies that may not have been completely melted. To remove this ice from the result, a water mask was generated using an image from the summer of that year (shown in Table 1). As the extent and location of thaw ponds and thermokarst lakes are expected to be variable and vulnerable to modification from climate change and supra-permafrost groundwater flow (You et al. 2017; Fraser et al. 2018), separate water masks were generated for each year (rather than just one), to most closely approximate the position of surface water bodies at that time. Use of the MDSII to detect marl water also assists with the removal of non-permanent water features that may not have been removed by the water mask. To generate the water mask, the Normalized Difference Water Index (NDWI) algorithm was applied to the summer top of atmosphere conversion:
The result of this algorithm is a range of values from −1 to 1, where positive values are considered to be water. These values were extracted and converted to polygons. As some variability also exists within the stage of rivers and lakes, the polygon mask is grown outwards by 1.5 pixels (45 m) to at least partially account for ice that may otherwise have been classified as an icing when it was actually surface water ice (Morse and Wolfe 2015). The MDSII result was converted from raster format to polygon format, and the water mask was erased from the MDSII polygons to leave only land-fast ice (henceforth considered to be icings).
After locating the icings, thermal imagery was used in an effort to distinguish spring and ground type icings by their thermal signatures. It was hypothesized that when thermal effects of elevation, vegetation, and water bodies were removed from the scene, warm anomalies of the snow pack would represent places where spring icings exist (because they continually discharge water), and that cold anomalies would represent areas where ground icings exist (because they have normally reached peak formation earlier in the winter). Sass et al. (2014) demonstrate this concept for a study area in the Prairie Parkland of northern Alberta, also using Landsat imagery. In the work of Sass et al. (2014), strong correlation was found between known spring locations and warm zones (presumably discharge zones) within the snowpack.
A thermal image, taken from the late winter of each year, was used to compute land surface temperature using the series of equations below and also illustrated in Fig. 6. For Landsat 4–5 TM, band 6 was used by default. For Landsat 8 OLI, band 10 was used. Although Landsat 8 OLI has two thermal bands (bands 10 and 11), the selection of band 10 for this work was arbitrary as only relative land surface temperature needed to be achieved. When calculating absolute land surface temperature is necessary, a split window algorithm that considers both bands 10 and 11 may be used instead (e.g., Du et al. 2015).
First, digital numbers were converted to spectral radiance:
Rλ = Spectral radiance (W·m−2·sr−1·μm−1)
Rmult = Sensor radiance multiplier (gain coefficient) (W·m−2·sr−1·μm−1)
Radd = Sensor radiance add (bias coefficient) (W·m−2·sr−1·μm−1)
DN = unaltered digital number (Chander et al. 2009)
Then to temperature in degrees Kelvin (to achieve at-satellite brightness temperature):
TU = Temperature (°K)
K2 = Thermal constant 2 (°K) (shown in Table 2)
K1 = Thermal constant 1 (mW·cm−2·sr−1·μm−1) (shown in Table 2)
Ru = Spectral radiance (W·m−2·sr−1·μm−1) (Wukelic et al. 1989)
Table 2.
Table 2. Thermal constants and bias/gain coefficients used for thermal anomaly calculations.

Note: TM, Thematic Mapper; OLI, Operational Land Imager.

To account for thermal interference related to vegetation density, at-satellite brightness temperatures were converted to land surface temperatures using Normalized Difference Vegetation Indices (NDVI) from the summer of each year. First, the NDVI was used to compute a proportion of vegetation (Pv, a scaled NDVI) where NDVImin is the lowest NDVI value in the scene and NDVImax is the highest NDVI value in the scene:
Then the proportion of vegetation was used to calculate the emissivity of vegetation (e):
At-satellite brightness temperature and the emissivity values were then used in conjunction to determine the land surface temperature (LST) where ω (wavelength of emitted radiance) is the median wavelength of the thermal band:
LST = Land surface temperature (°C)
BT = At satellite brightness temperature (°C)
e = Emissivity of vegetation (from eq. 7)
ω = Wavelength of emitted radiance (μm)
h = Planck’s constant (6.6261 × 10−34 kg·m2·s−1)
c = Velocity of light (2.9998 × 108 m·s−1)
s = Boltzmann constant (1.38 × 10−23 J·K−1) (Avdan and Jovanovska 2016)
Finally, the normality of the temperature distribution was verified using a Kolmogorov–Smirnov normality test and the LST values were standardized to z-scores to allow temporal comparison over each selected year.
After the methodology was successfully performed for the CMV, the icing extraction portion was replicated for just the Bogg Creek watershed for the year 2017 using higher resolution (5 m) RapidEye-3 optical imagery obtained from Planet Labs (image acquisition details shown in Table 1). Bogg Creek, shown in Fig. 1, is a sub-watershed of the CMV that is the focus area for shale oil and gas extraction. As such, it has been used as a field study location for the last several years. It forms part of a larger project whose goal is to establish baseline hydrologic monitoring conditions for the watershed in advance of oil and gas extraction.
Thermal anomalies were not obtained with RapidEye-3 imagery as a thermal band is not available on this sensor; only icings were obtained. As RapidEye-3 also does not have a SWIR, the NIR band was used for the algorithms required to discriminate icings. Although the NDSI and MDSII algorithms shown above utilize the SWIR band, several studies have also used the NIR band to distinguish snow from other land cover types (Bühler et al. 2015; Korzeniowska et al. 2017). Therefore, it is expected that use of the NIR band for this application is reasonable.

In situ methods

During the summer of 2018 and 2019, field work was completed in this watershed as part of the project described above. The data collected in the vicinity of TC-1 (Fig. 7) during these field campaigns were utilized to verify the groundwater discharge locations predicted by the remote method. This field data includes the following:
High-resolution thermal infrared imagery collected during a low-altitude aerial survey
A low-altitude (75–300 m) aerial survey was conducted across parts of the watershed to capture thermal infrared imagery of potential groundwater discharge. This survey utilized a hand-held FLIR Model T650sc thermal infrared camera used in conjunction with a Sony HDR-PJ430V video camera to record infrared and visual images of ground and water surfaces simultaneously, using the methods described by Conant (2009) and Rudolph (2019). It was assumed that during this time of the year (late August), groundwater temperatures would be colder than atmospheric, surface water, and ground surface temperatures. This temperature contrast between surface water features, vegetation, and emergent groundwater provides identification of colder thermal “anomalies” (potential springs) that may be less observable in the visible light spectrum. The survey was conducted along various transects of lake and stream edges, and along roads and seismic lines. A GPS tracker recording positions every second during the survey allowed for the geographic location of the anomaly to be determined later through post-processing. Not all anomalies were necessarily indicative of discharging groundwater; therefore, additional field verification was needed.
Vertical hydraulic gradients
Field verification of groundwater discharge conditions was achieved using the portable PPX48 PushPoint Sampler (MHE Products) to measure vertical hydraulic gradients and to collect groundwater samples. At a length of 122 cm, these samplers consist of 0.635 cm inner diameter stainless steel tubing, with a fine-tipped drive-point and slotted screen at one end, and a small welded handle and sampling port on the other. These portable and lightweight samplers act as a temporary piezometer that can be pushed by hand through soft sediments to take water and then removed afterwards. Small diameter (<0.635 cm) water-level tapes can be used within the PushPoint Sampler to measure groundwater levels or clear polymer tubing can be attached to the sampling port if water is under artesian conditions.
Groundwater geochemistry
At each potential spring location, ponded surface water was collected for major ion concentrations. As well, groundwater below the spring was collected from a depth of about 100 cm using a PushPoint sampler. Samplers were attached to 60 mL plastic syringes via 0.635 cm diameter clear polymer tubing, and samples drawn into the syringe by suction. Electroconductivity for each sample was determined with a handheld Oakton TDS/EC CON 10 m. Samples collected for metals were filtered within 6 h through 0.45 μm syringe filters and preserved with 1 mL 1:4 nitric acid. Collected samples were kept in coolers with icepacks and then refrigerated. Samples were sent within 3 days of collection to ALS Edmonton for analysis of major ions.
Fig. 7.
Fig. 7. Priority monitoring locations in the Bogg Creek watershed, given by center coordinates of icings. TC1, TC2, and TC3 represent the locations shown with thermal imaging in Fig. 13. Base map is a RapidEye-3 composite sourced from Planet Team 2017. Map generated using ArcGIS 10. [Colour online.]


Temperature distribution

Using the methods described in the previous section, temperature distributions and standardized temperature distributions (z-scores) were generated for years 2004, 2009, and 2016. An example from 2016 is shown in Figs. 8a and 8b. The z-scores were reclassified into integer intervals where lower negative values represent colder anomalies, higher positive values represent warmer anomalies, and 0 is the mean temperature of the scene. This was compared with both the distribution of icings by count and by areal coverage, using a χ2 goodness of fit test, as shown in Fig. 9. It was concluded that there was no significant difference (α = 0.05) in the observed distribution and the expected distribution (the distribution of the previous year) in any of the cases. It was also observed that the vast majority of icings, by both count and coverage, occur in weak anomalies (anomalies that are near the mean temperature) and that there is greater variation in both the discrete icing count and size between the years of 2016 and 2009 than between the years of 2009 and 2004. Because the vast majority of icings (by size and count) occurred near the mean temperature zones, spring-sourced and ground-sourced icings could not be definitively discriminated. This is explored further in the Discussion.
Fig. 8.
Fig. 8. Temperature distribution (a) and z-score result (b), shown as an example from 2016. Base maps are Landsat 8 Operational Land Imager images. Map generated using ArcGIS 10. [Colour online.]
Fig. 9.
Fig. 9. Comparison of thermal anomalies with discrete icing count for 2016 (a) and 2009 (b) and with the areal icing coverage for 2016 (c) and 2009 (d). Expected values are given by the icing distribution for the previously analyzed year. [Colour online.]

Icing recurrence in the CMV

Icings were delineated for each of the years 2004, 2009, and 2016, as shown in Fig. 10. Then, each of the two-year combinations were intersected to determine the percentage of recurring icings. All three years were also intersected to determine which icings would be most likely to recur again in the future. It was found that approximately 12.5% of the icings recur in all three years, and therefore are more likely to be of spring type. These icings represent the most promising field monitoring locations within the CMV. Figure 11 shows the icing overlap distribution for each year combination. Although the overlap range for all three years represents only the minority of icings, it is noted that the entire area covered by icings for each year is very consistent: 35.1 million square metres in 2004, 34.8 in 2009, and 32.9 in 2016 (σ = 0.97 m2).
Fig. 10.
Fig. 10. Example of the final icing result, shown for 2016, overlain on a Landsat 5 true colour (RGB) composite. Selected subset shows an enlarged view of the icing distribution. Map generated using ArcGIS 10. [Colour online.]
Fig. 11.
Fig. 11. Icing overlap (recurrence) for each possible year combination as a percent of total icing area. [Colour online.]

Icings in the Bogg Creek watershed

Following the analysis done for the entire CMV, icings were examined in more detail within the target shale oil and gas extraction watershed, Bogg Creek, using higher resolution RapidEye-3 imagery obtained from Planet Labs. A recurring icing complex was identified in the northwest portion of the watershed, as shown in Fig. 7. The development of icings in this complex from 2004, 2009, 2016, and 2017 is shown in Fig. 12. The icing complex is much less expansive in 2016 and 2017 than it is in 2004 and 2009. Nonetheless, the existence of the complex provides strong evidence to support the presence of spring-sourced groundwater discharge. In addition to the icings in this complex, a few other recurring icings within Bogg Creek were identified (also shown in Fig. 7).
Fig. 12.
Fig. 12. Progression of the Bogg Creek Icing Complex from 2004 to 2017 (icings shown in pink). Base map is a RapidEye-3 composite sourced from Planet Team 2017. Map generated using ArcGIS 10. [Colour online.]

Field verification in the Bogg Creek watershed

In August of 2018 and 2019, a field campaign was established in the Bogg Creek watershed and several of the icing locations given in Fig. 7 were visited with a high-resolution thermal camera. Additionally, field measurements were taken in the vicinity of TC-1. The following results support the hypothesis that recurrent icings are sourced from groundwater springs:
Thermal infrared imagery
Although high-resolution thermal camera imagery could not be collected during winter (due to poor flying conditions), very good agreement was found between the winter groundwater discharge locations given by the icings and the summer discharge locations given by the thermal camera. Selected summer thermal images are shown in Fig. 13. Using this imagery, cool zones are expected to be groundwater discharge points, as the land surface in the summer is warmer than the discharging groundwater.
Vertical hydraulic gradients
Artesian conditions were observed at two locations (samples GL1 and GL2) near TC-1 during data collection in 2018, with water flowing through the sampling port of the PushPoint Samplers. An upward vertical gradient was measured at TC-1 of about 0.04. This confirmed that these locations were in fact zones of groundwater discharge (springs). Upon returning to the general location a year later in 2019, artesian conditions were still noted to be occurring in the area (sample GL3 in Table 3).
Electroconductivity in the vicinity of TC-1 was high compared with other supra-permafrost groundwater (mean of 1678 μS·cm−1 compared with 929 μS·cm−1 for supra-permafrost groundwater) and surface water in the watershed (mean 206 μS·cm−1). Examining a piper plot (Fig. 14) of the major ions reveals a unique chemistry dominated by Ca/Na and HCO3/Cl with very little SO4. Geochemical facies of this groundwater, shown in Table 3, seems to suggest that it consists of a mixture of supra-permafrost and sub-permafrost groundwaters. However, the exact origin of these waters cannot yet be determined as Na and Cl may occur naturally in some of the unconsolidated overburden (data not shown) and is not necessarily unique to sub-permafrost groundwater. Nonetheless, this geochemical signature does suggest that subsurface-sourced water unique from other groundwater and surface water is discharging in the vicinity of TC-1.
Fig. 13.
Fig. 13. Thermal imagery collected from helicopter flight for TC1 and TC2 (left) and for TC3 (right). Groundwater discharge points are delineated with green circles. Imagery collected by B. Conant, Jr. [Colour online.]
Table 3.
Table 3. Major ion concentrations (mg·L−1) for spring sites in the vicinity of TC-1 (summer icing footprint).

Note: GL1 and GL2 were collected in 2018 whereas GL3 was collected in 2019. GL1, GL2, and GL3 refer to three sampling locations on the TC-1 icing footprint.

Fig. 14.
Fig. 14. Piper plot of site-wide endmembers, including supra-permafrost groundwater (asterisk), and sub-permafrost groundwater (X). Spring waters (green circles), collected from the GL sites described in Table 3, plot primarily within the Ca-HCO3 facies but display high concentrations of Na and Cl with very little SO4, making them unique within the watershed. Spring water appear to fall between supra-permafrost and sub-permafrost endmembers, although its origin cannot yet be definitively determined. [Colour online.]


Distinguishing spring and ground type icings

As aforementioned, it was hypothesized that icings would be segregated into either strong warm (2, 3, 4), or strong cold anomalies (−2, −3, −4) and that this property could be used to distinguish spring from ground icings. This discernment is important when establishing field monitoring sites as different mechanisms and physical processes within these two types of icings may affect the way that contaminants are moved to the surface and the way that icings interact with the flow system. This hypothesis forms its basis from the previously reviewed literature on icing phenomena, which indicates that ground icings may exhaust their supply of groundwater at some point during the winter and that spring icings have a continuous supply through the winter, reaching peak formation in the late spring (Carey 1973; Kane 1981; Yoshikawa et al. 2007). Weak anomalies (1, −1) are considered to represent either unaffected snowpack — snow that does not contain either active or inactive icings — or icings for which there is not strong enough evidence to conclude that they are one particular type.
As the majority of icings did fall within weak warm and cold anomalies, this suggests that the presence of icings, whether discharging groundwater or not, does not strongly affect the LST in the majority of cases for this study area. The icings falling in the weak warm category (1 > z > 0) may be discharging enough warm groundwater to fall above the mean temperature, but not enough to appear in the strong warm category. Similarly, the icings falling in the weak cold category (−1 < z < 0) may affect the LST enough to fall below the mean, but not enough to appear in the strong cold category. In theory, it is possible that the mean of the data (the 0-level z-score) could distinguish ground from spring icings; however, the results provide insufficient evidence to support this conclusion. If an icing were classified as spring type because it fell just above the mean z-score of 0, it may in fact be a ground icing that is still discharging warm water in the late winter or that has a low ice content. The same argument may be made for an icing classified as ground type because it fell just below the mean. Perhaps this is in fact a spring type icing that is discharging a small amount of groundwater on that particular day. These examples are to say that because the weak anomalies represent values closest to the mean, their contrast is not sufficient to discriminate icing types. Additionally, if the icings are distinguished from one another using the entire data set, there is room for a large degree of error in classifying those icings that lie in the weak anomalies. The icings that do fall into strong warm and cold anomalies may be more definitively classified as either spring or ground type based on current definitions of their physical occurrence, but because the majority do not fall into these strong anomalies, it is concluded that thermal anomalies may not be an ideal variable for discriminating icing type at this resolution. The spatial resolution of the Landsat thermal band is only 120 m (re-sampled to 30 m) and therefore may be too coarse to detect temperature changes resulting from groundwater discharge. This is especially true when the discharge locations are small and discrete, as they seem to be in this region. Fine-resolution thermal imagery may be more successful at discriminating spring and ground type icings in this region.
It should also be noted that there are other factors besides discharging groundwater that may play a role in determining LST. These include snow pack depth and density and differential thermal insulation (Pérez Díaz et al. 2015). Several measures were taken to try to account for variation in the expression of the anomalies. These included conversion to LST rather than brightness temperature, removal of water bodies from the distribution, selection of images with similar snow depth (as given by Environment Canada climate stations in Norman Wells, see Table 1), selection of cloud-free images, selection of images taken at the same time of day, and temperature standardization.

Recurring icings and implications

The fact that the minority of icings observed in the study region (12.5% by area) recur in each of the three observed years suggests that the majority of icings are intermittent or temporary. Although it is possible for ground type icings to recur in the same location year to year, it is expected that these recurring icings are more likely spring type as there is greater evidence to suggest that springs yield recurring icings. Upon examination of the three-year recurring icings in greater detail, no relationship was found between either extent or count with the strong anomalies. This result is not surprising given that the majority of icings from all years coincide with weak anomalies. Given these results, no definitive conclusions are drawn regarding the type of icings that three-year recurrences are likely to be. This does not mean that the three-year recurring icings are not of spring type, but it does suggest a large degree of variability in the hydrogeological regime, and that spring and ground type icings may not simply be distinguished based on the amount of water they are discharging in the late winter. Nonetheless, these recurring icings do provide promising locations for field monitoring, where further in situ work may be able to classify their type more definitively. It has been demonstrated by field observations that ground icings contain more organic material (appearing brown in colour) than do spring icings (Carey 1973), and also that carbonate precipitates are found in spring-derived icings (Hall 1980). These variables warrant further investigation in the CMV and in Bogg Creek. In this work, the observations that (1) 87.5% of the icings by area do not recur in all years, and (2) upon visual inspection, the distribution of anomalies is not consistent, further support the conclusion that these icings are intermittent and occur in different places on an annual basis.
In 2016, warm icings cover far more area than cold icings when compared with 2009 or 2004 (74.9%, 5.9%, and 10.7% for years 2016, 2009, and 2004, respectively). This may indicate that more groundwater springs are becoming active in the spring months as permafrost continues to thaw. Several studies in northern Canada and Alaska observe an increase in the amount of groundwater discharge as permafrost thaws and the active layer thickens. Walvoord et al. (2012) calibrated a model using data from the Yukon Flats Basin and found that projected future permafrost thaw results in increased river discharge from groundwater contributions, increased overall groundwater flux, and increased lateral reach of groundwater discharge in low-lying areas. Evans and Ge (2017), also using numerical modelling techniques, projected the amount of groundwater discharge that will occur on hillslopes following a 100-year warming period that induces permafrost thaw. They also observed significant increase in groundwater discharge following climate warming. These findings lend strong evidence to suggest that increasing proportions of warm icings in the CMV during late winter – early spring (when the surrounding landscape is still colder than discharging groundwater) could be, at the very least, partially resultant of permafrost thawing.

Icings in the Bogg Creek watershed

The occurrence of icings within Bogg Creek is minimal when compared with the occurrence of icings throughout the entire CMV. As aforementioned, there was one major icing hotspot (complex) identified in northwestern region of the watershed. The primary surficial texture class in this region is “fen”; therefore, it follows that the water table in this region is likely close to the surface and has strong potential for groundwater–surface interactions. Based on the recurrence of icings in this region, they are more likely of spring (permanent) type. Results of the in situ data collection support this conclusion. These include observations of discharging groundwater in 2018 and 2019, thermal anomalies located with an infrared camera, and unique geochemical characteristics that indicate groundwater here was unique from surface water and other supra-permafrost groundwater.
In 2004 and 2009, a substantially larger icing exists on the eastern side of the two lakes that is not present in 2016 and 2017 (Fig. 12). Although permafrost thaw often increases groundwater discharge, as previously discussed, complete thaw may result in drainage and depression of the water table thus restricting discharge. Permafrost thaw does tend to occur more readily in the immediate vicinity of lakes and water bodies as taliks beneath those water bodies provide pathways for groundwater flow that convects heat to the surrounding soil matrix (Scheidegger and Bense 2014). This is a possible cause of the decreased extent of the icing complex in 2016 and 2017. Precipitation patterns may also play a role in the occurrence of this icing complex; however, the relationship between icings and precipitation has not been clearly established in previous literature.

Characteristics of the hydrogeologic system in the CMV

The inconsistency in spatial icing distribution within the CMV does not appear to affect the total areal coverage of icings from year to year. The consistency in the total areal coverage of icings, calculated by summing the areas of each individual icing, suggests that the water available to supply icings is also consistent. Yoshikawa et al. (2007) concluded that the extent of icings monitored in the Brooks Range of Alaska is not nearly as sensitive to climate change as it is to source groundwater properties. It is suspected that this is also the case for icings in the CMV. As there have been numerous studies in the last few decades detailing dramatic climate changes in northern latitudes, yet total icing coverage remains consistent in the CMV, this lends evidence to the idea that climatic changes are not adversely affecting the volume of groundwater discharge during the winter. Rather, climatic changes affect the spatial distribution and size of icings. This information is valuable in that it provides an indication of the amount of groundwater that discharges to the surface during the winter (excluding that which discharges into surface water). This metric could be used in part to characterize the overall hydrologic system in this region, and to determine what proportion of water available in the system is derived from a subsurface source. Areal coverage on its own is not particularly useful for this purpose; therefore, an empirical equation first proposed by Sokolov in 1973 (Sokolov 1973) is suggested to compute the volume of icings:
where V = volume of icings, b and n = aufeis growth coefficients, and F = the area of icing coverage, derived from satellite imagery.
Empirical coefficients for this equation were determined for an aufeis field in the Brooks Range of Alaska by Hall and Roswell (1981):
In situ studies would be useful to refine these coefficients for more precise use in the CMV, and (or) to assess the applicability of the Brooks Range coefficients for use in this study area. An understanding of the volume of water that discharges throughout the winter would be valuable in gaining a more in-depth picture of the entire hydrogeologic system, and to aid with numerical modelling efforts.

Conclusions and recommendations

The locations, sizes, and recurrence of areal features that are hypothesized to be icings within the CMV and the Bogg Creek watershed were determined for years 2004, 2009, 2016, and 2017 using a combination of Landsat 4–5 TM, Landsat 8 OLI, and RapidEye-3 imagery. Additionally, thermal anomalies were calculated using the Landsat thermal band for years 2004, 2009, and 2016. Following these results, several important conclusions regarding icing occurrence were made:
The amount of winter groundwater discharging to the land surface in this region is stable from year to year, even though the spatial distribution of icings is not.
Ground and spring icings may not function in the same way across all regions; ground icings may still discharge water late in the winter, and spring icings may undergo periods of little or no discharge.
Due to variability in the distribution of recurring icings within thermal anomalies, and the large overall proportion of icings in weak anomalies for all years, coarse resolution thermal imagery should not be used to definitively distinguish ground and spring type icings. It should be used only to provide suggestions as to the mechanisms governing their formation. Thermal camera imaging aboard a helicopter, by contrast, was very effective at verifying predicated icing locations during summer.
Priority icing/groundwater discharge monitoring locations were established for the Bogg Creek Watershed and several were verified and examined in further detail using in situ methods. These priority monitoring locations are annually recurring icings and are likely spring sourced.
This work demonstrated a remote-sensing based method of identifying expected groundwater discharge zones. Moving forward, monitoring the icing distribution in this region with even higher resolution imagery, from unmanned aerial vehicles, for example, would be of great value. Further exploration of a combination of remote sensing and in situ methods conducted during winter (such as isotopic and geochemical methods) to resolve the source of the groundwater forming icings depicted in this work would assist in the overall understanding of the hydrologic regime. Additionally, monitoring selected icings and their growth throughout the winter season would allow for the establishment of growth coefficients that could more precisely estimate the volume of water that reaches the ground surface during winter.


The authors gratefully acknowledge funding from the Natural Sciences and Engineering Research Council of Canada, Global Water Futures through the Northern Water Futures project and the Government of the Northwest Territories Environmental Studies Research Fund. We are highly appreciative of the technical and logistical advice provided by Andrew Applejohn and Bruce Hanna (Government of the Northwest Territories). Technical support and field work facilitation were also graciously provided by the Sahtu Renewable Resources Board and we are thankful to Sahtu communities for reviewing and approving our research licenses (ARI Scientific Research License Nos. 16345 and 16581). The work would not have been possible without the logistical and technical support of personnel from Husky Energy. Many thanks are also owed to Brewster Conant, Jr. and Aaron Vandenhoff (Department of Earth and Environmental Sciences, University of Waterloo) for their assistance with field work and data interpretation and to Planet Labs ( for the provision of the RapidEye-3 imagery. Finally, we thank the anonymous reviewers of the Canadian Journal of Earth Sciences for their time in providing an in-depth review of this work.


AMEC. 2014. Central Mackenzie Valley Subsurface Groundwater Baseline Study. Indigenous and Northern Affairs Canada and the North, Northern Petroleum Resources. Retrieved from
Avdan U. and Jovanovska G. 2016. Algorithm for automated mapping of land surface temperature using LANDSAT 8 satellite data. Journal of Sensors, 2016(2): 1–8.
Bühler Y., Meier L., and Ginzler C. 2015 Potential of operational high spatial resolution near-infrared remote sensing instruments for snow surface type mapping. IEEE Geoscience and Remote Sensing Letters, 12(4): 821–825.
Burgess, M.M., and Smith, S.L. 2000. Shallow ground temperatures. In The physical environment of the Mackenzie Valley, Northwest Territories: a base line for the assessment of environmental change. Edited by L.D. Dyke and G.R. Brooks. Geological Survey of Canada Bulletin 547, pp. 89–103.
Bush, E., and Lemmen, D.S. 2019. Canada’s changing climate report. Government of Canada, Ottawa, Ont.
Cao Y. and Liu C. 2006. Normalized difference snow index simulation for snow-cover mapping in forest by geosail model. Chinese Geographical Science, 16(2): 171–175.
Carey, K.L. 1973. Icings developed from surface and ground water. Corps of Engineers, U.S. Army Cold Regions Research and Engineering Laboratory.
Carlson T.N. and Ripley D.A. 1997. On the relation between NDVI, fractional vegetation cover, and leaf area index. Remote Sensing of Environment, 62(3): 241–252.
Chander G., Markham B.L., and Helder D.L. 2009. Summary of current radiometric calibration coefficients for Landsat MSS, TM, ETM+, and EO-1 ALI sensors. Remote Sensing of Environment, 113(5): 893–903.
Conant, B., Jr. 2009. Infrared thermography to identify groundwater discharges for fisheries habitat assessments in the Northwest Territories. Conference proceeding, Groundwater Summit 2009.
Connon R., Devoie É., Hayashi M., Veness T., and Quinton W. 2018. The influence of shallow taliks on permafrost thaw and active layer dynamics in subarctic Canada. Journal of Geophysical Research: Earth Surface, 123(2): 281–297.
Côté, M.M., Duchesne, C., Wright, J.F., and Ednie, M. 2013. Digital compilation of the surficial sediments of the Mackenzie Valley corridor, Yukon Coastal Plain, and the Tuktoyaktuk Peninsula. Geological Survey of Canada, Open File 7289.
Cuenca J. and Sobrino J.A. 2004. Experimental measurements for studying angular and spectral variation of thermal infrared emissivity. Applied Optics, 43(23): 4598.
Derksen, C., Burgess, D., Duguay, C., Howell, S., Mudryk, L., Smith, S. et al. 2018. Changes in snow, ice, and permafrost across Canada. Chapter 5. In Canada’s changing climate report. Edited by E. Bush and D.S. Lemmen. Government of Canada, Ottawa, Ont. pp. 194–260.
Du C., Huazhong R., Qin Q., Meng J., and Zhao S. 2015. A practical split-window algorithm for estimating land surface temperature from Landsat 8 data. Remote Sensing, 7(1): 647–665.
Environment Canada. 2018. Government of Canada Historical Data. Available from
Environmental Protection Agency (EPA). 2016. Hydraulic fracturing for oil and gas: impacts from the hydraulic fracturing water cycle on drinking water resources in the United States. EPA/600/R-16/236F, Final Report. U.S. Environmental Protection Agency, Washington, D.C.
Evans S.G. and Ge S. 2017. Contrasting hydrogeologic responses to warming in permafrost and seasonally frozen ground hillslopes. Geophysical Research Letters, 44(4): 1803–1813.
Fallas, K.M., and MacNaughton, R.B. 2013. Geology, Norman Wells (Southeast) Northwest Territories, Geological Survey of Canada, Canadian Geoscience Map 100.
Fallas, K.M., Hadlari, T., and MacLean, B.C. 2013. Geology, Carcajou Canyon (northeast), Northwest Territories. Geological Survey of Canada. Canadian Geoscience Map 95.
Fraser R.H., Kokelj S.V., Lantz T.C., McFarlane-Winchester M., Olthof I., and Lacelle D. 2018. Climate sensitivity of high arctic permafrost terrain demonstrated by widespread ice-wedge thermokarst on Banks Island. Remote Sensing, 10(6): 954.
Hall D.K. 1980. Mineral precipitation in North Slope river icings. Arctic, 33(2): 343–348.
Hall D.K. and Roswell C. 1981. The origin of water feeding icings on the eastern North Slope of Alaska. Polar Record, 20(128): 433–438.
Hall D.K., Riggs G.A., and Salomonson V.V. 1995. Development of methods for mapping global snow cover using moderate resolution imaging spectroradiometer data. Remote Sensing of Environment, 54(2): 127–140.
Hayes, B.R., and Dunn, J. 2012. Deep Subsurface Saline Aquifer Characterization, Central Mackenzie Valley, Northwest Territories. NWT Open File 2012-07.
Jepsen S.M., Walvoord M.A., Voss C., and Rover J. 2015. Effect of permafrost thaw on the dynamics of lakes recharged by ice-jam floods: case study of Yukon Flats, Alaska. Hydrological Processes, 30(11): 1782–1795.
Kane D.L. 1981. Physical mechanics of aufeis growth. Canadian Journal of Civil Engineering, 8(2): 186–195.
Korzeniowska K., Bühler Y., Marty M., and Korup O. 2017 Regional snow-avalanche detection using object-based image analysis of near-infrared aerial imagery. Natural Hazards Earth System Sciences, 17(10): 1823–1836.
McFeeters S.K. 1996. The use of the Normalized Difference Water Index (NDWI) in the delineation of open water features. International Journal of Remote Sensing, 17(7): 1425–1432.
Morse P.D. and Wolfe S.A. 2015. Geological and meteorological controls on icing (aufeis) dynamics (1985 to 2014) in subarctic Canada. Journal of Geophysical Research F: Earth Surface, 120(9): 1670–1686.
NWT Government 2012. Climate observations in the Northwest Territories (1957–2012). Retrieved from
Pérez Díaz C.L., Lakhankar T., Romanov P., Khanbilvardi R., and Yu Y. 2015. Evaluation of VIIRS land surface temperature using CREST-SAFE air, snow surface, and soil temperature data. Geosciences, 5(4): 334–360.
Planet Team. 2017. Planet Application Program Interface: in Space for Life on Earth. San Francisco, Calif. Available from
Pollard, W.H., and van Everdingen, R.O. 1992. Formation of seasonal ice bodies. In Periglacial geomorphology. Edited by J.C. Dixon and A.D. Abrahams. John Wiley & Sons, Chichester, U.K. pp. 281–304.
Rouse, W.R., Bello, R.L., and Lafleur, P.M. 1997. The Low Arctic and Subarctic. In The surface climates of Canada. Edited by W.G. Bailey, T. Oke, and R. Rouse. McGill Queen’s University Press. pp. 198–211.
Rudolph, D.L. 2019. Regional hydrologic and ecologic characterization and baseline assessment of remote northern Canadian terrain in advance of shale oil and gas development. Second Annual Report to NWT ESRF Management Board.
Rudolph, D.L., Lotimer, T., and Barker, J.F. 2016. Final report baseline hydrogeological evaluation of Central Mackenzie Valley oil and gas exploration areas, Sahtu Region, Northwest Territories. Available from
Salomonson V.V. and Appel I. 2004. Estimating fractional snow cover from MODIS using the normalized difference snow index. Remote Sensing of Environment, 89(3): 351–360.
Sass G.Z., Creed I.F., Riddell J., and Bayley S.E. 2014. Regional-scale mapping of groundwater discharge zones using thermal satellite imagery. Hydrological Processes, 28(23): 5662–5673.
Scheidegger J.M. and Bense V.F. 2014. Impacts of glacially recharged groundwater flow systems on talik evolution. Journal of Geophysical Research: Earth Surface, 119(4): 758–778.
Shusan L., Benson C., Shapiro L., and Dean K. 1997. Aufeis in the Ivishak River, Alaska, mapped from satellite radar interferometry. Remote Sensing of Environment, 60(2): 131–139.
Sokolov, A.A. 1973. Experience on hydrologic substantiation of projected reservoirs. In Man_made lakes: their problems and environmental effects. Vol. 17. Edited by W.C. Ackermann, G.F. White, E.B. Worthington, and J.L. Ivens. pp. 264–271. Geophysical Monograph Series.
Walvoord M.A., Voss C.I., and Wellman T.P. 2012. Influence of permafrost distribution on groundwater flow in the context of climate-driven permafrost thaw: example from Yukon Flats Basin, Alaska, United States. Water Resources Research, 48(7): 1–17.
Woo M.-K. and Xia Z. 1995. Suprapermafrost groundwater seepage in gravelly terrain, Resolute, NWT, Canada. Permafrost and Periglacial Processes, 6(1), 57–72.
Wukelic G.E., Gibbons D.E., Martucci L.M., and Foote H.P. 1989. Radiometric calibration of Landsat Thematic Mapper thermal band. Remote Sensing of Environment, 28(C): 339–347.
Yoshikawa K., Hinzman L.D., and Kane D.L. 2007. Spring and aufeis (icing) hydrology in Brooks Range, Alaska. Journal of Geophysical Research: Biogeosciences, 112(4): 1–14.
You Y., Yu Q., Pan X., Wang X., Guo L., and Wu Q. 2017. Thermal effects of lateral supra-permafrost water flow around a thermokarst lake on the Qinghai-Tibet Plateau: thermal effects of lateral supra-permafrost water flow. Hydrological Processes, 31(13): 2429–2437.
Zhang, X., Flato, G., Kirchmeier-Young, M., Vincent, L., Wan, H., Wang, X. et al. 2019. Changes in temperature and precipitation across Canada. In Canada’s changing climate report. Edited by E. Bush and D.S. Lemmen. Government of Canada, Ottawa, Ontario. pp. 112–193.

Information & Authors


Published In

cover image Canadian Journal of Earth Sciences
Canadian Journal of Earth Sciences
Volume 58Number 2February 2021
Pages: 105 - 121


Received: 17 September 2019
Accepted: 11 May 2020
Accepted manuscript online: 23 June 2020
Version of record online: 23 June 2020

Key Words

  1. icings
  2. aufeis
  3. permafrost
  4. climate change
  5. groundwater discharge


  1. aufeis
  2. pergélisol
  3. changement climatique
  4. émergence d’eau souterraine



Brittney K. Glass [email protected]
Department of Earth and Environmental Sciences, University of Waterloo, 200 University Avenue West, Waterloo, ON N2L 3G1, Canada.
David L. Rudolph
Department of Earth and Environmental Sciences, University of Waterloo, 200 University Avenue West, Waterloo, ON N2L 3G1, Canada.
Claude Duguay
Department of Geography and Environmental Management, University of Waterloo, 200 University Avenue West, Waterloo, ON N2L 3G1, Canada.
Andrew Wicke
Department of Earth and Environmental Sciences, University of Waterloo, 200 University Avenue West, Waterloo, ON N2L 3G1, Canada.


Copyright remains with the author(s) or their institution(s). This work is licensed under a Creative Commons Attribution 4.0 International License (CC BY 4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original author(s) and source are credited.
A correction was made to the e-First version of this paper on 11 December 2020 prior to the final issue publication. The current online and print versions are identical and both contain the correction.

Metrics & Citations


Other Metrics


Cite As

Export Citations

If you have the appropriate software installed, you can download article citation data to the citation manager of your choice. Simply select your manager software from the list below and click Download.

Cited by

1. Groundwater vulnerability in the Yukon and Northwest Territories, Canada
2. Distribution Of Icings In The Northern (Russian) Part Of The Selenga River Basin And Their Role In The Functioning Of Ecosystems And Impact On Settlements
3. Thaw-induced impacts on land and water in discontinuous permafrost: A review of the Taiga Plains and Taiga Shield, northwestern Canada
4. Complex Vulnerabilities of the Water and Aquatic Carbon Cycles to Permafrost Thaw

View Options

View options


View PDF

Get Access

Login options

Check if you access through your login credentials or your institution to get full access on this article.


Click on the button below to subscribe to Canadian Journal of Earth Sciences

Purchase options

Purchase this article to get full access to it.

Restore your content access

Enter your email address to restore your content access:

Note: This functionality works only for purchases done as a guest. If you already have an account, log in to access the content to which you are entitled.





Share Options


Share the article link

Share on social media

Cookies Notification

We use cookies to improve your website experience. To learn about our use of cookies and how you can manage your cookie settings, please see our Cookie Policy.