Postglacial response of Arctic Ocean gas hydrates to climatic amelioration

Serov, Pavel; Vadakkepuliyambatta, Sunil; Mienert, Jürgen; Patton, Henry; Portnov, Alexey; Silyakova, Anna; Panieri, Giuliana; Carroll, Michael L.; Carroll, JoLynn; Andreassen, Karin; Hubbard, Alun

doi:10.1073/pnas.1619288114

Earth, Atmospheric, and Planetary Sciences

Postglacial response of Arctic Ocean gas hydrates to climatic amelioration

Pavel Serov [email protected], Sunil Vadakkepuliyambatta https://orcid.org/0000-0003-1118-1321, Jürgen Mienert, +7 , Henry Patton, Alexey Portnov, Anna Silyakova, Giuliana Panieri, Michael L. Carroll, JoLynn Carroll, Karin Andreassen, and Alun HubbardAuthors Info & Affiliations

Edited by Mark H. Thiemens, University of California, San Diego, La Jolla, CA, and approved May 1, 2017 (received for review November 22, 2016)

June 5, 2017

114 (24) 6215-6220

https://doi.org/10.1073/pnas.1619288114

PDF/EPUB

Significance

Shallow Arctic Ocean gas hydrate reservoirs experienced distinct episodes of subglacial growth and subsequent dissociation that modulated methane release over millennial timescales.

Abstract

Seafloor methane release due to the thermal dissociation of gas hydrates is pervasive across the continental margins of the Arctic Ocean. Furthermore, there is increasing awareness that shallow hydrate-related methane seeps have appeared due to enhanced warming of Arctic Ocean bottom water during the last century. Although it has been argued that a gas hydrate gun could trigger abrupt climate change, the processes and rates of subsurface/atmospheric natural gas exchange remain uncertain. Here we investigate the dynamics between gas hydrate stability and environmental changes from the height of the last glaciation through to the present day. Using geophysical observations from offshore Svalbard to constrain a coupled ice sheet/gas hydrate model, we identify distinct phases of subglacial methane sequestration and subsequent release on ice sheet retreat that led to the formation of a suite of seafloor domes. Reconstructing the evolution of this dome field, we find that incursions of warm Atlantic bottom water forced rapid gas hydrate dissociation and enhanced methane emissions during the penultimate Heinrich event, the Bølling and Allerød interstadials, and the Holocene optimum. Our results highlight the complex interplay between the cryosphere, geosphere, and atmosphere over the last 30,000 y that led to extensive changes in subseafloor carbon storage that forced distinct episodes of methane release due to natural climate variability well before recent anthropogenic warming.

Marine surveys of the Arctic Ocean continental shelf and slope are continuously disclosing new seafloor methane seeps associated with gas hydrate reservoirs (1–3). Gas hydrates are crystalline solids that consist of methane trapped in a lattice of hydrogen-bonded molecules of water (4). Due to their extensive distribution throughout the Arctic and elsewhere, hydrates are an integral part of a dynamic global carbon cycle (5, 6) where methane and heavier gases (i.e., ethane/propane) are sequestered and released over time. Under stable—high-pressure/low-temperature—conditions, gas hydrates constitute a potentially massive natural subseafloor carbon sink and storage capacitor. However, even under stable conditions, some ongoing methane seepage is likely to occur due to preferential fluid migration from deep, thermogenic hydrocarbon reservoirs or due to methanogenesis within organic-rich marine sediments. Despite this finding, under warming and/or depressurization, hydrate dissociation can drive large-scale natural gas release with potentially profound impacts. Abrupt episodes of methane emissions from the seafloor may attain the atmosphere (6) and thereby become a potent feedback for abrupt climate change (5, 7). Methane released into the water column also affects its geochemical signature and pH due to aerobic oxidation, leading to enhanced levels of carbon dioxide (8). However, moderate methane release is regulated by, and is also the basis for, marine chemosynthetic ecosystems that thrive in the vicinity of venting gas seeps, with new extremophiles continually discovered (9–11). Gas hydrates also sculpt and influence seafloor morphology with methane-derived carbonate crusts and pavements formed above gas venting systems, and, furthermore, hydrate dissociation within sediments has been linked to megascale submarine landslides (12), pockmarks (13), craters (14), and gas dome structures (15).

Gas and water that constitute a hydrate crystalline solid within the pore space of sediment remain stable within a gas hydrate stability zone (GHSZ) that is a function of bottom water temperature, subbottom geothermal gradient, hydrostatic and lithostatic pressure, pore water salinity, and the specific composition of the natural gas concerned. Generally, the GHSZ increases in thickness with greater water depth (4). In contrast to other Arctic regions, where gas hydrates remain stable to 300 m below sea level (mbsl), or even shallower in subsea permafrost regions (16), the modern GHSZ along the southwestern Svalbard margin appears deeper at 370 mbsl to 390 mbsl. Here, the relatively warm, ∼2.7 ^oC northward flowing West Spitsbergen Current exerts strong control on the spatial extent and thickness of the GHSZ. It has been argued that recent warming of this current has triggered active recession of the upper GHSZ, thereby promoting enhanced methane seepage (17, 18). An alternative hypothesis suggests that seasonal variations in bottom water temperature drive fluctuations of gas hydrate decomposition and transient methane release (19). To date, gas hydrates have not been observed at or close to the upper edge of the hydrate stability zone offshore of Svalbard. Due to the largely unknown composition of gas in marine sediments coupled with a paucity of cores and actual hydrate samples, previous estimates for the GHSZ (17, 19) are based on theoretical considerations alone, which may be at odds with the actual hydrate stability conditions at the seabed. Thus, the fate of gas hydrates on the Svalbard margin in response to past, ongoing, and future oceanic warming remains unclear.

