Study area

In 2008–2014, we conducted four marine expeditions and four drilling campaigns in the study area between 70–74° N and 129–131° E with focus on the near-shore Laptev Sea southeast of the Lena Delta, the Buor-Khaya Bay (BKB, between 70–74° N and 129–131° E), and the Dmitry Laptev Strait (DLS, between 72.5–73.5° N and 138–143° E, Fig. 1). In drilling campaigns, we investigated the thermal regime, geomorphology, lithology and geocryology of sediment cores extracted from drilled boreholes and sediments sampled along the drilling transect (Fig. 2). We also performed few geoelectrical surveys, results of which were validated by direct measurement of electrical resistivity of recovered sediments. In marine expeditions, we collected conductivity-temperature-depth (CTD) data, performed high-resolution sub-bottom profiling, sonar-derived imagery and visual observations (using an autonomous underwater vehicle) of geomorphological features of the seafloor (subsea thermokarst, ice scouring and pockmarks) associated with gas releases.

Figure 1: Study area bathymetry and position of the investigated sites. (a) Red and blue rectangles mark study areas, where drilling was conducted in 2011–2014, position of the sites investigated in marine expeditions, data from which are presented in Figs 4, 5, 6, 7, 8, 9 are shown as black triangles (2D sites) and red lines (transects); two black crosses in the blue rectangle show position of the drilling transect conducted in 2012–2014 (shown enlarged in b,c); (b) position of the boreholes drilled in March 2011–2013; (c) enlarged position of the drilling transect performed the northern part of MI in 2012–2014. Full size image

Figure 2: Geomorphological structure of the sediments along the drilled transect in the near-shore zone of the ESAS. (a) Cross section from Cape North, MI towards Cape Muostakh, Bykovsky Peninsula with marked positions of the re-drilled boreholes. Boreholes 4D-14 (former 301), 4D-13 (former 303), 3D-14 (former 304) and 2D-13 (former 305) are re-drilled boreholes; borehole 4D-12 is the borehole drilled for the first time. Red spikes identify changes in the position of the IBPT in drilled boreholes during the last 31–32 years; purple spike in borehole 4D-12 show change in the position of the IBPT since the time of inundation. The top of the red/purple spike marks the former top of the IBPT; the low end of the spike marks the recently observed IBPT; numbers below each borehole show the length of recovered sediment cores. (a) 1—marks of drilled boreholes; 2—marks of sediment sampling sites; 3—plant remains within sediments; 4—gravel; 5—pebbles; 6—frozen ground in the coastal IC; 7—change in position of the IBPT during the last 31–32 years (boreholes 4D-14, 4D-13, 3D-14 and 2D-13); 8—change in position of the IBPT in borehole 4D-12 (since inundation); 9—wedge-like structure of ice in the coastal IC; 10—medium sand; 11—fine sand; 12—silty sand; 13—sandy silt; 14—clayey silt; 15—silty clay; 16—clay; (b) shelf-ward dynamics in grain-size distribution in sediments collected at 5 cm depth at selected sites along the transect (weight-% of major size fractions): 17—change in contribution of sands (>0.1 mm); 18—change in contribution of silt (0.1–0.01 mm); 19—change in contribution of clay (<0.01 mm); shift in sedimentology from predominant contribution of large-medium grain-size fractions (20) to medium-fine grain-size fractions (21). Full size image

To assess modern rates of downward permafrost degradation, we re-drilled four boreholes first drilled in 1982–83 (ref. 20) and also drilled one new borehole northwest of Muostakh Island (MI); MI represents remains of the Bykovsky Peninsula submerged <0.5 kyr ago11. Water depths vary from 2.5 to 3.4 m. Bottom sediments are predominantly silty sands. The seafloor here represents an abrasion-accumulative terrace of the former coastal plain formed via inundation about 2.7–3 kyr ago21. MI IC thickness is 31 m, 10 m of which extends below sea level21. At the time of deposition, the IC temperature varied from −25 to −28 °C; the current MI permafrost temperature is −10.4 °C (ref. 23). Sand and silt are the predominant sediment types observed in the area; ice-poor sands of Pliocene-early Pleistocene age underlay the IC deposits20. Mechanical denudation processes dominate lithogenesis in this area; ice transport, bottom erosion, and sedimentary matter resuspension and redistribution are widespread34. This area has been strongly affected by global warming17,23, and exhibits a high coastal erosion rate, which has approximately doubled during the last 62 years23,35. During the last 5 kyr, after the glacio-eustatic sea level reached its highstand at the ∼30 m isobaths in the studied area, the shallower near-shore zone submergence has been occurring primarily via erosional processes36,37.

