Digital Commons @ Michigan Tech Digital Commons @ Michigan Tech Geophysical Observations of Taliks Below Drained Lake Basins on Geophysical Observations of Taliks Below Drained Lake Basins on the Arctic Coastal Plain of Alaska the Arctic Coastal Plain of Alaska

Lakes and drained lake basins (DLBs) together cover up to ∼ 80% of the western Arctic Coastal Plain of Alaska. The formation and drainage of lakes in this continuous permafrost region drive spatial and temporal landscape dynamics. Postdrainage processes including vegetation succession and permafrost aggradation have implications for hydrology, carbon cycling, and landscape evolution. Here, we used surface nuclear magnetic resonance (NMR) and transient electromagnetic (TEM) measurements in conjunction with thermal modeling to investigate permafrost aggradation beneath eight DLBs on the western Arctic Coastal Plain of Alaska. We also surveyed two primary surface sites that served as nonlake affected control sites. Approximate timing of lake drainage was estimated based on historical aerial imagery. We interpreted the presence of taliks based on either unfrozen water estimated with surface NMR and/or TEM resistivities in DLBs compared to measurements on primary surface sites and borehole resistivity logs. Our results show evidence of taliks below several DLBs that drained before and after 1949 (oldest imagery). We observed depths to the top of taliks between 9 and 45 m. Thermal modeling and geophysical observations agree about the presence and extent of taliks at sites that drained after 1949. Lake drainage events will likely become more frequent in the future due to climate change and our modeling results suggest that warmer and wetter conditions will limit permafrost aggradation in DLBs. Our observations provide useful information to predict future evolution of permafrost in DLBs and its implications for the water and carbon cycles in the Arctic.