Here, we present the discovery of intensive cold seep activity clustered on the apexes of several ∼500-m-wide gas hydrate-bearing domes at 370 mbsl to 390 mbsl in Storfjordrenna, northwestern Barents Sea (Fig. 1). Such formations, close to the shallow termination of the GHSZ, have rarely been observed in the Arctic, and their origin has yet to be investigated. We refer to these domes as “gas hydrate pingos” (GHPs) because they are morphologically similar to ice-bearing onshore pingos (20) and their offshore counterparts (21, 22). Terrestrial and offshore pingos form in permafrost regions where water-saturated soils freeze and expand (20, 23). The primary difference between the permafrost-related mounds and the domes imaged here is that, instead of ice, GHPs are formed from methane-derived authigenic carbonates and gas hydrates, which render them susceptible to changes in their ambient temperature and pressure environment.

Fig. 1.

Gas leakage system at Storfjordrenna. Compilation of observations, including seabed topography from high-resolution multibeam data (10-m grid cell), 2D seismic cross-section (300 Hz), and single-beam echosounder data (38 kHz) tracing streams of gas bubbles (gas flares) in the water column. Different colors within gas flares indicate backscattering strength of the reflected acoustic signals (red is the highest values, light green is the lowest). Vertical trains of large dots in the water column show locations of the water samples and, by colors, concentrations of dissolved methane measured. Red and gray marks at the seafloor indicate coring sites with and without gas hydrates, respectively. (*Inset*) Red arrow shows West Spitsbergen Current (warm Atlantic Water); blue arrow shows East Spitsbergen Current (cold Polar Water).

The wider Barents Sea region experienced profound subglacial temperature, pressure, and isostatic variations during the last glacial cycle (24–26). A cooling climate ∼35,000 y ago initiated the growth of the marine-based Barents−Kara Sea ice sheet, providing extensive high-pressure/low-temperature subglacial conditions across the continental shelf off Svalbard (13). Analysis of sediment cores from the region reveal that the ice sheet advanced across the shelf at ∼27,000 calendar (cal.) y B.P. and was at its maximum extent at the shelf break west of Svalbard by ∼24,000 cal. y B.P. (27). After a prolonged period of relative stability, deglaciation commenced rapidly from ∼20,000 cal. y B.P. (28, 29) onward. Hemipelagic mud present in a sediment core from Storfjordenna, some 12 km south of our GHP site, constrains local deglaciation to around 19,000 cal. y B.P. (30). The receding ice sheet left a series of grounding zone wedges and several generations of plow marks, indicating alternating phases of standstill and active, calving retreat (29). Concurrent with and promoting deglaciation, ambient Arctic water of ∼1.5 °C encroached onto the shelf (30). Marine sediment δ¹⁸O records reveal that, during Heinrich event 1 (H1, 15,000–13,000 cal. y B.P.), Bølling and Allerød interstadials (13,000–11,000 cal. y B.P.), and the Holocene Optimum (9,000–8,000 cal. y B.P.), Atlantic bottom water—on average 3 °C warmer—displaced the cooler ambient Arctic water body that was present immediately after deglaciation (30, 31).

Storfjordrenna’s complex environmental history, and that of the wider Barents Sea shelf, raises several important questions in relation to gas hydrate storage and decomposition. Did the GHPs develop as a result of deglaciation or due to more recent ocean warming? How did the GHSZ respond to ice sheet retreat and the subsequent marine incursion of Arctic waters? When and for how long did stable gas hydrates exist during glaciation? How thick were they? Addressing these questions requires a quantitative and unified understanding of the interaction between the ice sheet, ocean, and subsurface methane hydrate reservoir over timescales spanning the last glaciation into the near future. To this end, we characterize the newly discovered site at Storfjordrenna, along with the documented recovery of gas hydrate from the Svalbard−Barents Sea shelf, to provide boundary conditions for a time-dependent coupled ice sheet/GHSZ model that describes the evolution and dynamics of the glacial and subglacial gas hydrate systems in this sector.

Gas Hydrate Pingos and Methane Venting

Five discrete GHPs were geophysically imaged within a 2.5 km² area on the flank of the glacially eroded cross-shelf Storfjordrenna (Fig. 1 and Fig. S1). All of the GHPs have subcircular or elongated shapes with diameters of 280 m to 450 m and heights of 8 m to 10 m. Their existence within a ground zone of vigorous paleo-ice stream activity evidenced by multiple megascale glacial lineations indicates formation after the last ice sheet retreated from the area.

Fig. S1.

Map of the study area with locations of samples and seismic lines. (*Inset*) Red arrow shows West Spitsbergen Current (warm Atlantic water), and blue arrow shows East Spitsbergen Current (cold polar water).