Lithology and thermal regime of sediment cores

Five boreholes (4D-14, 4D-13, 3D-14, 2D-13 and 4D-12) were drilled along the transect located 300–2,900 m from Cape North, MI and extending towards Cape Muostakh of the Bykovsky Peninsula (Fig. 2a). Borehole 4D-14 was drilled 300 m from Cape North, MI; a 9.5 m sediment core was recovered. An IBPT was identified at 8.6 m below sea level (b.s.l.). The lithological structure of the thawed sediments in borehole 4D-14 was uniform (Supplementary Fig. 1) and consisted of light grey/grey coarse sand with grain-size fraction >0.1 mm composing 99% of sediments (Fig. 2b). The cryostructure of the frozen part of the 4D-14 sediment core was massive, which could represent the ice wedge of the submerged coastal IC (Supplementary Figs 2 and 3).

Borehole 4D-13 was drilled 600 m from Cape North, MI. The sediment core recovered from the borehole was 21 m long and was composed of five units (Supplementary Fig. 1). The uppermost 6 m is characterized by partly bedded brown coarse-grained silty sand, which is underlain by 1.8 m of sandy silt followed by 10.6 m of partly bedded medium-grained silty sand, which included plant and wood remains. Below, there was a 1.1 m layer of clayey silt with organic inclusions. The lowermost 1.5 m layer was composed of partly bedded silty sand containing plant and wood remains. The IBPT in borehole 4D-13 was identified at 11.4 m b.s.l.; the temperature of the unfrozen part of the sediment varied from −1 to 0.1 °C (Fig. 3a), and was interpreted as thawed/cryotic (thawed is unfrozen at temperatures >0 °C; cryotic is unfrozen at temperatures <0 °C). Cryology of the 4D-13 sediment core below the IBPT was characterized as mostly non-structured ice within frozen sediments, which in places included fine lenticular ice, predominantly within silty sand (Supplementary Figs 2 and 3).

Figure 3: Thermal regime of subsea permafrost based on results of direct observations. (a) Data show that along the re-drilled transect, sediment temperatures in boreholes 4D-13, 2D-13 and 4D-12 varied from −2 to +1 °C. (b) Examples of thermal state of sediment cores recovered from a few boreholes in 2011–2014; the results show that sediment temperatures varied from −3 to +1 °C, increasing downwards in some cores. Full size image

Borehole 3D-14 was drilled 850 m from Cape North, MI; the recovered sediment core was 17.6 m long. The lithological structure of the sediment core was heterogeneous (Supplementary Fig. 1). The uppermost 1 m of sediments consisted of coarse/medium grey sand; this layer was underlain with 3.8 m of bedded silty sand interlayered with two thin layers (∼0.2 m) of medium-grained sand. Below, there was a 1 m layer of sandy silt followed by 2.4 m of medium-grained sand; the next unit (3.4 m) consisted of partially bedded sandy silt with organic inclusions and inclusions of plant remains. Below, there was a 3.6 m layer of silty sand followed by a 2 m layer of clayey silt containing organic inclusions. The IBPT was identified at 12.8 m b.s.l. The cryostructure of the sediment core below the IBPT was identified as mostly consisting of inclined lenticular ice; 0.5–2 mm ice lenses were included in the sediments (Supplementary Figs 2 and 3).