Introduction
Lakes and drained lake basins (DLBs) combined are estimated to cover up to ∼80% of the western Arctic Coastal Plain of Alaska (∼30,000 km 2 ) Hinkel et al., 2005;Jones & Arp, 2015). There are a variety of lake types in the Arctic, but the most common are thermokarst lakes in lowland regions with ice-rich permafrost Kling, 2009) that form due to permafrost thaw and surface subsidence. Deeper lakes developed in permafrost terrain are often underlain by layers or bodies of perennially unfrozen ground below the lake bed known as a talik (van Everdingen, 1998). Arctic lakes can persist for thousands of years, but, due to ongoing margin expansion and other landscape changes, they eventually drain laterally to create a mosaic of extant lakes and DLBs (Hinkel et al. 2007;Mackay, 1992).
Lakes can serve as one of the few nonclimatic thawing agents in continuous permafrost terrain. In a lake, if the water column does not typically freeze to lake bottom (floating ice) during winter, permafrost will thaw below the lake to form a talik. Because taliks can extend to depths of hundreds of meters through the permafrost (Mackay, 1997), direct measurements are challenging and, to our knowledge, there are only a few borehole studies that investigated sublake taliks (e.g., Brewer, 1958;Brewer et al., 1993;Burn & Smith, 1990;Heslop et al., 2015;Mackay, 1997;Roy-Léveillée & Burn, 2017). Thermal modeling indicates that a lake talik thickness between 28 and 53 m will form over 3,000 years in the Arctic Coastal Plain of Alaska given lakebed temperatures between 1 and 3 °C (Ling & Zhang, 2003). The timing of talik refreeze depends on a combination of factors including the climate, talik thickness, lithology, and temperature of the surrounding permafrost (Ling & Zhang, 2004). In general, remnant taliks are estimated to refreeze over several decades through a combination of top-down and bottom-up permafrost aggradation (Ling & Zhang, 2004;Mackay, 1992). Ling and Zhang's (2004) thermal modeling indicates that following lake drainage, a talik can refreeze in 40-157 years, where top-down refreezing can be 6-8 times faster than bottom-up. Other thermal modeling results suggest that under future warmer climate scenarios, the number of floating ice lakes underlain by taliks may increase on the Arctic Coastal Plain of Alaska (Matell et al., 2013). Considering the broad diversity of DLB locations, sizes, and ages, and the scarcity of direct measurements of taliks, we are motivated to acquire a better understanding of where taliks exist below DLBs and the rates of talik refreezing given regional climatic conditions over the last several decades.
DLBs are a viable target for geophysical observations. Variations in geophysical properties depend on many parameters including lithology, porosity, ice content, unfrozen pore water content, pore water salinity, pore pressure, and temperature (Minsley et al., 2015;Yoshikawa et al., 2006). Thawed sediments are assumed to have electrical resistivities lower than frozen sediments due to the lower resistivity of liquid water compared to ice (Ross et al., 2007). In Arctic coastal environments, both thawed and frozen sediments can present low electrical resistivities due to saline pore water (Creighton et al. 2018;Keating et al., 2018;Overduin et al., 2012;Ross et al., 2007) that may make the discrimination of frozen, unfrozen or partially unfrozen pore water ambiguous. Electrical and electromagnetic methods, such as the electrical resistivity tomography (ERT) and transient electromagnetic (TEM), respectively, measure a physical property that is indirectly related to the presence of water; i.e., the measured electrical resistivity of the subsurface layers represents the bulk resistivity of the sediments or rocks plus the frozen or unfrozen pore water (Minsley et al., 2015), which can complicate interpretation. Surface nuclear magnetic resonance (NMR) is the only noninvasive geophysical method that can directly detect unfrozen water as it measures the response of hydrogen protons on liquid water (Behroozmand et al., 2015;Walsh, 2008). However, surface NMR has a limitation to detect water in very small pore spaces or partially unfrozen pore water, which can produce a signal too short to be measured.
Electrical and electromagnetic geophysical methods have been used to image taliks in permafrost areas based on the electrical resistivity distribution of the subsurface (Creighton et al., 2018;Minsley et al., 2012Minsley et al., , 2015Yoshikawa et al., 2006;Yoshikawa & Hinzman, 2003;You et al., 2017). Surface NMR can be used in combination with electrical resistivity measurements to support the interpretation of frozen versus unfrozen pore water on permafrost environments (Creighton et al., 2018;Keating et al., 2018;Parsekian et al., 2013Parsekian et al., , 2019. Parsekian et al. (2013) combined measurements and synthetic forward modeling to demonstrate that surface NMR can detect different types of thawed layers in permafrost. For example, Creighton et al. (2018) and Parsekian et al. (2019) used surface NMR and TEM methods to identify taliks below thermokarst lakes on the Arctic Coastal Plain of Alaska. Recent studies in the Arctic have combined surface NMR and ground-penetrating radar (GPR) to estimate glacial water in an ice sheet firn in Greenland (Legchenko et al., 2018), and identify frozen and thawed sediments below river aufeis (icings) on floodplain of the North Slope of Alaska (Terry et al., 2020).
Other geophysical methods have been used to successfully image taliks. For example, Schwamborn et al. (2002) used seismic reflection to determine a talik thickness of 95 m below a thermokarst lake in continuous permafrost of Siberia. However, to our knowledge, it is the only published study that used seismic reflection to image a lake talik. GPR and ERT were combined to image a shallow (<25 m depth) open talik (talik that penetrates permafrost completely) below a thermokarst lake in the discontinuous permafrost zone of Alaska (Yoshikawa & Hinzman, 2003) and Qinghai-Tibet Plateau, China (You et al., 2017), which can connect lakes to the subpermafrost groundwater. Another relevant example is the study by Yoshikawa et al. (2006) that combined seven different geophysical surveys, including GPR, ERT, seismic refraction, very low frequency electromagnetic, surface NMR, airborne and ground-based TEM, to investigate their abilities to distinguish permafrost, massive ground ice, and talik, at three pingos in the discontinuous permafrost zone of Alaska. Their results suggest that, depending on the depth of investigation, GPR and ERT are the most useful methods for detecting massive ground ice. However, the ERT electrodes require a good ground contact, which is difficult in frozen ground.
Despite the importance of lake taliks for carbon and water cycling in the Arctic, there are still significant observational knowledge gaps related to talik formation and the refreezing rate after lake drainage. Understanding how a talik refreezes after lake drainage may help to predict its contribution to the permafrost carbon feedback. Given the ubiquity of lakes and DLBs in the Arctic Coastal Plain of Alaska and the paucity of information about postdrainage permafrost aggradation, we aim to answer the following questions: (1) What is the rate of permafrost aggradation below DLBs with variable drainage histories? (2) What is the geophysical character of a talik and how does it change through the refreezing process? To answer these questions and test our predictions of talik refreezing following lake drainage, we used surface NMR and TEM measurements in conjunction with talik refreeze modeling to investigate permafrost aggradation beneath DLBs on the western Arctic Coastal Plain of Alaska ( Figure 1). To our knowledge, there are very few geophysical measurements of DLB taliks, and there are none using surface NMR alone or combined with TEM. The combination of these two noninvasive geophysical methods is ideal to achieve our objective to distinguish frozen versus thawed subsurface layers.

Study Sites
We investigated eight DLBs and two primary surface sites (terrain that has not been eroded) in the Teshekpuk and Inigok regions on the Arctic Coastal Plain of Alaska (Figure 1), all within the continuous permafrost zone. Based on surface characteristics including relief and vegetation, primary surface sites are assumed to have not experienced lake formation and drainage since at least the beginning of the Holocene (Brown, 1965;Kanevskiy et al., 2013), and they were used here as experimental control sites and for comparison with DLB measurements. To aid our geophysical survey interpretations, one borehole in each region was used to infer lithology, resistivity, and temperature logs (Clow, 2014;Gryc, 1988).
The Teshekpuk region has low topographic relief and predominantly Holocene ice-rich silty soils Kanevskiy et al., 2013), where DLBs occupy about 62% of the landscape and lakes 23% (Jones & Arp, 2015). At the Drew Point borehole (70.880°, −153.904°), the permafrost thickness is ∼320 m (Clow, 2014), the average permafrost temperature is ∼−5 °C (Clow, 2014), and the average temperature at 1.2 m depth was −8.0 °C between 1999 and 2010 (Farquharson et al., 2016). The Inigok region is characterized by late Pleistocene and early Holocene deposits of aeolian sand (Carter et al., 1981;Galloway & Carter, 1993) and relatively low ground ice content (Kanevskiy et al., 2013). At the Inigok North borehole (70.265°, −152.766°), the permafrost thickness is ∼286 m and the average permafrost temperature is ∼−5 °C (Clow, 2014). At another borehole site in the region (69.990°, −153.094°), the average temperature at 1.2 m depth was −5.9 °C between 1999 and 2010 (Farquharson et al., 2016). One of the DLBs in the Inigok region contains a pingo, which is associated with permafrost aggradation consisting of an ice-cored mound formed from closed system freezing of pore water in the talik following lake drainage, and upward deformation of the frozen overburden (Mackay, 1973(Mackay, , 1998. The coordinates around the center and the drainage history for the studied DLBs are shown in Table 1. The drainage periods were determined based on historical aerial imagery , where the oldest photos are from 1949. The mechanisms and degrees of drainage (partial or complete) were inferred based on either historical aerial imagery or field observations . Peatball Bluff (70.708°, −153.972°) and Inigok Camp (70.003°, −153.085°) are the primary surface sites (underlain by marine silt and aeolian sand deposits, respectively) ( Figure 1).

Methods
We used surface nuclear magnetic resonance (NMR) and transient electromagnetic (TEM) measurements in conjunction with talik refreeze modeling to investigate remnant taliks beneath DLBs on the western Arctic Coastal Plain of Alaska. Both surface NMR and TEM measurements consisted of a one-dimensional investigation of subsurface geophysical properties. The measurements were performed around the center of each DLB (coordinates in Table 1).

Transient Electromagnetic
The TEM method uses a transmitter wire loop on the Earth's surface to produce an electromagnetic field that interacts with the subsurface and measures the decay of the induced electromagnetic field through time with a receiver wire loop. The measured induced field is inverted to estimate a vertical subsurface distribution of electrical resistivity. Geophysical inversion is an optimization process that uses the field data to estimate the subsurface distribution of a physical property (Aster et al., 2018).
The TEM depth of investigation (DOI) is partially controlled by the transmitter loop area and the subsurface resistivity structure, where the more resistive the substrate, the deeper the DOI (Christiansen et al., 2006). Electromagnetic noise sources were not present at our remote study sites. A full explanation of the underlying physical principles of the TEM method can be found in Fitterman and Stewart (1986) and Christiansen et al. (2006 TEM data were acquired with the ABEM WalkTEM instrument (Guideline Geo, Stockholm, Sweden) at six sites in the Teshekpuk region and at four sites in the Inigok region in April of 2016, 2017, and 2019 ( Figure 1 and Table 2), when the ground is snow covered and near-surface ground temperatures are at their annual minimum. All soundings used a nested central loop configuration that consisted of three loops: (1) a circular transmitter loop with an area of 1,600 m 2 (single turn 160 m long wire), (2) a square receiver loop with an effective area of 200 m 2 (40 m long wire, 10 m × 10 m with two turns internally), and 3) a second square receiver loop with an effective area of 5 m 2 (40 m long wire, 0.5 m × 0.5 m with 20 turns internally).
The TEM raw data were processed and inverted with the SPIA software (AGS, Aarhus, Denmark). First, each dataset was preprocessed by removing unreliable early and late time gates to increase the signal-tonoise ratio (SNR). Next, the data were inverted to estimate the electrical resistivities of the subsurface. The inversion is based on a one-dimensional vertical resistivity constraint using a smooth model with 20 layers of fixed thicknesses. The data residual (δ) for all the inversion model results varied between 0.4 and 0.9 (Table 2). For all TEM acquisitions, the DOI was deeper than 100 m.
We used the modified Archie Equation (Archie, 1942;Ruffet et al., 1995) to estimate volumetric water content (VWC) using the TEM resistivities: where  is the bulk resistivity measured with TEM,  w the pore water resistivity,  s the surface resistivity (related to surface conduction effect), Φ the porosity, and m the cementation exponent. Assuming a full saturation,  Φ VWC , the equation (1) becomes:

Surface Nuclear Magnetic Resonance
Surface NMR is the only noninvasive surface-based geophysical method that directly detects unfrozen water belowground. It uses a wire loop on the Earth's surface to produce an electromagnetic excitation of subsurface hydrogen protons in unfrozen water aligned at equilibrium with the geomagnetic field. Subsequently, the resonating magnetic field generated when the protons relax back to the alignment with the geomagnetic field is measured with the same loop. The excitation pulse frequency is tuned to the Larmor frequency (f L ) of hydrogen. The f L is defined locally for each measurement and it is the product of the geomagnetic field total intensity and the hydrogen gyromagnetic ratio. In this work, all surface NMR measurements are related to unfrozen water because the relaxation time of ice is too short to be detected with surface NMR. Moreover, surface NMR has an inherent limitation in its ability to detect water in very small pore spaces, such as in claystone, or partially unfrozen pore water, which results in relaxation times too short to be measurable.
The free induction decay (FID) is the most common surface NMR pulse sequence, which uses a single excitation pulse. The pulse moment is the amplitude and duration of the transmitted excitation pulse. Longer pulse moments can increase the DOI but limit the ability to resolve short relaxation times. The surface NMR DOI is partially controlled by the loop diameter (∅), the excitation pulse duration (t P ), and the subsurface electrical properties. Electrically conductive materials attenuate the signal, decreasing the DOI. Moreover, the surface NMR measurements may be affected by gradients in the geomagnetic field, which can induce the protons to oscillate at a slightly different f L . The Geophysical Institute of the University of Alaska has several stations across the state to monitor the geomagnetic field activity and, based on their publicly available archive (https://www.gi.alaska.edu/monitors/magnetometer/archive), there were no significant gradients on the geomagnetic field during the surface NMR data acquisitions. Therefore, spatial gradients can be neglected (Parsekian et al., 2019). Magnetic field gradients as a function of depth are unknown but assumed to be small given that the sedimentary geologic substrate has a low magnetic mineral content at our study sites (Gryc, 1988).
The exponential relaxation of hydrogen protons is inverted to estimate the vertical distribution of volumetric water content (VWC) that is related to the signal amplitude, and the relaxation pulse duration (  2 T ), which is proportional to the mean pore size (Hertrich, 2008). More details about the surface NMR method can be found in Walsh (2008) and Behroozmand et al. (2015).
The surface NMR FID data were acquired with the GMR instrument (Vista Clara, Mukilteo, WA, USA) at four sites in the Inigok region in April 2015 and at three sites in the Teshekpuk region in April 2016, with 75 and 90 m circle loop diameters (∅), respectively ( Table 2). The instrument dead time (interval between the end of the excitation pulse and the start of data recording) is 5.5 ms. Due to the lack of nearby electromagnetic noise sources, it was not necessary to use noise canceling loops. We used a t P of 20 ms for most of the sites, except for Imakruak Basin where 10 ms was used to increase the SNR, and for Three Creatures Basin where 40 ms was measured to compare with 20 ms. The data were acquired with either 4 or 8 stacks. Before each acquisition, we used a Geometrics G816 magnetometer to determine the geomagnetic field magnitude at five locations within each loop to estimate the f L (Table 2). After each measurement started, we checked the peak of the spectrum to ensure that we were transmitting within <0.5 Hz of on-resonance (Grombacher & Knight, 2015;Walbrecker et al., 2011). The geomagnetic field inclination values were calculated with the NOAA IGRF model based on the location of each surface NMR measurement and, in general, were about 80.8° at the Teshekpuk region and 80.4° at the Inigok region.
The surface NMR raw data were processed and inverted with the GMR processing software (Vista Clara, Mukilteo, WA, USA; Walsh, 2008). First, each dataset was processed using a 200 Hz bandpass filter resulting in dead times between 15.7 and 15.9 ms. Bad stacks were discarded resulting in either 2 or 8 stacks per site. Then the data were inverted using both single-exponential and multiexponential fitting algorithms (Walsh, 2008). We inverted the data considering the relaxation during and after the excitation pulse (Grombacher et al., 2017;Walbrecker et al., 2009). The inversion kernels were calculated using the local geomagnetic field magnitude and inclination angle. For all sites except Blustery Basin, we also acquired TEM data and the measured subsurface average resistivity  (Table 3) was used in the surface NMR inversion (Parsekian et al., 2013). The DOI and the depth uncertainties were based on resolution matrices (Müller-Petke & Yaramanci, 2008) estimated during the inversion using singular value decomposition, which considers the regularization parameters (signal amplitude and noise level) and the measurement configuration (∅,  , and pulse moments). The maximum DOI was 50 m ( Table 2). The GMR software does not provide the data residual for the inversion model results. However, we have high confidence that our results are consistent because previous works have used the same surface NMR equipment on the North Slope of Alaska (Creighton et al., 2018;Parsekian et al., 2019;Terry et al., 2020) and the same processing software (Terry et al., 2020).

Talik Refreeze Modeling
We modeled top-down talik refreezing and permafrost aggradation beneath DLBs using the Geophysical Institute Permafrost Laboratory version 2 (GIPL2) model (Jafarov et al., 2012;Marchenko et al., 2008). This model estimates the transient subsurface temperature based on one-dimensional heat conduction and phase change. It includes the thermal effect of snow and the subsurface is divided into layers with different porosities and thermal properties (frozen and thawed heat capacities and thermal conductivities). The sediments and sedimentary rocks are assumed to be fully saturated, therefore, the volumetric porosity is assumed to the same as the frozen or unfrozen volumetric water content. GIPL2 has been successfully applied to model permafrost temperature in Alaska (Arp et al., 2016;Jafarov et al., 2012;Nicolsky et al., 2009Nicolsky et al., , 2017Romanovsky et al., 2007).
We chose the Derksen and Imakruak DLBs at the Teshekpuk region for the talik refreeze modeling because we have evidence of when they drained based on historical aerial imagery , and the geophysical results suggest the presence of remnant taliks for both basins. Derksen basin drained partially between 1966 and 1972, therefore the shallow water temperature was considered in the model as the upper boundary condition. Imakruak basin drained partially between 1955 and 1966, and then completely drained in 1999. In this case, we considered surface temperatures to be influenced by shallow water between 1960 and 1999 and air temperature and snow on the ground as a climatic forcing after 1999.
Daily mean air temperature and snow depth model inputs were based on publicly available (http://climate. gi.alaska.edu/acis_data) historical data from Utqiaġvik (formerly known as Barrow), located about 120 km northwest of the Teshekpuk region. For partial drainage, we estimated the shallow water temperature using a linear regression based on the observed water temperature measurements for Derksen Basin for the years 2007-2009 (Arp et al., 2016) and the daily mean air temperature records. We assumed an initial talik thickness of 70 m with a temperature of 4 °C from the surface to 50 m depth and 0 °C between 50 and 70 m depth, and −2 °C at 80 m depth. For the bottom boundary condition, we assumed a constant temperature gradient of 1 °C per 100 m, which is consistent with the permafrost thickness ∼320 m estimated at Drew Point borehole (Clow, 2014; Figure 1). Because we do not know exactly when the lake drained, our simulations started on July first, which is a reasonable time for a lake drainage to occur on the Arctic Coastal Plain of Alaska due to the timing of snow and lake ice melt (Jones & Arp, 2015). Moreover, we assumed thermal equilibrium in the ground temperature regime right before lake drainage, which is a reasonable assumption considering that the lakes would have formed thousands of years before drainage (Bockheim et al., 2004;Hinkel et al., 2003).
The subsurface properties were defined based on the lithologic description from the Drew Point borehole (Gryc, 1988) and field observations made from shallow (<1 m) permafrost cores. We assumed a surficial layer of lake sediments with 3 m thickness and porosity of 0.50. We also assumed a porosity of 0.20 for the sandstone (below 24 m depth) based on the values reported for sandstone of the Colville Group on the Umiat region on the North Slope of Alaska (Fox et al., 1979). Because the unconsolidated sediments (between 3 and 24 m depth) are formed of varied sediment textures, including sand, silt, and clay with different proportions within this layer, we performed a sensitivity analysis by varying the porosity of the unconsolidated sediments layer between 0.30 and 0.50 (Manger, 1963). We also performed a sensitivity analysis by varying the thermal properties according to the range of values reported in existing literature (Nicolsky et al., 2017;Roy-Léveillée & Burn, 2017): thawed thermal conductivity (K t ) between 0.9 and 1.6 W m −1 °C −1 , frozen thermal conductivity (K f ) between 1.5 and 2.4 W m −1 °C −1 , thawed heat capacity (C t ) between 3.10 × 10 6 and 3.27 × 10 6 J m −3 °C −1 , and frozen heat capacity (C f ) between 1.50 × 10 6 and 2.17 × 10 6 J m −3 °C −1 . The lake sediments, unconsolidated sediments, and sandstone layers were assumed to have the same thermal properties. Figure 2. Three categories of surface nuclear magnetic resonance (NMR) signal detection: (a) very low signal-to-noise ratio (SNR) (3:1) with limited or no evidence of exponential relaxation; (b) low SNR (10:1) with evidence of exponential relaxation between 0.5 and 2.0 A·s on the pulse moment axis and 0 and 50 ms in time; and (c) high SNR (30:1) with evidence of exponential relaxation. For each category, the left panel shows the stacked signal amplitude (nV) as a function of pulse moment (A·s), which is a proxy for depth, the right top panel is the average signal amplitude (nV) as a function of time (ms), which reflects the total water content, and the right bottom panel is the signal frequency spectrum. FFT, fast Fourier transform.

Geophysical Surveys
The surface NMR data quality was divided into three categories of signal detection ( Figure 2): (1) very low SNR (∼3:1) with limited or no evidence of NMR relaxation and no peak on the Larmor frequency of water in the spectrum (Figure 2a), therefore, we interpret this as no liquid water detected; (2) low SNR (∼10:1) with evidence of exponential relaxation and an identifiable low peak at the Larmor frequency of water (Figure 2b), indicating the presence of low, but real liquid VWC, and insufficient information to quantitatively interpret VWC distribution with depth; and (3) high SNR (∼30:1) with strong evidence of exponential relaxation and clear peak on Larmor frequency of water (Figure 2c), indicating ample subsurface liquid water and sufficient information to interpret VWC distribution with depth. It is worth noting again that surface NMR may not detect water held in very small pore spaces or partially unfrozen pore water, which may lead to very short relaxation times. RANGEL ET AL.    (Figure 3d). The left panels show the NMR signal time-series and that a higher average signal amplitude leads to a higher SNR since our study sites present low noise levels. The right panels show the signal frequency spectrum, where a peak on the Larmor frequency is observed when the SNR is at least 10:1 (Figures 3a, 3c, and 3d), and no peak when the SNR is very low (Figure 3b). This comparison shows that in general, a shorter t P can lead to a higher SNR. Longer t P can reach deeper DOI, however, it can result in lower SNR due to shorter relaxation times being obscured by the transmitting pulse.
We interpret the presence of remnant taliks based on either unfrozen water content estimated with surface NMR or TEM resistivities in DLBs compared to measurements on primary surface sites and borehole resistivity logs. For each study region, we have borehole lithologic descriptions (Figures 4a and 5a) (Gryc, 1988) along with TEM resistivities; (c, d) Peatball Bluff; (e) Bean Basin; (f) Schmutz Basin; (g, h) Imakruak Basin; and (i, j) Derksen Basin. US = unconsolidated sediments, and SS = sandstone of the Colville Group. DOI, depth of investigation; SE, single-exponential relaxation; MEA, multiexponential relaxation after excitation pulse; MED, multiexponential relaxation during and after excitation pulse. The gray shading zone on the resistivity profiles represents the ± standard deviations and the vertical line at 10 Ωm [log 10 (10) = 1] represents the interpreted minimum resistivity for permafrost. The horizontal lines (e-j) represent the interpreted depth to the top of talik. the interpretations. We assume that borehole sites (Figures 4b and 5b) and primary surface sites (Figures 4c  and 5c) represent permafrost with zero or very low liquid water content, and interpret that permafrost resistivity is >10 Ωm. At the Drew Point borehole, the TEM resistivity sounding is compared with the resistivity log ( Figure 4b) and both datasets reveal consistent results >10 Ωm.
At Peatball Bluff, even though we observed a thin zone with resistivities <10 Ωm between 5 and 10 m depth (2 and −3 m elevation, Figure 4c  At Inigok Camp, the resistivities were >10 Ωm (Figure 5c), suggesting only permafrost down to the depth of investigation (>100 m). Additionally, even though the surface NMR result (Figure 5d) presents a maximum VWC of 0.055 m 3 m −3 , it shows an irregular distribution with depth and signal category 1. Therefore, we suspect that these VWC values are spurious due to the very low SNR and no peak on the Larmor frequency of water. At Blustery Basin, only a surface NMR measurement was conducted, and it has a signal category 1. Although we observed a maximum VWC of 0.035 m 3 m −3 (Figure 5e), we interpret that no liquid water was detected, suggesting that the remnant talik beneath Blustery is possibly completely refrozen. For Deep Basin, we observed resistivity <10 Ωm (Figure 5f) below ∼25 m depth (20 m elevation), suggesting the possibility of a remnant talik. However, this basin has an NMR signal category 1 (Figure 5g), suggesting that no liquid water was detected. We suspect that the pore water is partially frozen, which generated a signal that is below the limit of the surface NMR detection. Perhaps only TEM was able to distinguish the partially unfrozen pore water. For Three Creatures Basin, we interpret a remnant talik from ∼22 m depth (26 m elevation) based on observed resistivity <10 Ωm ( Figure 5h) and VWC of 0.118 m 3 m −3 (Figure 5i) corresponding to signal category 3. Lastly, at Pingo Basin (Figure 5j), only a TEM measurement was conducted and we observed resistivities <10 Ωm below ∼40 m depth (3 m elevation), suggesting a remnant talik at this depth.
For Derksen (Figures 4i and 4j) and Three Creatures (Figures 5h and 5i) basins, the surface NMR results present a significant VWC shallower than the transition where the resistivity is <10 Ωm, which suggests that the depth to the top of talik may be shallower than interpreted. However, since both surface NMR and TEM methods present some uncertainty in depth, our interpretation is based on both results and we assumed the depth to the top of talik where the resistivity transition was observed. On the other hand, it is also possible that, due to differences in porosity and salinity in Derksen and Three Creatures basins, they may have a higher resistivity threshold for the frozen-unfrozen transition compared to the primary surface sites (Peatball Bluff and Inigok Camp) and the borehole sites (Drew Point and Inigok North). If we had TEM and surface NMR data for more DLBs, we could potentially calibrate a new resistivity threshold.
For Inigok Camp ( Figure 5d) and Blustery Basin (Figure 5e), although the surface NMR results present a significant VWC (0.055 m 3 m −3 and 0.035 m 3 m −3 , respectively), we presume that those VWC values are not real because both sites have a signal category 1, i.e., a very low SNR and no peak on the Larmor frequency of water. Moreover, the inversion artifacts observed for Inigok Camp (Figure 5d) and Blustery Basin (Figure 5e) results can be attributed to some instrument noise during the data acquisition. Deep Basin also presents a surface NMR signal with category 1 (Figure 4g) and no water was detected, therefore, our interpretation for this basin relied on the TEM results only (Figure 4f). Finally, although Imakruak Basin result presents a very low VWC (0.015 m 3 m −3 , Figure 4g), we know that this water signal is real because it has a visible peak on the Larmor frequency (Figures 2b and 3a). Table 3 includes a summary of the DLBs geophysical results, including TEM average linear resistivity (  ) until ∼80 m depth, surface NMR VWC MAX , and signal category. Table 3 also includes the interpreted depth to the top of taliks (Figures 4 and 5) and estimations of top-down talik refreezing rates, which are calculated based on the depth to the top of talik (m) and the years between the drainage period inferred from historical aerial imagery (Table 1) and the geophysical data acquisition (Table 2). It is possible to observe that, in general, a complete lake drainage results in a faster refreezing rate.  to depths where water was detected. By comparing to the lithologic profiles (Figures 4a and 5a), the top light gray zone corresponds to the unconsolidated sediments layer, which has a maximum VWC of 0.098 and 0.118 m 3 m −3 , and average  2 T of 68 and 60 ms for Derksen and Three Creatures basins, respectively. The lower darker gray zone probably corresponds to the sandstone layer, which has a maximum VWC of 0.009 and 0.035 m 3 m −3 , and average  2 T of 64 and 79 ms for Derksen and Three Creatures basins, respectively. Compared to the sandstone, the unconsolidated sediments have a higher VWC and a lower difference on the average  2 T between the two basins. The sandstone has a smaller average pore size than the unconsolidated sediments, therefore, we would expect a lower  2 T for the sandstone. However, the sandstone average  2 T is slightly higher than the unconsolidated sediments. This can be attributed to the low VWC of the sandstone, which leads to higher uncertainty on  2 T inversion results. For the three signal categories defined in Figure 2, we plotted surface NMR inverted VWC versus TEM inverted resistivities (Figure 7), limited to coincident observations within the surface NMR DOI (maximum 50 m depth). The plotted surface NMR VWC is the average value among the single-exponential (SE), multiexponential after excitation pulse (MEA), and multiexponential during and after excitation pulse (MED) inversion results (Figures 4 and 5) for a given depth and the uncertainty was estimated based on the minimum and maximum VWC. We used equation (2) to estimate VWC with the TEM resistivities, where  is the measured TEM resistivities and m was assumed as 1.3, which is a typical value for unconsolidated sands (Archie, 1942). Relatively high pore water salinity (Collett & Bird, 1988) and low clay content (Black, 1964) were observed for late Quaternary sediments on the Arctic Coastal Plain of Alaska, therefore, the surface conductivity effect may be very low (Ruffet et al., 1995). Based on the surface conductivity for sandstones with low clay content between 13.6 × 10 −4 and 32.6 × 10 −4 S m −1 (or 735.3 and 306.7 Ωm) reported by Walker et al. (2014) and Glover (2016), we assumed  s of 700 Ωm. Due to the high salinity, the pore water controls the bulk resistivity, and it can be up to ten times more conductive than the bulk resistivity (Glover et al., 1994). Therefore, because the lowest measured bulk resistivity is ∼3 Ωm, we assumed  w of 0.3 Ωm. The corresponding root mean square errors (RMSE) are shown in Figure 7. Despite our limited information to determine m and  s , these plots show a weak relationship between the inverted VWC and resistivities.
This may be explained in part due to the difficulty of calibrating an accurate electrical transform model that applies to pore spaces that may be filled both with solid (ice) and brine liquid phase. The low VWC associated with resistivities <10 Ωm may be related to water in very small pore spaces or partially unfrozen pore water. Ice formation expels dissolved solids and, consequently, the remaining unfrozen pore water can present higher salinity, which results in lower resistivities. As previously mentioned, water in very small RANGEL ET AL.  Figure 2; 2 Drainage degree at the location where the geophysical data was acquired; 3 Approximate depth to the top of talik around the center of each DLB, interpreted based on the geophysical results (Figures 4 and 5); 4 Years between the drainage period (Table 1) and the geophysical data acquisition (Table 2); 5 Top-down talik refreezing rates, which is the 3 talik top depth (m) divided by the 4 drainage time (years). No talik was detected below Blustery Basin, which is interpreted as refrozen.

Table 3
Geophysical Results and Interpretation pore spaces or partially unfrozen pore water can be a limitation for surface NMR detection, which are characterized by relaxation times too short to be measured. Figure 8 presents the top-down thermal modeling results of permafrost aggradation and talik refreeze compared to the depth to the top of the talik interpreted from the geophysical data. We plotted the −2°C, −1°C, and 0 °C isotherms for unconsolidated sediments with an average VWC of 0.40 and average thermal properties values for all layers (K t = 1.25 W m −1 °C −1 , K f = 1.95 W m −1 °C −1 , C t = 3.185 × 10 6 J m −3 °C −1 , and C f = 1.835 × 10 6 J m −3 °C −1 ). The gray shaded zones represent the possible range of the depth to the top of the talik (0 °C isotherm) by varying the VWC of the unconsolidated sediments between 0.30 and 0.50. The darker gray zone was estimated using average thermal properties values, and the lighter gray zone was estimated by varying the thermal properties (K t = 0.9-1.6 W m −1 °C −1 , K f = 1.5-2.4 W m −1 °C −1 , C t = 3.10 × 10 6 to 3.27 × 10 6 J m −3 °C −1 , and C f = 1.50 × 10 6 to 2.17 × 10 6 J m −3 °C −1 ). For Derksen Basin (Figure 4a), the results do not change within the discretization of the model (1 m layer thickness between 5 and 16 m depth) when varying the thermal properties values. A higher VWC results in a shallower refreezing depth due to the thermal properties of frozen or unfrozen water (Romanovsky & Osterkamp, 2000). In terms of thermal properties, a combination of lower thermal conductivity and higher heat capacity results in shallower refreezing depth.

Thermal Modeling
For Derksen Basin (Figure 8a), the geophysical results suggest the depth to the top of the talik was 8.4-9.6 m in 2016, which is within the range of possible depths of 8-11 m obtained from modeling. For Imakruak Basin (Figure 8b), the geophysical results suggest the depth to the top of the talik was 19-25 m in 2016, and 20-38 m for the model results. It is possible to observe in Figure 8b that the talik refroze at a faster rate after complete drainage occurred in 1999. The model and geophysical results are in agreement and provide two lines of evidence of a remnant talik for Derksen and Imakruak basins.
For Imakruak, the −2 and −1 °C isotherms reach deeper depths and have a similar trend as the 0 °C isotherm. On the other hand, for Derksen, the −2 and −1 °C isotherms are shallower and were influenced by the warmer conditions between around 1980 and 2005.

Drainage Age and Remnant Talik
Based on historical aerial imagery , we have evidence that DLBs in the Teshekpuk region drained after 1949, which is the year of the oldest available imagery, and the DLBs in the Inigok region drained before 1949. The results suggest the presence of a remnant talik for all the studied DLBs in the Te-RANGEL ET AL.

10.1029/2020JB020889
14 of 21 Error bars indicate the ± mean uncertainty for each method. The gray curves represent the VWC calculated with the modified Archie equation (2) (Ruffet et al., 1995). RMSE, root mean square error.
shekpuk region (Derksen, Imakruak, Bean, and Schmutz), which is consistent with their recent drainage. For the DLBs in the Inigok region, the results suggest the presence of a talik beneath Three Creatures, Deep, and Pingo basins, and no talik was detected at Blustery Basin. Based on the hypothesized talik refreezing time on the order of decades (Ling & Zhang, 2004;Mackay, 1992), our geophysical results suggest that: (i) Three Creatures Basin probably drained shortly before 1949 because the surface NMR result showed a relatively high VWC, similar to Derksen Basin; (ii) Deep Basin drained somewhat before Three Creatures Basin since our results suggest a deeper top of remnant talik and low VWC was detected with surface NMR (Table 3); (iii) Pingo Basin probably drained before Deep Basin because the results suggest a deeper top of remnant talik; and (iv) Blustery is the oldest basin and probably drained hundreds of years before present since our results suggest that its talik is probably completely refrozen. Three Creatures and Deep basins drained partially and then completely where the geophysical data were acquired, and Pingo and Blustery basins drained completely ( Table 3). The drainage degree plays a role on the top of talik refreezing depth and can be a confounding factor in the interpretation of the drainage timing. Other surface characteristics can help to constrain it, such as icewedge formation and vegetation succession, and also radiocarbon dating of the lowermost organics in the postdrainage sediments (Hinkel et al., 2003).

Geophysical Detection of Taliks
Based on borehole resistivity logs and TEM results on primary surface sites, we interpret that permafrost resistivity is >10 Ωm in our study sites. This interpretation is consistent with the minimum permafrost resistivity observed by Collet and Bird (1993) based on 156 borehole logs from the Prudhoe Bay and Kuparuk River oil fields, and by Overduin et al. (2012) based on electrical resistivity measurements around Utqiaġvik, all in the Arctic Coastal Plain of Alaska. Compared to surface NMR and other geophysical methods (e.g., GPR and ERT), TEM soundings can achieve deeper DOI (>100 m in our case). Minsley et al. (2012) used airborne frequency domain EM to image taliks in discontinuous permafrost at the Yukon Flats, Alaska, reaching depths of ∼100 m. Their results show higher resistivities compared to our results due to differences in equipment and lithology. Creighton et al. (2018) used TEM surveys to determine the talik thickness below Peatball Lake in the Teshekpuk region, suggesting that during the 1,400 years of the lake lifespan (based on radiocarbon dating), Peatball Lakes's talik reached ∼91 m depth. They also investigated the sensitivity of TEM to measure talik thickness based on forward modeling and inversion of synthetic data, showing that TEM can detect a talik thickness up to 130 m depth at Peatball Lake. Our results do not show clear evidence of the deepest extent of the remnant talik. However, the TEM resistivities increase to values close to 10 Ωm below ∼65 m depth for Derksen (−64 m elevation), Bean (−63 m elevation), and Pingo basins (−22 m elevation), and >10 Ωm ∼ 60 m depth (−57 m elevation) for Schmutz Basin, which may be related to the depth of the base of the talik and permafrost interface. The resistivity profile for Schmutz Basin (Figure 4f) present values close to 10 Ωm between −9 and −18 m elevation, which may be related to partially unfrozen pore water. However, we are cautious about interpreting it as evidence for a closed talik because it is right at the resistivity threshold for permafrost.
Our surface NMR measurements have a shallow DOI (<50 m) due to the low subsurface resistivity and the high geomagnetic field inclination on our study sites (>80°) (Parsekian et al., 2019). The surface NMR result  (Figure 4a) with average volumetric water content (VWC) of 0.40 and average thermal properties values (K t = 1.25 W m −1 °C −1 , K f = 1.95 W m −1 °C −1 , C t = 3.185 × 10 6 J m −3 °C −1 , and C f = 1.835 × 10 6 J m −3 °C −1 ). The gray shaded zones represent the range of possible 0 °C isotherms assuming unconsolidated sediments with VWC between 0.30 and 0.50. The darker gray zone was estimated using average thermal properties values, and the lighter gray zone was estimated by varying the thermal properties (K t = 0.9 to 0.6 W m −1 °C −1 , K f = 1.5-2.4 W m −1 °C −1 , C t = 3 × 10 × 10 6 to 3.27 × 10 6 J m −3 °C −1 , and C f = 1.50 × 10 6 to 2.17 × 10 6 J m −3 °C −1 ). For Derksen Basin (a), the results do not change within the precision of the model when varying the thermal properties values. The shading was included only for 0 °C isotherm to preserve the clarity of the figure. The black dots are the depth to the top of the talik based on the geophysical results and the error bars indicate the ± mean uncertainty estimated with surface nuclear magnetic resonance (NMR) inverted resolution matrices.
at Blustery Basin suggests that the remnant talik is possibly completely refrozen. A TEM survey at this site would help to confirm whether this is the case. The results for Imakruak (Figure 4h), Derksen (Figure 4j), and Three Creatures (Figure 5i) are similar to the expected result simulated for a closed talik by Parsekian et al. (2013), which corroborates our interpretation of a remnant talik for those DLBs. The results for the primary surface sites, Peatball Bluff ( Figure 4d) and Inigok Camp (Figure 5d), and Blustery Basin (Figure 5e) are similar to their simulated result for permafrost not impacted by lake development at the surface. Terry et al. (2020) used surface NMR to detect a talik with VWC ∼0.10 m 3 m −3 beneath the Kuparuk River aufeis on the North Slope of Alaska, similar to our observed values for the remnant taliks below Derksen and Three Creatures basins, and they also observed VWC < 0.03 m 3 m −3 for permafrost, which is similar to our primary surface sites and Blustery Basin results.
The analysis of our surface NMR dataset shows its advantages and limitations to detect unfrozen water found in remnant taliks beneath DLBs. The relaxation pulse duration (  2 T ) for the category 3 results (Figure 6) show that the unconsolidated sediments present a consistent average  2 T values (60-68 ms) due to its relatively high VWC (0.098-0.118 m 3 m −3 ), while the sandstone layer demonstrate a wider range of  2 T values (64-79 ms) and low VWC (0.009-0.035 m 3 m −3 ). The sandstone layer may have smaller pore sizes than the unconsolidated sediments, which can result in short relaxation times likely below the surface NMR detection limit. However, our sandstone  2 T values are consistent with previously observed ranges between 64 ms (Hein et al., 2017) and 80 ms (Legchenko et al., 2002). Yet, the observed resistivities of <10 Ωm for the sandstone layer support the interpretation that its pore water is probably partially or completely unfrozen. The plotted surface NMR VWC versus TEM resistivities compared to the calculated VWC using the Archie Equation (Figure 7), highlight the unambiguous ability of surface NMR to detect unfrozen water when the SNR is sufficiently high (>10:1). The high RMSE values (0.070-0.082 m 3 m −3 ) can be explained because the VWC was low (0.009-0.035 m 3 m −3 ) in the sandstone layer (below 25-30 m depth), but the TEM resistivities were also low (3-6 Ωm), leading to a high VWC (>0.100 m 3 m −3 ) calculated with the equation (2). Moreover, the parameter m can be sufficiently different for the unconsolidated sediments and sandstone layers (Archie, 1942), and our generalized values do not fit both well.
In general, geophysical methods, including electrical, electromagnetic, and seismic, present limitations to distinguish between frozen and unfrozen materials in permafrost environments (Yoshikawa et al., 2006), mainly when pore water salinity plays an important role (Ross et al., 2007). Our results show that surface NMR and TEM are a valuable combination to detect remnant taliks below DLBs. Both methods enabled relatively rapid measurements (<1 h each) even under Arctic climatic conditions.

Geophysical Properties of Partially Refrozen Materials
In general, the permafrost in our study sites has low resistivity due to the relatively high concentration of salt deposited during marine transgressive events through the geologic history of the Arctic Coastal Plain of Alaska (Brigham-Grette & Hopkins, 1995;Dinter et al., 1990;Gryc, 1988). Salinity and pressure can depress the freezing point of water, and permafrost contains partially unfrozen pore water at temperatures below 0 °C (Gryc, 1988;Kleinberg & Griffin 2005;Keating et al., 2018;Mackay, 1997), which enables electrolytic conduction (Ross et al., 2007;Yoshikawa et al., 2006). The temperature logs from Drew Point and Inigok North boreholes showed temperatures <−5 °C down to 100 m depth (Clow, 2014). Collett and Bird (1988) observed a maximum pore water salinity of 19 ppt and pressure of 22 kPa/m down to 1,500 m depth in boreholes at the Prudhoe Bay and Kuparuk River oil fields and suggested that the freezing point depression is no more than ∼2 °C. Thus, the freezing point depression due to salinity or pore pressure may be negligible for subsurface material with temperatures <−2 °C. However, it can be an important factor for a remnant talik that is still undergoing refreezing with temperatures >−2 °C and, consequently, may contain partially unfrozen pore water. Figure 8a shows how dynamic the −2 and −1 °C isotherms can be during the talik refreezing process after lake drainage. Furthermore, Figure 8b shows that the depth to the top of a talik can be shallower if the freezing point depression is 2 °C. Recent studies based on laboratory measurements of the resistivity response during freezing of porous media showed that above 0 °C, the bulk resistivity is controlled by the pore water resistivity which increases with decreasing temperature, and below 0 °C, the bulk resistivity depends on the unfrozen water content and its resistivity (Herring et al., 2019;Ming et al., 2020). This evidence suggests that TEM measurements can detect partially unfrozen saline pore water in remnant taliks, since electrical current can flow through very small pore spaces. Kleinberg & Griffin (2005) conducted NMR measurements in cores up to 428 m depth from the Arctic Coastal Plain of Alaska, to estimate the unfrozen water content in sediments with temperatures below 0 °C. Their results show that sediments with total porosity of ∼0.37 can have liquid water-filled porosity of 0.04-0.14 at −3 °C and 2 T <5 ms. Because the surface NMR instrument dead time is 5.5 ms and  2 T ≤ 2 T (Behroozmand et al., 2015), the partially unfrozen pore water signal may relax too fast to be detected by surface NMR, which is consistent with our observations.

Thermal Modeling
The top-down thermal modeling of permafrost aggradation and talik refreeze provides a second line of evidence, complementary to our geophysical observational data, that a remnant talik is present below Derksen and Imakruak basins ( Figure 8). Historical aerial imagery shows that Imakruak drained a decade earlier than Derksen, and the results suggest that the top of the talik at Imakruak is about two times deeper (Figure 8). Considering that both basins are in the Teshekpuk region, i.e., they have a similar climate history and lithology, this difference can be attributed to the fact that Derksen drained partially and Imakruak drained completely in 1999. To enable an objective comparison, the modeling was done independently of the geophysical results. Even considering the best information available in our thermal model, it still depends on some assumptions such as: (1) estimates of the volumetric water content, (2) constant volumetric water content through time, and (3) similar climate record as at Utqiaġvik, located ∼120 km northwest. When all other variables remain constant, the modeling results show that a lower volumetric water content results in a deeper talik due to the influence of frozen or unfrozen water content on thermal conductivity and latent heat effects (Romanovsky & Osterkamp, 2000).
Our modeling results suggest that after 40 years since lake drainage, the depth to the top of the talik is 8-11 m for Derksen Basin, and 16-26 m for Imakruak. Ling and Zhang's (2004) thermal modeling of talik refreeze after a complete lake drainage on the Arctic Coastal Plain of Alaska suggest a top-down refreeze of 23-25 m in 40 years, which is consistent with our results considering that Derksen Basin drained partially and Imakruak Basin drained completely only after 30 years of the simulation period.
Even though the TEM resistivity and surface NMR VWC measurements are sensitive to unfrozen VWC, not necessarily to temperature, and a significant amount of pore water can be partially unfrozen between −2 and 0 °C (Romanovsky & Osterkamp, 2000), our model results still agree very well with the geophysical results. For Imakruak (Figure 8b), the −2 °C isotherm is around the same depth as the resistivity is <10 Ωm (Figure 4g) and the detected VWC is very low for this basin (Figure 4h), indicating that the pore water is mostly frozen. For Derksen Basin (Figure 8a), the −2 °C isotherm was at ∼4.5 m depth in 2016, therefore, it is possible that there is some unfrozen VWC between 4.5 and 9 m depth that was not detected by the geophysical measurements.
Our observations show that a remnant talik refreezing below young DLBs are easier to interpret and model, but it is more challenging in old DLBs, suggesting more complexity than predicted by numerical simulations (Ling & Zhang, 2004) because some processes are very difficult to account for in the thermal modeling. For example, talik evolution depends on the lake ice regime because a talik develops when the lake water is deeper than the winter ice thickness (floating ice) and, throughout its lifespan, a lake can have different ice regimes depending on the climate conditions (Arp et al., 2016). Pore water salinity is also a difficult parameter to be considered and it is an important component of the permafrost on the Arctic coastal plains that influences permafrost degradation and aggradation (MacKay, 1997;Mackay & Burn 2002a, 2002b. Future modeling efforts of permafrost dynamics should consider incorporating these components and associated geophysical constraints.

Implications of Permafrost Aggradation Dynamics After Lake Drainage
Our geophysical observations and thermal modeling results may provide useful information to predict future evolution of permafrost aggradation beneath DLBs across the Arctic. Lake drainage in the Arctic has a significant influence on permafrost dynamics (Hopkins, 1949) and developing techniques to better understand talik refreeze rates is critical to our understanding of Arctic landscape dynamics. For example, if lake drainage frequency increases in the region, it could result in more active permafrost aggradation, which would likely decrease the emissions of permafrost carbon to the atmosphere ( van Huissteden et al., 2011). On the other hand, if lake drainage frequency decreases, lake talik development could increase and probably lead to an increase in permafrost carbon emissions (Walter et al., 2007). Lake drainage frequency is already increasing in the boundary between discontinuous and continuous permafrost zones (Lantz & Turner, 2015;Nitze et al., 2018), and may increase in continuous permafrost zone in the future . Future warmer and wetter conditions in the Arctic (Bintanja & Selten, 2014) may stimulate the mechanisms leading to lake drainage, including lake expansion , bank overtopping (Jones & Arp 2015), and ice-wedge degradation (Liljedahl et al., 2016), and also increase the impacts on the water availability , tundra habitat, vegetation productivity (Lara et al., 2018), and snowdam outburst floods .

Conclusions
Geophysical methods provide a nondestructive and indirect way to investigate taliks in permafrost regions. This study provides new observations of permafrost aggradation below DLBs in the continuous permafrost zone of Alaska by combining geophysical measurements and thermal modeling to independently determine the depth to the top of remnant taliks. Our results show evidence of remnant taliks below several DLBs that have drained before and after 1949 (oldest aerial imagery). Partially unfrozen pore water represents a limitation for surface NMR detection due to the very short water signal and only TEM was able to detect thaw in small pores because electrical current can flow as long as there is continuity in the pore-filling fluid. Both geophysical and thermal modeling results showed a talik refreezing rate consistent with previous modeling results on the western Arctic Coastal Plain of Alaska. Complete lake drainage results in faster permafrost aggradation rates than a partial drainage due to the thermal properties of remnant lake water in the basin. Our findings may help to constrain talik refreezing models and provide valuable information to predict future impacts of climate warming on permafrost aggradation after lake drainage and subsequently both water and carbon cycling in the Arctic.