Hydroacoustic observations reveal that four out of five GHPs persistently emit natural gas (Fig. 1). Gas bubbles, represented by hydroacoustic anomalies within the water column, emerge and concentrate from the topographic summits of the GHPs. The area was surveyed three times in May, July, and October 2015, and the observed gas flares were continuous, with many of them rising to at least ∼200 mbsl and the largest ones rising to 20 mbsl. Undetected bubbles—those smaller than the resonance frequency of the echosounder signal (∼0.6 and 0.8 mm at 20 and 200 m water depth, respectively)—cannot be discounted from breaking the ocean surface without being traced (32). Despite the lower hydrostatic pressure, rising methane bubbles gradually shrink by diffusion into the ambient water column (33, 34). Geochemical analysis indicates that gas seepage supplies the water column with up to 130 mL/L of dissolved methane (Fig. 1), some ∼40 times higher than the ambient concentration. A towed camera vehicle equipped with a methane sensor surveyed 0.5 m to 2.0 m above seabed and traced concentrated plumes of dissolved methane associated with the GHPs. The location of the methane plumes along with the gas release sites coincides with the GHP summits, confirming that persistent, focused methane expulsion is closely linked to the specific morphology of each GHP (Fig. 1).

Cores acquired from the GHPs reveal that gas hydrate-bearing hemipelagic sediments with abundant carbonate concretions (Fig. S2) are present in distinct layers below the seafloor (40 cm to 70 cm and 90 cm to 120 cm below the seafloor in GHP’s summits; 120 cm to 130 cm and 205 cm to 220 cm below seafloor in the GHPs’ flanks). Outside the GHPs, sediments do not contain carbonate inclusions, indicating reduced or absent influx of methane. In a pattern identical to gas expulsion and flares, gas hydrates appear exclusively within the topographic highs, comprising multiple layers with different textures that include disseminated, massive, and layered hydrates that occur at various depths within the cores. Some sediment layers exhibit a liquefied, soupy material due to the dissociation of hydrates that is typically observed in recovered sediment cores where the temperature and pressure conditions have changed greatly on opening (35, 36). For all GHPs cores, the released gas was predominantly methane with an unambiguous thermogenic (i.e., depleted) isotopic signature (δ¹³C_average = −47, n = 8; δD_average = −177; n = 8) with additional low admixtures of higher methane homologs (C1/C2−C3_average = 111.3; n = 87) (Table S1).

Fig. S2.

Seafloor images. (A) Blocks of authigenic carbonates colonized by sea anemones. (B) Streams of gas bubbles rising from a seafloor with patches of bacterial mats. (C) Bacteria colonizing tube worms and carbonate outcrops. (D) Carbonate concretion from 50-cm subbottom depth.

Table S1.

Molecular and isotopic composition of gas

Station	Latitude	Longitude	Water depth, m	Sampling device	Interval, cm below seabed	Methane concentration, µM	C1/(C2−C3)	δ¹³C, ‰ vs. vPDB	δD, ‰ vs. vSMOW
CAGE 15-2 HH 904 MC	76° 06.4290 N	15° 58.0500 E	379	Multicorer	5–7	13.00	20.37
					7–9	10.72	16.62
					9–11	11.63	16.41
					11–13	24.64	17.35
					13–15	27.36	17.65
					15–17	31.64	18.02
					17–19	51.20	16.81
					19–21	48.00	16.85
					23–25	79.15	17.81
					25–27	61.57	16.97
					27–29	518.85	22.47
CAGE 15-2 HH 911 GC	76° 06.4150 N	15° 58.0500 E	379	Gravity corer	5–7	149.54	143.88
					15–17	45.81	17.99
					25–27	50.88	13.86
					35–37	323.54	71.20
					45–47	530.87	24.25
					55–57	10,609.44	460.06
					65–67	6,935.08	384.91	−48.4	−184
					75–77	4,802.57	199.96
CAGE 15-2 HH 912 GC	76° 06.4017 N	15° 58.1174 E	379	Gravity corer	5–7	60.45	142.15
					15–17	172.45	120.86
					25–27	117.05	51.09
					35–37	631.76	47.05
					45–47	800.92	55.30
					55–57	1,384.82	76.19
					65–67	2,032.30	97.21
					75–77	2,673.00	158.21
					85–87	12,675.55	623.84	−50.5	−130
					95–97	7,378.60	466.78
CAGE 15-2 HH 913 GC	76° 06.9095 N	16° 00.1476 E	381	Gravity corer	0–2	227.46	67.23
					10–12	378.01	57.16
					20–22	301.45	45.94
					30–32	224.76	27.08
					40–42	7,863.33	129.85
					50–52	7,061.01	99.86
					60–62	3,064.76	78.23
					80–82	3,228.96	150.15
					90–92	8,697.49	226.06
					110–112	1,490.17	177.46
					120–122	12,125.21	377.53
					150–152	14,441.39	364.33
					185–187	27,856.60	655.89	−44.3	−189
CAGE 15-2 HH 914 GC	76° 06.3429 N	16° 02.0971 E	381	Gravity corer	0–2	22.50	56.16
					10–12	109.64	29.95
					20–22	88.63	21.16
					40–42	323.77	26.39
					60–62	4,236.12	203.99
					80–82	5,911.36	239.38
					100–102	7,756.12	231.91
					110–112	3,884.55	139.74
					150–152	3,930.04	147.70
					180–182	3,195.42	146.47
					200–202	5,968.10	238.26	−44.7	−190