Borehole 2D-13 was drilled 2,500 m from Cape North, MI; a 30.4 m-long sediment core was recovered from the borehole (Supplementary Fig. 1). The uppermost 4.5 m of the sediment core consisted of medium- to fine-grained silty sand with lenses of coarse-grained sand and gravel included in the sediments. Below was 10.3 m of dark grey, greenish or black silty clay, which included lenses of plant remains and organic inclusions. The next unit was composed of 2.5 m of medium- to fine-grained sand, dark grey in colour; this unit was followed by 11.5 m of clayey silt, dark grey or greenish, which included organic inclusions, plant detritus and wood remains. The lowermost layer consisted of 1.6 m of silty sand, dark grey in colour. The IBPT was identified at 19.3 m b.s.l.; the cryostructure below the IBPT was ice-rich and predominantly non-structured with rare inclusions of lenticular ice and small lenses <0.1 mm in size; between these two layers we identified a layer (0.5 m) of cryotic sediments. The temperature of unfrozen sediments varied from 0 to 1 °C; these sediments were identified as thawed sediments interlayered with cryotic sediments (Fig. 3a; Supplementary Figs 2 and 3).

Borehole 4D-12 was drilled 2,900 m from Cape North, MI; the sediment core recovered from the borehole was 28 m long (Supplementary Fig. 1). The lithology of the upper 20 m of this core consists of multiple layers of silty clay and clay interspersed with thin layers (<1 m) of clayey sand (at 9.5 m b.s.l.) and clayey silt (at 10.6 m b.s.l.). Below, we identified a 4.5 m layer of medium-grained sand followed by 0.6–0.9 m layers of silty clay, clay, sandy clay and clayey sand. All layers included plant and wood remains and gravel. The IBPT was identified at 26.4 m b.s.l.; the cryostructure of the sediment core below the IBPT was massive and, at some depths, was interlayered by non-structured ice. The temperature of unfrozen sediments varied from −0.5 to 0.7 °C, and sediments were interpreted as remaining in a cryotic/thawed state (Supplementary Figs 2 and 3). The thermal regime in the three re-drilled boreholes presented in Fig. 3a is compared with the thermal state of the few other boreholes drilled in the near-shore ESAS area in 2011–2014 (Fig. 3b; other results from these boreholes are not presented in this paper). Volumetric ice content of different cryostructures identified in the drilled boreholes is presented in Supplementary Table 1.

Rates of permafrost degradation

To elucidate modern rates of downward permafrost degradation, we compared the position of the IBPT in four re-drilled cores with the position of the IBPT in four boreholes first drilled in 1982–83 at the same locations. The boreholes were localized using a shore-based theodolite technique (see Methods). The IBPT positions observed in 1982–1983 at sites 301, 303, 304 and 305 were at 3.3–4.2 m, 5.8–7 m, 8.3–8.6 m and 16–16.8 m b.s.l., correspondingly (Table 1). In 2012–2013, the IBPT positions were identified at 8.6 m, 11.4 m, 12.8 m, and 19.3 m b.s.l. at sites 4D-14 (former 301), 4D-13 (former 303), 3D-14 (former 304), and 2D-13 (former 305), correspondingly. IBPT deepening during the last 31–32 years varied from 9.3 to 18.3 cm year−1 with a mean rate of 14±3.1 cm (mean±s.e.m.) per year during the last 31–32 years. Position of the IBPT between drilled boreholes was established based on results of geoelectric investigation performed along the drilling transect. Sediment temperature, lithology and geochemistry were used to convert recorded electrical resistivity to possible ice content of sediments along the drilled transect (Supplementary Fig. 4; Supplementary Table 2).

Table 1 Rates of ice-bonded permafrost table deepening in the study area. Full size table