CAGE 15-2 HH 918 GC	76° 07.0080 N	16° 02.6314 E	371	Gravity corer	0–2	10.18	8.19
					20–22	157.62	10.45
					40–42	380.68	23.65
					60–62	420.32	26.45
					80–82	9,241.07	192.71	−51.8	−161
CAGE 15-2 HH 920 GC	76° 06.6998 N	16° 00.6489 E	386	Gravity corer	1–3	3.52	18.39
					20–22	8.05	38.43
					40–42	61.91	12.23
					60–62	68.89	13.18
					80–82	101.85	11.36
					100–102	126.04	10.58
					120–122	141.72	11.20
					140–142	331.44	18.74
					160–162	2,768.97	120.72
					180–182	12,281.65	232.46
					200–202	8,965.01	379.38
					220–222	7,754.05	608.42
					250–252	11,298.93	201.12	−44.5	−194
CAGE 15-2 HH 921 GC	76° 06.5342 N	16° 00.1685 E	393	Gravity corer	368–370	5,834.18	14.67	−48.5	−183
CAGE 15-2 HH 938 MC	76° 06.4410 N	15° 58.1590 E	379	Multicorer	0–2	4.19	10.46
					8–11	13.21	17.76
					18–20	31.24	19.22
CAGE 15-2 HH 939 GC	76° 06.4158 N	15° 57.4570 E	387	Gravity corer	60–62	8.09	12.38
					80–82	14.93	19.83
					140–142	45.73	21.21
					160–162	4.71	18.04
					180–182	15.15	20.32
					200–202	80.81	39.52
					240–242	14.50	18.00
					280–282	59.78	21.39	−44.0	−185
					300–302	85.49	27.36
					320–322	27.96	25.55
					340–342	11.00	57.82
					350–352	36.91	20.58

vPDB, Vienna Pee Dee Belemnite standard; vSMOW, Vienna standard mean ocean water.

High-resolution seismic data reveal deep-rooted (>150 m below seafloor) subvertical amplitude masking zones underlying each of the GHPs that we interpret as chimneys through which thermogenic-derived gas migrates upward (Fig. 1). High-amplitude reflectors around these gas chimneys indicate local accumulations of free gas, gas hydrate, or authigenic carbonates. Given these seismic data and the regional geological setting, we infer that an existing fault system within the upper Paleocene−Eocene and Pliocene−Pleistocene sedimentary rocks provides high-permeability zones for upward thermogenic gas migration from underlying hydrocarbon-rich Triassic−Jurassic formations (37–39).

Glacial History and Evolution of Gas Hydrate Stability

During the Last Glacial Maximum, a large grounded ice stream occupied Storfjordrenna and drained ice from a major accumulation center over southern Svalbard (40). Empirically constrained ice flow modeling reveals that grounded ice entered Storfjordrenna ∼35,500 y ago with the onset of glaciation, overriding today’s GHP site (25). Within the next 2,000 y (by ∼33,500 y ago), the Storfjordrenna ice stream had advanced to the shelf break (Fig. 2A). During the next ∼10,000 y, the ice stream was relatively stable, experiencing only minor fluctuations of the ice front, which was mostly pinned to the continental shelf edge. Ice was between 900 m and 1,000 m thick above the GHP site at this time, whereas subglacial temperatures fluctuated between −0.5 °C (the pressure-dependent melting point of ice) and −6 °C dependent on the ice stream configuration (Fig. 2B). Ice sheet retreat commenced around 22,500 y ago in line with global climate amelioration. The active ice stream retreated from the GHP site ∼21,000 y ago, until it attained a stable position 40 km farther upstream around 18,000 y ago. Under continued atmospheric and ocean warming coupled with ongoing eustatic sea level rise, the ice stream retreated back to inner Storfjorden by 14,500 y ago.

Fig. 2.

Evolution of the Storfjordrenna ice stream and postglacial oceanographic changes. (A) Time-lapse setting of the ice stream along the line indicated in *Inset*. GHPs are not to vertical and horizontal scale. (B) Changes of the ice and GHSZ thickness, bottom temperature, sea level (59), and isostatically adjusted seabed at GHP site throughout the last 37,000 y.

By coupling the glacial evolution with a transient gas hydrate model (Methods), we underscore the tight spatial and temporal relationship between GHSZ depth in Storfjordrenna and ice sheet dynamics (Figs. 3 and 4). The GHSZ model essentially solves the conductive heat flux equations based on ambient pressure and thermal conditions provided by ice and/or ocean above and geothermal inputs from below. A more sophisticated multiphase fluid flow model could be adopted, but robust application of such a model requires accurate definition of a wide range of input parameters related to sediment and fluid properties, along with their evolution over time. Given the complex geological and environmental history at our study site, including significant episodes of glacial erosion/isostasy, compaction of sediments, subglacial and marine deposition, and formation and melting of subglacial permafrost, all of these conditions would require accurate parameterization in a multiphase model. Hence, given the available data and the environmental complexity of the study site, we reason that application of a GHSZ model is a more pragmatic and robust approach in this instance.