To assess rates of permafrost degradation before 1982, we calculated the time since inundation by applying an empirical modelling approach based on modern bathymetry and observational time series of coastal erosion. Following Bauch et al.36, we assume that by 5 kyr ago, sea level in this region achieved its highstand at the 30 m isobath. After that, coastal zone inundation occurred due to the erosion-based processes of thermo-denudation and thermo-abrasion37, implying that between two subsequent isobaths, the speed at which water depth adjusted to existing sea level was proportional to coastal erosion rates (Supplementary Figs 5 and 6). The mean coastal erosion rate was set based on observational data collected at regular MI stations, where coastal erosion has been monitored over the last few decades11,35. These time series documented that from 1983 to 2013 the coastline position shifted inland by ∼180 m, establishing a mean annual coastal erosion rate of 6 m year−1. We calculated that IBPT deepening rates from the time of inundation until 1982–1983 varied from 3.9 to 8.5 cm year−1 with a mean rate of 5.7±2.8 cm year−1, less than half as much as IBPT deepening rates observed during the last 31–32 years (Table 1). Calculations of pre-1982 rates of permafrost degradation are more uncertain than those observed during the last three decades, and reflect the current state of knowledge on this topic (Supplementary Discussion).

We found in our study that the further we travelled from the coast, the deeper we identified the position of the IBPT; at the same time, rates of permafrost degradation during the last 30 years were higher in boreholes 4D-14, 4D-13 and 3D-14 located closer to the coast. Indeed, the position of the IBPT changed from 8.6 m b.s.l. in borehole 4D-14, which was 300 m from Cape North, MI to 19.3 m b.s.l. in borehole 2D-13, which is 2,500 m from the coast; rates of permafrost deepening decreased from 18.3 cm year−1 in borehole 4D-13 to 9.3 cm year−1 in borehole 2D-13. The position of the IBPT in the drilled boreholes often coincided with the boundary between coarse- to medium-grained sediments (sands-silts) and medium- to fine-grained sediments (silt-clay). A specific feature of the lithology of the boreholes drilled closer to the coast (4D-14, 4D-13 and 3D-14) is that the uppermost sediment layers (5 cm from the top) in these boreholes were composed of sands, which represent remains of ICs of Pleistocene age; thus, in these boreholes either no Holocene-age accumulations occurred or erosion-driven processes balanced them. In the uppermost layers, sands with grain size >0.1 mm composed 70 to 99% of the total sediment weight, while fractions of silt (0.1–0.01 mm) and clay (<0.01 mm) was <30%. This distribution of grain size was different from that in boreholes further from the coastline (2D-13 and 4D-12) in which the fraction of finer sediments with grain size <0.1 mm increased to 50–95% (Fig. 2b).

Specific features of seafloor morphology

We used high-resolution reflection seismic and sonar-derived imagery to identify specific seafloor morphology features that could reflect ongoing permafrost-degradation processes such as thermokarst. When on-land permafrost thaws, it leaves behind a regularly shaped hummocky landscape that resembles a polygonal pattern4,7. In the study area, we frequently observed such a pattern on the seafloor (Fig. 4). Because the pattern was not buried by Holocene sediments, we suggest that the seafloor appearance reflects modern thermokarst processes, which keep developing after permafrost submergence. However, a polygonal pattern characteristic of thermokarst, as observed in the near-shore area of the ESAS, could also represent a relic feature and remain unburied by sediments due to the predominantly erosional character of sedimentological processes in some near-shore areas of the ESAS.

Figure 4: Two types of seafloor polygonal structure observed in the study area. (a) A high-resolution side-scan sonar image showing a type of polygonal structure, which represents a regularly shaped depression surrounded by a raised shaft; (b) a high-resolution side-scan sonar image showing a type of polygonal structure characterized by a regular pattern consisting of raised oval-like shapes. Full size image

Another permafrost destabilization mechanism, specific to the shelf system, is ice scouring, which mechanically disintegrates the upper frozen and/or thawed sediment layers38,39. Ice scouring was observed as a ubiquitous morphological feature, not only of the shallow part of the ESAS but also offshore (Figs 1a and 5a–c). Ice scouring is evidenced by a long, linear, relatively-straight furrow a few tens of metres wide extending for many tens of kilometres40. In the ESAS, ice scouring penetrated as much as 10 m into the sediments, and where surface sediments were underlain with free gas, strong ebullition to the water column through the scours was observed (Fig. 5d). Ice scouring likely provides an important mechanism, not only for unroofing shallow gas accumulations (Fig. 6a), but also for allowing gases to avoid anaerobic oxidation due to removal/disintegration of the uppermost sediment layers, where anaerobic oxidation of CH 4 occurs within the sulphate-reduction zone.