Fig. 3.

Evolution of GHSZ in outer Storfjordrenna throughout the last glacial cycle. Blue line indicates contours of the ice sheet. Red dashed line shows location of GHP site. GHPs are not to vertical and horizontal scale.

Fig. 4.

Growth and collapse of GHSZ in outer Storfjordrenna. Isostatic movements, subsurface temperature distribution, GHSZ, thickness of ice, and permafrost resulted from our modeling. Gas chimneys and faults are not to vertical and horizontal scale. (A) Setting during the Last Glacial Maximum: ∼200-m-thick GHSZ, patches of subglacial permafrost. (B) GHSZ-free shelf during the H1. Seabed gas efflux is unhampered. (C) Continuous GHSZ on the shelf by the end of the Younger Dryas Interstadial. Gas chimneys intersect with ∼40-m GHSZ.

The GHP site and adjacent shelf were outside of the GHSZ until the onset of ice sheet advance 35,000 y ago (Figs. 2B and 3). Cold subglacial temperatures (−2 °C_mean) combined with high overburden pressures in excess of 8 MPa (equating to 900-m overburden of ice), established a ∼200-m-thick subglacial GHSZ at the GHP site that was sustained for 13,500 y. Around 30,000 y ago, the subglacial GHSZ merged with the subseafloor GHSZ on the continental slope, forming a continuous gas hydrate field across the entire region (Fig. 3). Throughout this glacial episode, the thickness of the subglacial GHSZ varied by around 20%, dependent on ice thickness, basal temperatures, and concomitant overburden pressure (Fig. 2B). After final deglaciation, the impact of an inherited glacio-isostatic depression of ∼85 m at the GHP site promoted the preservation of a 100-m-thick GHSZ up until around 15,500 y ago (Figs. 2B and Fig. 3). Eventually, however, inflowing warm Atlantic Water at 4.0 °C to 5.5 °C associated with the H1 event and the Bølling−Allerød interstadials (30) combined with ongoing isostatic rebound destabilized any remnants of the GHSZ from the area (Fig. 3). Northern Hemisphere cooling during the Younger Dryas stadial at ∼12,000 y ago and the incursion of the cold East Spitsbergen Current (30) initiated a second phase of gas hydrate formation with a ∼60-m-thick GHSZ established across the shelf that once again connected with the persistent offshore GHSZ beneath the continental slope (Fig. 3).

Analogous to H1, the Holocene optimum was likewise associated with an intrusion of warm, ∼4 °C Atlantic Water from outer Storfjordrenna, and led to a further episode of gas hydrate destabilization (Fig. 3). From 8,000 y onward, a steady transition to modern oceanographic conditions, with bottom water temperatures experiencing a steady decline from 4.0 °C to 2.0 °C, somewhat surprisingly promoted moderate gas hydrate growth at the GHP site up to the present. Today, Storfjordrenna hosts two competing water masses: warm and saline Atlantic Water and Arctic Water that is cold and fresh, the interplay of which yields strong seasonal fluctuations in bottom water temperature from 0.5 °C to 2.0 °C, dependent on prevailing synoptic conditions (41). Annual bottom water temperatures observed since the 1950s (42) have, however, remained steady, and thus gas hydrates in the area have remained stable (assuming a similar gas composition to that at the GHP site; Fig. S3).

Fig. S3.

Bottom water temperatures and limits of gas hydrate conditions. (A) Water temperature profiles obtained by CAGE in 2015. (B) Historical water temperature data. (*Insert*) Location of considered temperature profiles. Color of the lines indicates year of measurement.

Varying Methane Leakage Activity

Through synthesis of direct observations with hybrid ice sheet/GHSZ modeling, we demonstrate that an extensive, well-developed subglacial gas hydrate system formed across outer Storfjordrenna during the Last Glacial Maximum. This hydrate system subsequently experienced repeated cycles of reemergence/dissociation during the Late Glacial and Holocene periods driven by changes in oceanographic conditions and gradual glacio-isostatic recovery. Due to its episodic nature, the changes in the GHSZ forced distinct phases of seafloor methane expulsion. During phases when the seafloor was within the hydrate stability envelope, gas hydrate growth incorporated existing natural gas, partially filling sediment pore space and thereby reducing its permeability to ascending fluid flow. Conversely, during phases of gas hydrate decomposition, seafloor gas emissions were amplified due to hydrate-bound gas release and free gas venting from deeper thermogenic reservoirs.

The occurrence of discrete layers of methane-derived authigenic carbonates in shallow sediment cores acquired from the GHPs supports our inference of distinct phases of enhanced methane release since deglaciation. Increased methane flux induces anaerobic oxidation of methane near the seafloor, which produces excess HCO³⁻, thereby enhancing authigenic carbonate precipitation (43, 44). Hence, high methane seepage activity associated with conditions of hydrate dissociation is favorable for carbonate precipitation.

The Barents Sea ice sheet covered the West Svalbard shelf for over 13,500 y, driving continuous gas entrapment in and beneath a thick and extensive subglacial GHSZ. On regional deglaciation, the corresponding abrupt increase in temperature and decreased pressure conditions triggered a period of thinning and shrinkage of the GHSZ (Figs. 2B and 3). Reduced pressure and warmer bottom waters resulted in the complete disappearance of GHSZ within <5,000 y after the ice sheet retreated from shallow regions of the seafloor. Throughout the postglacial period, a ∼5-m-thick section of hemipelagic sediments containing present gas hydrates and authigenic carbonates was deposited across the seafloor (28, 30). Driven by the pronounced warming of bottom water to 5.5 °C (30) from 15,500 y ago onward, any remnant GHSZ collapsed, thereby releasing gas hydrates that had accumulated for more than 18,000 y. Decomposition of gas hydrates caused pore volume expansion and activated large-scale release of formerly hydrate-bound methane that was vented through gas chimneys in the seafloor. Laboratory experiments and numerical simulations of seabed gas dome growth indicate that buoyancy forces and the corresponding enhanced pressure from upwelling methane confined within a gas chimney are sufficient to create seabed domes of a few hundred meters in diameter (45–47). We propose that it was this excess pressure-related doming that initiated the growth of the GHPs around 15,500 y ago.

Corresponding to the Younger Dryas, a ∼1,000-y episode of oceanic cooling stimulated extensive GHSZ regrowth and the cessation of methane seepage across the shelf (Figs. 3 and 4C). Gas hydrate heaving, a process analogous to frost heave under permafrost conditions, also would have contributed to sediment upheaval within GHPs at this time. Successive fracturing of sediments caused by excess pore pressure would have led to cracks that eventually fill with hydrates (48), thereby leading to further GHP volume expansion.

A rapid recession of the GHSZ took place, associated with a warming period of bottom water at the Holocene Optimum (Fig. 3). From ∼6,500 y ago onward, oceanographic conditions were broadly comparable to those today, with Arctic-derived bottom waters (<2 °C) prevailing in outer Storfjordenna. These cooler oceanic conditions gradually led to the establishment of a new GHSZ up to 60 m thick that has persisted through to the present day (Fig. 3). Our analysis demonstrates that complex changes in temperature and pressure conditions led to episodic gas hydrate formation in outer Storfjordrenna, which strongly modulated seafloor methane release and the formation of authigenic carbonates and GHPs during the Late Pleistocene and Holocene.

Besides Storfjordeanna, several glacial troughs with depths in excess of 350 mbsl have been eroded into the Barents and Kara Sea shelf (Fig. 4A). These troughs and associated deeper shelf areas must have developed extensive GHSZ during the last glaciation that subsequently experienced episodic phases of collapse and reemergence driven by changing subglacial, isostatic, and oceanographic conditions (49–51). Given the abundance of hydrocarbon provinces within these formerly glaciated margins, we propose that the GHPs we document here could be more common and extensive across the Arctic, where submarine gas hydrate systems exist. Recent surveys off West Greenland support this proposition, where hydrate-bearing seafloor features appear to be associated with deep gas migration channels (52). Furthermore, across the East Greenland shelf, δ¹³C records in benthic and planktonic foraminifera indicate at least three methane release episodes since deglaciation related to dissociating hydrates (53). It is also likely that many GHPs that reside outside of the present-day GHSZ have collapsed, forming large depressions, a phenomenon that has been widely reported in previously glaciated trough systems in the Arctic (13, 14, 54).

Despite considerable seafloor methane seepage from formerly glaciated Arctic shelves, the actual flux of methane that attains the atmosphere remains unconstrained. Recent studies show that a broad gas seepage area extending along the Northwestern Barents sea from 74° to 79° contributed only 0.07% to net atmospheric methane (3). This finding resonates with recent airborne measurements revealing a distinct absence of high atmospheric methane concentration during the summer (55). The role of the water column in critically regulating methane transfer to the atmosphere is not fully understood, and it remains unclear as to whether oceanic methane degradation has limits where large and abrupt fluxes of seafloor release could overcome filter systems, thereby forcing a potent atmospheric feedback, as has previously been proposed.

Earth has experienced a wide range of climate extremes over its geological history (56), such as the Permian−Triassic catastrophe 252 million years ago (57) when both carbon dioxide and methane were released on a massive scale into the atmosphere. It has recently been proposed that such an event was reinforced by global-scale hydrate dissociation and methane release triggered by initial global warming after a prolonged “Snowball Earth” glacial episode (58). Our inferences regarding a glacial gas hydrate capacitor are worth consideration when investigating the causes of past episodes of global-scale gas methane release evident in the geological record.

Despite the growing number of seep-related features that have been recently discovered across the seafloor of the Arctic, shallow gas hydrate systems remain poorly understood and documented, particularly where they have undergone a complex environmental history. This study reveals that abrupt changes in pressure and temperature conditions associated with the interplay of grounded ice, postglacial isostatic rebound, and influx of variable ocean currents all critically modulate the GHSZ, thereby driving distinct episodes of natural gas storage and release. To date, these processes have not been well described or quantified, and any attempts to understand the past and determine the future impact of Arctic methane emissions on global climate need to comprehensively account for them.

Methods

Description of methods of seismic and hydroacoustic data acquisition and processing, sediment sampling, and geochemical analyses are provided in SI Methods. SI Methods also contains an extensive description of the ice sheet model and the conductive heat flux model of GHSZ.

SI Methods