Figure 5: Ice scouring and associated features observed in the ESAS. (a) A high-resolution side-scan sonar image of an ice scour observed in the seafloor; (b) a vertical profile of the ice scour as it occurs on the high-resolution sub-bottom profile image in one location. The depth of the depression created by ice scour reaches ∼4 m; (c) a high-resolution side-scan sonar image demonstrating wide occurrence of ice scours in the seafloor in the far offshore area of the ESAS; (d) a high-resolution sub-bottom profile image showing multiple vertical profiles of the ice scours and bubble plumes propagating from the seafloor in the area of investigation. Full size image

Figure 6: Specific morphological features of the seabed and the seafloor observed in the study area. (a) Three-dimensional coloured view of sub-bottom morphology based on interpretation of high-resolution sub-bottom profile data (vertical profile) combined with side-scan backscatter data; colour variations refer to the amplitude response of the seabed (maximum shown in red, minimum shown in blue). Morphological features within sediments with highest amplitudes are shown as bright spots (red-orange colours). As seen from the image, bubble releases to the water occur where ice sours unroof shallow gas accumulations in the seabed; b.s.f. refers to below seafloor; (b) a view of multi-beam sonar backscatter data used to detect typical pockmarks in the seafloor in the study area where submerged thaw lake was identified as shown in Fig. 7. Full size image

The fate of submerged taliks

The fate of gas-migration pathways such as thaw-lake taliks that formed before ESAS inundation was assessed by repeated geophysical surveys in the Dmitry Laptev Strait (DLS) specifically in regions of high dissolved CH 4 concentrations, called ‘hot spots’41 (here between 72.5–73.5° N and 138–142° E, Fig. 1a). Such hot spots suggest elevated sediment permeability for ascending gas, which could only be possible if deep taliks serve as gas-migration pathways within permafrost. Despite the fact that this area was submerged relatively recently, <<1 kyr ago, which is insufficient time for upwardly developing taliks to penetrate through the entire permafrost body and reach the permafrost top7, significant variations in position of the IBPT (>120 m) in the recovered sediment cores were observed during the drilling campaign performed in the DLS in the 1980s (ref. 5).

In our investigation of DLS, interpretation of high-resolution seismic imagery allowed identification of a presumed submerged thaw-lake basin beneath the seafloor with a central washed-out zone (Fig. 7). Repeated observations along the same transect 3 years later revealed a gas front moving upwards, because the top gas front boundary changed from ∼3 m below the seafloor in 2008 to the surface of the seafloor in 2011. Hydro-acoustical images of gas bubbles released from doming surface sediments were recorded, and increased concentrations of dissolved CH 4 in the water column (≤120 nM) were measured. We also observed visible morphological signs of recent gas releases (pockmarks) within the site (Fig. 6b). On the basis of our data, we suggest that submerged thaw-lake taliks may not freeze; instead, they may keep developing, creating pathways for ascending gas. Such pathways propagating throughout the entire permafrost body could allow vigorous CH 4 venting, if over-pressured gas reservoirs beneath are unroofed.

Figure 7: High-resolution sub-bottom profile images of gas plumes propagating through the centre of a submerged lake basin in 2008 and 2011. (a) Column-like acoustic anomalies (blanked areas) interpreted as gas plumes are shown as areas with no colour and blue colour. Acoustic anomalies highlighted in blue correspond with the centre position of the submerged lake basin and interpreted as gas propagating through the submerged thaw-lake talik. In 2008, the top boundary of this acoustic anomaly was located ∼5 m below the seafloor; (b) in 2011, the position of the top boundary of this acoustic anomaly was reaching the seafloor, causing doming of the surface layer of sediments. During the survey in 2011, bubble plumes of CH 4 releasing to the water column from the seafloor were observed at this site; multiple pockmarks in the seafloor were also documented within the site as shown in Fig. 6b. Full size image