We acquired data during Centre for Arctic Gas Hydrates, Environment and Climate (CAGE) 15-2, CAGE 15-5, and CAGE 15-6 cruises onboard RV Helmer Hanssen in May−June, July, and October 2015. Some of the multibeam data were collected in 2013 using the same Simrad EM 300 multibeam echosounder system operated at a 30-kHz sound wave frequency and 135° angular coverage.

Single-beam echosounder data were acquired using a Simrad EK 60 system operating at frequencies of 18, 38, and 120 kHz and beam angle of 6.5°. For gas flare mapping, we used a 38-kHz echosounder signal. This system effectively detects streams of gas bubbles in the water column due to high acoustic impedance contrast between water and free gas (33). A persistent issue involved in hydroacoustic mapping of gas bubbles is that some of the bubble streams deviate from the vertical due to currents and, thus, leave the echosounder’s footprint (59). However, reliable flare length estimates were obtained using a dense net of survey lines with 80-m spacing.

During the CAGE 15-2 cruise, we deployed a towed vehicle assembled on the multicorer frame and equipped with a multidisciplinary instrumentation in support of oceanography (MISO) 16-megapixel digital still color camera, CTD (conductivity, temperature, density) (SBE 25plus) real-time system, providing both altimetry and depth data, sensor of dissolved methane HydroC Plus, and a six-Niskin bottle array. The camera took high-resolution images every 10 s, transmitting low-resolution replicate images through the data cable in real time, giving a possibility of visually guided sampling. We performed seven tow camera dives with towed vehicle sea bottom tracks of >12.8 km and a total survey duration of 10.3 h.

For 2D reflection seismic surveys, we used two mini generator−injector air guns (980 cm³ total air gun volume); this generates a high-frequency signal, which is suitable for studying the shallow subsurface with 5 m vertical resolution. The shooting rate was 3 s, and the sampling rate was 0.25 ms.

For sediment sampling, apart from the tow camera multicore system, we used a gravity core with a 6-m-long core barrel, holding an inner plastic liner of 10 cm diameter. The total weight of the sampling device is ∼2 tons, which was enough to retrieve up to 369-cm-long sedimentary cores. We visually described and subsampled all cores onboard. For subtracting a gas phase from the bottom sediment matrix, we applied a conventional headspace gas sampling technique. Bulk sediments (5 mL) were transferred into 20-mL headspace glass vials, containing 5 mL of 1-molar NaOH solution and two glass beads, capped with rubber septa, sealed with aluminum crimp caps and shaken afterward. We performed gas analyses onboard using a ThermoScientific Trace 1310 gas chromatograph equipped with a flame ionization detector (GC-FID) and a ThermoScentific TG-BOND alumnia (Na₂SO₄) 30 m × 0.53 mm × 10 μm column.

Water samples from the tow camera system and CTD rosette were transferred into 100-mL glass bottles, which were capped and sealed. Before GC-FID analyses, we added 5 mL of laboratory pure N₂ gas to each of the vials and shook them for 2 min.

Ice-Sheet Evolution and Subglacial Temperature Modeling.

Modeled data relating to the thermodynamic evolution of the Barents−Kara Sea Ice Sheet through the last glacial cycle were implemented from an optimum ice sheet reconstruction previously developed for the Late Weichselian Eurasian ice sheet complex (25). Briefly, the 3D thermomechanical ice sheet model applied is a first-order approximation of the Stokes equations (60–62). The approach adopted to solve the stress/strain field equates to the L1L2 classification of higher-order models and includes longitudinal (membrane) stresses (63). Surface mass balance is derived using a distributed degree day calculation based on a reference seasonal climatology from mean (1950−2000) precipitation and temperature patterns (WorldClim, www.worldclim.org) (63). A further mass balance term for sublimation is incorporated across the High Arctic. Calving losses at marine terminating margins are coupled to changing sea level using a standard empirical function relating calving flux to ice thickness, water depth, and calving front geometry (64). Basal sliding is determined by a Weertman sliding law (65), adjusted using an exponential decay function to initiate basal motion at subglacial temperatures ≤0.75 K below the pressure melting point. A spatially variable distribution of geothermal heat flux is applied, interpolated onto the model domain from core measurements sourced from the Global Heat Flow Database (66). However, geothermal heat flux values around the western Svalbard margin are largely dominated by low mean continental shield values close to 50 mW⋅m⁻². Isostatic loading is implemented using an elastic lithosphere/relaxed asthenosphere scheme (67), which provides a computationally pragmatic solution in the absence of a full spherical Earth model.

The optimum experiment presented was hot-started from an initial modeled ice configuration at 37,230 y ago, a period of relatively warm interglacial climate conditions concurrent with minimal terrestrial-based ice cover. Ice extent, thickness, and the loaded topography at this time were determined from a previous Mid-Weichselian (marine isotope stage 4) experiment with sufficient spin-up time for the ice sheet and isostatic loading to attain a transient equilibrium with the forcing climate. From this inherited configuration, the model was integrated forward in time on a finite difference grid with a resolution of 10 km through perturbations in temperature using the North Greenland Ice Core Project (NGRIP) δ¹⁸O record (68), and global eustatic sea level (69). The deglaciation of the Barents and Kara seas within the optimum experiment has been tightly constrained using a range of geophysical data, including use of an independent glacial isostatic adjustment model to compare transient ice load distributions against observed relative sea level histories (24). Further lines of data used include radiocarbon dates for the timing of maximum extension and pace of retreat, cosmogenic exposure ages for probable surface profiles of the ice sheet, and geomorphological flow sets as indicators of switching ice flow patterns during ice retreat (25).