Another possible mechanism for preventing taliks from freezing and/or causing talik formation could be groundwater flow through coastal sediments, especially in the areas underlain with the faults. This could cause formation of so-called tectono-genic taliks42. These taliks could be identified by groundwater discharge that could be manifested as large point sources, which are temporally and spatially variable and could have a significant impact on the geochemical parameters of coastal waters43. Groundwater is usually terrestrially derived and enriched in naturally occurring radionuclides such as, ex224Ra and 223Ra, relative to seawater44. Releases of groundwater, including intra-permafrost water, might lead to formation of taliks within subsea permafrost, as we observed in one near-shore area, where a geoelectric survey showed spots of low electrical resistivity within the study area characterized by high resistivity of the sediment (Supplementary Fig. 7). High concentrations of ex224Ra and 223Ra (19 and 0.65 dpm/100 l−1 versus <0.1 and <0.05 dpm/100 l−1 in normal seawater) and thermohaline anomalies (higher temperature and lower salinity of seawater) were detected in the bottom water, which suggests groundwater discharge at this location. The location of this site coincides with the position of the fault45, which could have resulted in the development of a tectono-genetic talik.

Gas movement through sediments

We then investigated another part of the study area where possible development of a deep talik was suggested by modelling results16. We conducted high-resolution seismic surveys along the same transect in two subsequent years. The top boundary of the observed acoustic anomaly moved upwards by ∼5 m in the course of just one year (2011–2012) (Fig. 8). It was thus concluded that the observed acoustic anomaly could not be ice within sediments; if it were, the position of the top ice boundary would have remained the same over the course of one year. In the marine environment, low-amplitude seismic anomalies, referred to as washed-out or semi-blanked zones, are associated with the presence of gas in sediments28,46; in permafrost, these anomalies may also result from physical property variations or changes associated with talik development29,47. It was shown that gas usually accumulates in the central portion of the channels presented by unconsolidated sediments held within consolidated sediments28. In permafrost areas, unfrozen sediments within taliks could represent such unconsolidated (unfrozen) sediments within consolidated (frozen) sediments. It is thus logical to expect that gas would accumulate and propagate in the central part of the talik, which is also what was observed at the reported site (Fig. 7). Finally, the velocity at which gas appears to be ascending is consistent with Darcy’s law (∼7 m year−1) assuming typical values for marine sediment permeability to gas, gas viscosity and gas density, suggesting that gas is propagating within the thawed/cryotic sediments. This also excludes the presence of ice-bonded permafrost in the sediments.

Figure 8: Upward movement of gas front observed in the same area of investigation in two subsequent years. (a) A high-resolution sub-bottom profile image is interpreted based on colour variations referring to change from maximum seismic amplitudes (shown in red) to minimum amplitudes (shown in blue). Blanked (or washed-out, shown in white) zones are interpreted as attenuated and scattered by gas bubbles (acoustic turbidity)28,46. The upper boundary of the blanked zone indicate position of the gas front in 2012, ∼2 m below the seafloor; red arrows indicate areas of most pronounced upward movement, ≤5 m during 1 year; (b) a high-resolution sub-bottom profile image shows that the upper boundary of the blanked zones in 2011 is located ≤10+ m below the seafloor, in the same areas marked with red arrows in a; blanked or semi-blanked acoustic anomalies are interpreted as gas-charging sediments based on validation by drilling (see Fig. 9). Full size image

To distinguish between subsea permafrost physical properties such as water/ice and gas within sediments, we validated one such acoustic anomaly by drilling a borehole at the site where such an acoustic anomaly was observed and measuring the gas content above the acoustic anomaly and within it (Fig. 9a). Within the observed acoustic anomaly, we measured an increase in CH 4 concentration by two orders of magnitude (Fig. 9b). These data allowed us to interpret the observed acoustic anomaly as gas incorporated within sediments. Multiple pockmarks and bubbles propagating into the water column overlying the seafloor also point to gas releases occurring in the study area (Supplementary Fig. 8).