Transient Modeling of Subsurface Heat Flow and Gas Hydrate Stability.

The subsurface evolution of heat flow along the profile (Fig. 2A) over the last 37,000 y was reconstructed using a 2D finite difference heat flow model (70). The seafloor depths were sampled along the 2D profile (Fig. 1) every 100 m from the isostatically adjusted seafloor depths. A 1,201 × 400 cell temperature grid (100 m × 5 m cell dimensions) was generated, with upper boundary at the sea level and basal boundary 2 km below sea level. The initial temperature conditions were defined using thermal gradients from World Heat Flow Database (66) close to the study area. We produced an initial temperature grid assuming a 1D linear temperature gradient, with ocean bottom water temperatures adjusted to 37,000 y ago (30, 31). Thermal diffusivity is constrained by variable thermal conductivities along the profile from World Heat Flow Database (66), with average density of 1,900 kg⋅m⁻³ for the bulk sediments and an average specific heat capacity of 2,000 J⋅kg⁻¹⋅K⁻¹ (19). The 2D finite difference heat flow model (70) is run for 37,000 y, accounting for the environmental variations such as ice thickness, isostatic adjustments, sea level, and bottom water/ice temperatures over the period. We assume no significant in situ generation of heat by the formation and dissociation of gas hydrates. The 2D temperature profile was sampled at 100-y intervals to estimate the GHSZ thickness.

To identify the base of hydrate stability, we integrate results from our heat flow model with theoretical hydrate stability phase diagrams generated using the CSMHYD (publicly available executable for prediction of hydrate phase equilibrium) program (4), which uses an algorithm based on Gibbs energy minimization, accounting for different pressure and temperature conditions, the composition of gas forming hydrates, and the presence of inhibitors of hydrate formation (e.g., dissolved ions). Using the gas composition from gas hydrates, and a constant pore water salinity of 35‰, the CSMHYD program estimates the pressure at which hydrates are stable for any given temperature. With temperature constrained from the heat flow model, we estimate the hydrostatic pressure at each cell location using standard values for density of seawater (1,027 kg⋅m⁻³) and acceleration due to gravity (9.8 m⋅s⁻²). This pressure grid is then compared with theoretical predictions from the CSMHYD program to determine the depth hydrate stability at each cell location for each 100-y time interval. Thus, our model simulates transient behavior of an interface separating a zone where methane and water coexist (below) and zone of gas hydrate and water equilibrium (above). Such a model design allows us to predict the dynamic behavior of GHSZ and its collapse. Time lag between GHSZ retreat and seabed methane venting is an inevitable discrepancy of the model. Comparison of conductive heat flow models with alternative mass transport models indicates a possible time lag of hundreds to a thousand years, which would not affect our main conclusions (71). Application of mass transport model is not possible, due to largely unknown properties of the sediment matrix and its complex changes during the last 35,000 y.

Acknowledgments

We thank the crew of the RV Helmer Hanssen for their assistance in sediment and water sampling and acquisition of the seismic data. We also thank Kate Waghorn for processing 2D seismic lines. The research is part of the Centre for Arctic Gas Hydrate, Environment and Climate and was supported by the Research Council of Norway through its Centres of Excellence funding scheme Grant 223259.

Supporting Information

Supporting Information (PDF)

Download
998.28 KB

References

1

N Shakhova, et al., The East Siberian Arctic Shelf: Towards further assessment of permafrost-related methane fluxes and role of sea ice. Philos Trans A Math Phys Eng Sci 373, 20140451 (2015).

Featured Topics

Articles By Topic

Featured Topics

Articles By Topic

Featured Topic

Articles By Topic

Significance

Abstract

Sign up for PNAS alerts.

Gas Hydrate Pingos and Methane Venting

Glacial History and Evolution of Gas Hydrate Stability

Varying Methane Leakage Activity

Methods

SI Methods

Ice-Sheet Evolution and Subglacial Temperature Modeling.

Transient Modeling of Subsurface Heat Flow and Gas Hydrate Stability.

Acknowledgments

Supporting Information

References

Information

Published in

Classifications

Submission history

Keywords

Acknowledgments

Notes

Authors

Affiliations

Notes

Competing Interests

Metrics

Altmetrics

Citations

Cite this article

Cited by

View options

PDF format

Login options

Recommend to a librarian

Purchase options

Restore content access

Share

Share article link

Share on social media

Further reading in this issue

Response time in economic games reflects different types of decision conflict for prosocial and proself individuals

Viral and cellular mRNA-specific activators harness PABP and eIF4G to promote translation initiation downstream of cap binding

Auxin steers root cell expansion via apoplastic pH regulation in Arabidopsis thaliana

Multiple sclerosis and gut microbiota: Lachnospiraceae from the ileum of MS twins trigger MS-like disease in germfree transgenic mice—An unbiased functional study

Bodily maps of emotions

Evidence of a social evaluation penalty for using AI

Sign up for thePNAS Highlights newsletter