Methane production controls in a young thermokarst lake formed by abrupt permafrost thaw.

Methane (CH4 ) release to the atmosphere from thawing permafrost contributes significantly to global CH4 emissions. However, constraining the effects of thaw that control the production and emission of CH4 is needed to anticipate future Arctic emissions. Here are presented robust rate measurements of CH4 production and cycling in a region of rapidly degrading permafrost. Big Trail Lake, located in central Alaska, is a young, actively expanding thermokarst lake. The lake was investigated by taking two 1 m cores of sediment from different regions. Two independent methods of measuring microbial CH4  production, long term (CH4 accumulation) and short term (14 C tracer), produced similar average rates of 11 ± 3.5 and 9 ± 3.6 nmol cm-3  d-1 , respectively. The rates had small variations between the different lithological units, indicating homogeneous CH4 production despite heterogeneous lithology in the surface ~1 m of sediment. To estimate the total CH4 production, the CH4 production rates were multiplied through the 10-15 m deep talik (thaw bulb). This estimate suggests that CH4  production is higher than emission by a maximum factor of ~2, which is less than previous estimates. Stable and radioactive carbon isotope measurements showed that 50% of dissolved CH4 in the first meter was produced further below. Interestingly, labeled 14 C incubations with 2-14 C acetate and 14 C CO2 indicate that variations in the pathway used by microbes to produce CH4 depends on the age and type of organic matter in the sediment, but did not appear to influence the rates at which CH4  was produced. This study demonstrates that at least half of the CH4 produced by microbial breakdown of organic matter in actively expanding thermokarst is emitted to the atmosphere, and that the majority of this CH4 is produced in the deep sediment.

INTRODUCTION

Thawing permafrost causes previously frozen organic matter (OM) to be available for microbial degradation and emits additional greenhouse gases (Zimov et al., 1997, 2006). At the global scale, this may exacerbate global warming, in particular with regard to methane (CH4) emissions (Schuur et al., 2015). Increasingly, rapid abrupt thaw of large volumes of ice‐rich permafrost may be augmenting arctic CH4 emissions (Schneider von Deimling et al., 2015; Turetsky et al., 2020; Walter Anthony et al., 2018, 2021). To anticipate future Arctic CH4 emissions, it is necessary to understand processes that act at local scales to constrain CH4 production and emissions.

Historical and current carbon emissions follow the RCP 8.5 to <1% error which predicts that ~40–50% of circumpolar permafrost will abruptly thaw by 2100 (Schwalm et al., 2020). With abrupt thaw comes increased thermokarst lake formation, which enhances CH4 emissions (Walter Anthony et al., 2018). Abrupt thaw begins with the collapse of ice‐rich permafrost and the onset of thermokarst activity that lowers the ground surface relative to the water table. The ground surface eventually falls below the water table and a lake forms. As the thaw of permafrost progresses, the lake deepens and the margins of the lake collapse into the center, which causes sediment from the margins, which are sometimes older than the bottom sediments, to be deposited on top (Figure 1). The lake and talik then expand laterally. The large volume of recently thawed anoxic sediments emit CH4 produced by the microbial breakdown of OM (Desyatkin et al., 2009). Thermokarst lakes formed by the abrupt thaw of ice‐rich permafrost are becoming increasingly prevalent as the Arctic warms.

FIGURE 1

Schematic thermokarst lake formed by abrupt thaw. The collapse of ice‐rich permafrost causes ground subsidence and leads to lake formation. Erosion of older sediments exposed on the margins (forest soil, paleo‐active layer) as the lake expands is deposited on top of the old forest floor (which may appear as a peat horizon). As the taberal sediments (thawed permafrost) deepens, the collapse of permafrost continues and expands laterally while in‐situ permafrost thaw occurs in taberal sediments below. Aquatic plants and organic matter from the surrounding landscape are then deposited under lower sedimentation rates (recent sediments)

In surficial Arctic soils, CH4 production appears to vary by 3 to 4 orders of magnitude with wild fluctuations between different environments and studies (Treat et al., 2015). Measurements of CH4 production are typically obtained by incubating sediment at a set temperature that resembles as much as possible, the environment. Large variations are also associated with CH4 production rates in the talik underlying thermokarst lakes. For example, in thermokarst lakes in interior Alaska, the CH4 production ranged between 2 and 35 nmol cm−3 d−1 in Doughnut and Vault Lakes (Martinez‐Cruz et al., 2018), respectively, but values of up to 351 nmol cm−3 d−1 are also reported in Vault Lake (Heslop et al., 2015), while incubation conditions and temperature were similar in both studies (3–4°C). Other studies show that the sediments of glacial lake Otto and thermokarst lake Goldstream in Alaska have negligible CH4 production (Sepulveda‐Jauregui et al., 2018), while the top layer of sediments of a thermokarst lake in the North Slope of Alaska has CH4 production of around 1000 nmol cm−3 d−1 (de Jong et al., 2018). There is therefore considerable variability in the estimates of CH4 production in rapidly thawing permafrost regions.

A comprehensive understanding of the controls on CH4 production is essential for predicting future CH4 emissions. The CH4 production appears to correlate with OM concentration (Heslop et al., 2015) and OM biolability (Heslop et al., 2019), while it is also related to the depositional history and the resident microbial community (Holm et al., 2020). However, the ability to anticipate future CH4 emissions is limited by knowledge of OM fermentation in thermokarst sediments (Sepulveda‐Jauregui et al., 2015, 2018). In long‐term incubations of thawed permafrost, the proportion of CO2 and CH4 production appears to stabilize at 1:1 (Holm et al., 2020; Knoblauch et al., 2018; Walter Anthony et al., 2014) indicating complete fermentation of OM. However, the breakdown of OM in the methanogenic zone of the talik does not always tend toward this stoichiometry or molar ratio in lake emissions for several reasons: (1) iron oxides found in sediments increase CO2 production to the detriment of CH4 formation (Vigderovich et al., 2019), (2) OM degrades in quality and influences CH4 formation (Harwood et al., 1998), (3) CH4 oxidation re‐oxidizes CH4 to CO2 (Bastviken et al., 2002; Segarra et al., 2013; Elder et al., 2018; Martinez‐Cruz et al., 2018), and (4) groundwater transports the soluble CO2 from watershed soils more than the CH4 (Kling et al., 1991; Maberly et al., 2013; Weyhenmeyer et al., 2015).

Field observations in thawing permafrost show significant variations in the proportions of CO2 and CH4 produced (Hodgkins et al., 2014), which has an important control on the resulting warming scenarios (Knoblauch et al., 2018; Lokshina et al., 2019; Sepulveda‐Jauregui et al., 2015; Walter Anthony et al., 2014). To constrain CH4 production in regions of quickly degrading permafrost, a young and rapidly expanding thermokarst lake was studied by combining porewater geochemical profiles and two independent methods of measuring CH4 production. We constructed a CH4 budget of the lake and compared it with CH4 emission measurements to quantify the sources and sinks of CH4.

MATERIALS AND METHODS

Site characterization and sampling

Big Trail Lake (BTL) (64°55'9.47"N, 147°49'20.08"W) is situated 8 km north of Fairbanks Alaska in the Goldstream valley, at the junction between Goldstream Creek and Big Eldorado Creek, and downstream ~4 kilometers to the west of Goldstream Mine in Fox. BTL is an actively expanding thermokarst lake that currently has a surface of about 40,000 m2 and an average and maximum depth of 1.2 m and 4.1 m, respectively (Figure 2a). The lake has a record of thermokarst activity, including vigorous CH4 seeps (Walter Anthony et al., 2021). Degrading ice wedges in the lake bottom are likely driving the ground subsidence that led to the formation and expansion of the lake. Airborne photographic records indicate that BTL formed from a wetland sometime between 1949 and 1967 (Walter Anthony et al., 2018). Under increasing polar amplification, Goldstream Valley has experienced its three summers with the most precipitation between 2015 and 2020 in the 91‐year record (Douglas et al., 2020b).

FIGURE 2

(a) Bathymetry map of Big Trail Lake indicating the two coring sites with white dots. (b) and (c) Lithological units from sites S and D (respectively) color coded alongside the dry weight TOC % (white circles). Ages according to 14C dating (yr BP, calibrated) are presented on the right, located parallel to their sampling depth (site D only). The 14C dates are from TOC of bulk sediment, not macrofossils. The color scheme for the units is consistent throughout this manuscript. The photo associated with site D core 2 was taken at the time of sampling

The two sampling sites which were selected in BTL reflect the lake's changing shoreline history. Site S (shallow) was selected from the northern portion of the lake characterized by fen vegetation in 1949. Site D (deep), located in the western basin of the lake, was covered by boreal forest in 1949. Both locations had changed to open thermokarst lake by 1967, but currently, the water depth at site D (2.5 m) is over four times greater than in site S (0.6 m). This suggests greater ground ice melt and surface subsidence at site D compared to site S. Some possible explanations could include lower ground ice content at site S, more groundwater enhancement of thaw at site D, or more protection against thaw at site S, perhaps due to thicker ice‐free sediments above the ice wedges (Kessler et al., 2012).

Porewater and sediment geochemistry

Two sediment cores from each site were extracted in March 2020 by percussion coring following (Heslop et al., 2015). From the first core, porewater was extracted with Rhizons following (Dickens et al., 2007) within 24 h after sampling, and dissolved CH4 was extracted following (Hoehler et al., 2000) within 3 h after sampling. The second core was sampled for long‐term incubations and short‐term incubations (Supplemental Methods). Unless stated otherwise, the geochemical analyses were performed at Ben Gurion University of the Negev, Israel. The water content of the sediment was measured by weight difference from drying an aliquot of sediment at 60°C for 4 days. The CH4 concentrations were measured by GC‐FID (Thermo) with a Packed ShinCarbon column (precision of ±2 μM). Total organic carbon (TOC) was analyzed by pyrolysis with a Rock‐Eval Analyzer (Vinci). Acetate concentration was determined by two‐stage ion chromatography at Aarhus University following (Glombitza et al., 2014). Dissolved Fe(II) was measured using ferrozine assay (Stookey, 1970). Mn2+ and major cations (Na+, Ca2+, Mg2+) were measured by ICP‐OES on a Perkin‐Elmer optima 3300 with a precision of 2%. SO4 2− was measured by ion chromatography on a Compact IC Flex, Metrohm with detection limit of 6 µM. Hydrogen sulfide was measured as per Cline (1969). Nitrite (NO2 −) and Nitrate (NO3 −) were analyzed in the Interuniversity Institute for Marine Sciences in Eilat, Israel, using a colorimetric method (Grasshoff et al., 2009) by a Flow Injection Autoanalyzer (Lachat Instruments Model QuickChem 8500).

Measurements of the stable carbon isotope of CH4 were done by GC‐IRMS on Thermo Delta V. Isotopic measurements are expressed in δ (‰) = 103 ((Rsample/Rstandard) − 1), where R is 13C/12C and standards refer to the Vienna Pee Dee Belemnite (VPDB) with an analytical error of ±0.1 ‰. Dissolved inorganic carbon concentration and isotope ratio were also measured by gas‐source IRMS with precisions of ±0.1 mM and ±0.1‰, respectively. Radiocarbon measurements of sediment TOC and CH4 were performed at the A. E. Lalonde AMS Laboratory, University of Ottawa, Canada, on a 3 MV tandem accelerator mass spectrometer with preparation methods described in Crann et al. (2017) and Murseli et al. (2019). Magnetic susceptibility (MS) of core sediments followed the procedure of Kletetschka et al. (2021) at the University of Alaska Fairbanks. Further information regarding sampling procedures and measurements is available in Supplemental Methods.

Sediment incubations

Long‐ and short‐term sediment incubations were conducted at 4°C, which reflects average temperature in the sediment of thermokarst lakes from the region (Heslop et al., 2015). Long‐term sediment incubations were used to calculate the CH4 production rate based on the accumulation of CH4 in the headspace of the bottle over a timescale of months under anoxic conditions (see supplemental methods, Section 3.1), reported in units of nanomole per cubic centimeter of sediment per day. Short‐term sediment incubations amended with 14C‐labeled tracers were performed following Beulig et al. (2018) and conducted throughout time periods of hours to days. Briefly, 10 kBq of 14CH3COO− was amended to anaerobic sediment in a sealed glass bottle and incubated for a time interval. Recovery of 14C that was transferred in the time interval to CO2 (the dissolved inorganic carbon pool) and CH4 allowed the determination of the acetate turnover rate (ATR) where [Ac] is the acetate concentration from the porewater extraction (Figure S1b), φ is porosity and AAc, ACH4 and ACO2 are the activity of acetate, CH4 and CO2 respectively, at the end of the incubations. The acetate oxidation rate (AOR) was then calculated as

And the acetate fermentation rate (AFR) was calculated as

The CO2 reduction rate (CO2RR) was calculated from a separate incubation with labeled 14CO2, as where [DIC] is the dissolved inorganic carbon concentration.

ATR is the sum of AOR and AFR. A fraction of ATR, the AOR, is oxidized to CO2 while the balance, AFR, is fermented to CH4 and CO2. The CO2RR on the other hand is an independent measurement of a separate microbial process. The sum of AFR and CO2RR is the CH4 production rate, also expressed in nanomole per cubic centimeter of sediment per day.

One limitation of the 14C tracer incubation is the uncertainty on the measurement of acetate concentrations. Acetate concentration measured from the porewater was around 5 µM but had a minimum value of 2 µM and a maximum of 11 µM with no systematic trend (Figure S1b). Therefore, to calculate the acetate fermentation rate (AFR) and the acetate oxidation rate (AOR), the acetate concentration used was the average acetate concentration of 5 µM. This is around the expected acetate concentration in sediments that have active microbial communities in equilibrium (Glombitza et al., 2016; King, 1991) but is lower than the concentration found in some bogs and fens (Duddleston et al., 2002; Patzner et al., 2020).

Gas bubble analysis

Gas traps were deployed around BTL in the summer of 2019 below the surface to collect bubbles emanating from the lake bottom for isotopic analysis. The traps were left between 1 and 17 days and the gas collected subsequently. Stable isotope analysis of CH4 gas was performed at Florida State University on a Thermo Delta V and is expressed using the same notation and reference as above.

RESULTS

Lithology and diagenetic zones in BTL sediments

The lithology between the two sites was similar (Figure 2 and Figure S4) and was assumed as laterally continuous. The following unit numbers and description thus refer to both sites. The lithology is described from the bottom up.

Unit 1 is the deepest unit that was sampled. It is likely substantially composed of the original deposits from the Goldstream formation as described by Péwé (1975). This unit would have been part of the seasonal active layer before the onset of thermokarst talik formation (Figure 1) and may have been thawed year‐round during an earlier climate regime (e.g., Holocene Thermal Maximum, (Heslop et al., 2019)). The TOC content of Unit 1 was lower than the other units, 1–7%, and The TOC 14C age was 2980 ± 30 yr BP (Figure 2). Unit 1 stretched from 45 and 71 cm downwards at sites S and D, respectively. It was not possible to determine the lithology further down than the bottom of the Big Trail Lake core; however, similar thawed yedoma sediments were observed down to 16 m at nearby Goldstream Lake (Walter Anthony et al., 2021).

Unit 2 was a peat layer with TOC between 22% and 28% (Figure 2) with no visible clastic sediments (Figure S4). The 14C age of the TOC was 1730 ± 30 yr BP (Table 1), younger than the overlying units, which likely represents a mixture of ages. Unit 2 had pieces of wood, sphagnum moss and is interpreted as the Holocene forest floor (Figure 1), which was initially flooded at the onset of thermokarst. Therefore, it incorporated young wetland plants. Unit 2 thickness was 20 cm (site S) and 28 cm (site D).

TABLE 1

14C ages and Fraction modern (Fm) of CH4 bubbles released from boreholes after sediment cores were extracted from the lake, CH4 produced in sediment incubations, and sediment TOC

Sample	Year BP (calibrated)	±	Fraction modern	±	AMS Laboratory sample ID
CH4 borehole bubbles (site D)	10,760	45	0.262	0.001	UOC‐14406
CH4 borehole bubbles (site S)	6320	30	0.455	0.002	UOC‐14021
CH4 produced in incubation (site D, Unit 1, 70–90 cm)	3294	31	0.663	0.002	UOC‐14410
TOC Unit 4 (site D, 4 cm)	4332	58	0.583	0.004	UOC‐14401
TOC Unit 3 (site D, 38 cm)	2194	34	0.760	0.003	UOC‐14402
TOC Unit 2 (site D, 55 cm)	1727	32	0.806	0.003	UOC‐14403
TOC Unit 1 (site D, 73 cm)	2977	37	0.690	0.003	UOC‐14404

For sediment samples, the sampling depth is denoted in parentheses.

Unit 3 was an organic‐rich silt (Figure 2 and Figure S4). The top of the unit was darker than the bottom. The thickness varied between ~10 cm at site S and ~32 cm at site D. The TOC content of Unit 3 was 6–7% and its 14C age was 2200 ± 30 yr BP (Figure 2). The sediments from Unit 3 were likely derived from the Goldstream formation but contained more OM than Unit 1. This is likely a result of lake margin erosion, homogenization and deposition in the lake (Figure 1).

The shallowest unit, Unit 4, was a black, organic‐rich mud layer (gyttja) containing organic carbon recently deposited such as leaves and macrophytes (bryophytes and other submerged aquatic plants) in a silt matrix. The thickness was ~15 cm at site S and ~10 cm at site D. Total organic carbon (TOC) was 5–12% (Figure 2). The TOC of the bulk sediment from a depth of 4 cm at site D was dated by 14C to an age of 4340 ± 60 years before present (yr BP) (Table 1). The greater thickness at site S likely reflects the shallow depth of this site (0.6 m water depth) where a larger amount of macrophytes were observed compared to site D (2.5 m water depth). The older 14C age of Unit 4 compared to Unit 3 may be explained by older Pleistocene‐aged yedoma sediments mobilized from lake margins and redeposited in the lake, and/or incorporation into aquatic plant tissues 14C‐depleted CO2 originating from oxidation of 14C‐depleted CH4 that bubbled from deep within the talik.

Carbon isotopic composition of CH4 and organic matter

The CH4 concentrations in the porewater were at or above CH4 saturation except in the upper 10 cm of site S (Figure S1). Higher CH4 than saturation, in some samples, indicated that CH4 bubbles were present in the sediment at the time of sampling. The δ13CCH4 profile provided more information on the sources of this CH4 (Figure 3). Whereas in core S, the δ13CCH4 oscillated between −55 and −60 ‰ with only a slight increase between 20 cm and 50 cm, site D had a clear trajectory. In site D, δ13CCH4 decreased from −60 ‰ at 5 cm to −63 ‰ at 15 cm followed by a general increase from −63 ‰ at 20 cm to −53 ‰ at 65 cm. The measurements between 70 cm and 90 cm showed a general decrease back to −60 ‰. These results contrasted with the δ13CCH4 obtained from incubations of sediment from the various depths (Figure 3).

FIGURE 3

Isotopic composition of CH4 in BTL porewater and incubations collected from site S (a) and site D (b). White squares denote the porewater CH4, white circles denote the incubations. The error bar of CH4 produced in incubations is equal to one standard deviation (n = 2 separate incubations). Dashed line denotes the ebullition end member, which is the mean δ13CCH4 of a bubble survey in BTL (n = 30, σ = 6.3 ‰, Table S1). Background colors represent the lithological units described in Figure 2

The δ13CCH4 of measured in incubations represents the CH4 that is formed directly from the breakdown of OM at a specific depth horizon. In both sites, the incubation δ13CCH4 produced values around −52 to −55 ‰ in the lower part of the cores. In the first 20 cm of both sites, more negative δ13CCH4 values were obtained. The more negative δ13CCH4 values in the first 20 cm were most pronounced at site S, where the δ13CTOC showed a similar trend of the same magnitude (Figure S5). Both sites had less negative δ13CCH4 in the incubations than in the porewater profile.

Measurement of CH4 bubbles from the water column had a CH4 content of 78% ± 10%. The δ13CCH4 varied from −73.5 to −56.7 ‰. There was a relationship between the type of ebullition and the δ13CCH4. Vigorous ebullition sources, also called “Hotspots,” had the strongest depletions in 13C (−73.5 to −59.9 ‰), whereas weaker seep sources had less depletion of 13C (−67.2 to −56.7 ‰) with one outlier at −43.0 ‰ (Table S1). Vigorous ebullition sources may result from bubble migration from deeper in the talik, whereas weaker seep sources may result from shallower bubble migration (Walter Anthony et al., 2021; Walter et al., 2008). It was not possible to constrain the type of ebullition migrating into the shallow sediment. Thus, the average value of the bubbles (−64 ± 6 ‰) was selected as an endmember for calculations.

The CH4 recovered from bubbles outgassing from the lake bottom was older than CH4 produced during the incubation of sediment from 60–90 cm depth (Table 1). While the age of CH4 in bubbles collected at the lake surface (emitted from the boreholes immediately after core retrieval) was dated at 10,764 ± 44 and 6320 ± 32 yr BP, the CH4 recovered from bottles that were incubated had an age of 3300 ± 28 yr BP (Table 1). This age is similar to the 14C age of the TOC of the incubated sediment (Unit 1, 2984 ± 33 yr BP).

Methane production rates

Rates of CH4 production in the sediments of BTL varied between a minimum of 3.5 and maximum of 24.3 nmol cm−3 d−1 (Figure 4), at the mean in‐situ temperature (4°C). In the long‐term incubations, the highest CH4 production rates were observed in lithological Unit 2 (16 nmol cm−3 d−1), which was the organic‐rich, peaty layer and the lowest were observed in Unit 1 (5 nmol cm−3 d−1) which was the silt with low organic carbon content at the bottom of the profile. There were no significant differences in CH4 production rates of between site S and site D, suggesting CH4 production is fairly homogeneous laterally.

FIGURE 4

CH4 production in BTL sediment incubations, collected from site S (a) and site D (b), incubated at the mean in‐situ temperature (4°C). White circles denote the results of the long‐term CH4 accumulation experiments, whereas the diamonds denote the short‐term 14C tracer experiments. The error bar of long‐term CH4 accumulation is equal to one standard deviation (n = 3 separate incubations). Background colors represent the lithological units described in Figure 2

In the radiotracer short‐term incubations, the highest CH4 production rate was observed in Unit 4 (24 nmol cm−3 d−1, Figure 4). The rates of the lower lithological units showed relatively minor variations, while Unit 3 showed slightly lower CH4 production rates (mean values of 4.4, 7.2, and 8.8 nmol cm−3 d−1 for Units 3, 2, and 1, respectively). The values for the 14C short‐term incubation in Unit 4 were higher than the values for the long‐term incubations. Between the long‐term incubation and the short‐term 14C incubations, only minor differences were observed in the overall CH4 production. This was evident in the very similar mean CH4 production rates of the two methods for the 1 m sediment section (11 ± 3.5 and 9 ± 3.6 nmol cm−3 d−1, in the long‐term incubations and radiotracer short‐term incubations, respectively).

Methane production pathways

The total CH4 production measured with the short‐term tracer incubations were broken down into the acetate and CO2 components. First, the acetate turnover is discussed. In site D, the acetate oxidation rate (AOR) varied between 1.4 and 10.7 nmol cm−3 d−1 with the highest values obtained for Unit 4 (Figure 5). There were considerable differences in the AOR within the other units with a decreasing AOR trend in Unit 3 followed by an increase in AOR in Unit 2. However, there were no remarkable differences in AOR between Units 3–1 (Figure 5). Acetate fermentation rate (AFR) varied between 0.7 and 26.7 nmol cm−3 d−1 (Figure 5). As with AOR, a maximum was clearly evident in Unit 4 which may be driven by the presence of fresh labile OM input from the water column. Highly variable rates were observed in Unit 3 and 2. There were no significantly higher AFR in Units 3–1 (Figure 5). The higher TOC content of Unit 2 compared to Units 1 and 3 was therefore not a significant factor in controlling the conversion of acetate to CH4 or CO2 when rates were normalized to the volume of sediment. When the rates were expressed per gram of dry weight of soil, there were significant variations in production between the units (Table S2), but these were induced by the large density differences between the units that are made of either peat or silt (Table S2).

FIGURE 5

Depth‐dependent variations in the rate of OM degradation processes at site D. (a) White triangles are the acetate oxidation rate (AOR), white squares are the acetate fermentation rate (AFR). (b) White triangles are CO2 reduction rate (CO2RR). (c) White circles are the fraction of CO2RR relative to the total CH4 production (CO2RR + AFR), black dashed line is the expected fraction of CO2 reduction based on cellulose degradation to CH4 and CO2 in methanogenic environments. Background colors represent the lithological units described in Figure 2

In parallel to the 14C labeled acetate incubations, replicate samples from each depth were also used to measure the CO2 reduction rate (CO2RR) with labeled 14CO2. The CO2RR were between 0.4 and 6.7 nmol cm−3 d−1 (Figure 5). The lowest CO2RR was measured in Unit 4, where CO2 reduction accounted for only about 1% of the CH4 produced. In Units 3 and 2, the CO2RR were comparable but variable. Unit 1 had the highest CO2RR in the profile. Interestingly, the AFR and CO2RR appeared to balance each other to yield roughly even CH4 production rates throughout Units 3, 2, and 1 (Figure 5). Since AFR and CO2RR estimates were obtained independently and from different bottles, there were no possible experimental biases or within sample feedback that could account for this observation.

DISCUSSION

Abrupt thaw and thermokarst activity of BTL

According to historical aerial photography, BTL formation started between 1949 and 1967 (Walter Anthony et al., 2021). Present‐day thermokarst expansion of BTL suggests the following stages of lake formation that fit formation sequence of thermokarst lakes (Figure 1):

Stage 1: Prior to lake formation, the land was a forest/bog ecosystem with its ground ice still intact. The initial pre‐thermokarst forest floor was estimated to be about 6.15 ± 0.72 m above the present day lake surface, based on aerial photography on the earliest lake formation date (72 yr BP). Mean summer precipitation of 18.5 cm yr−1 thaws the active layer ~0.7 ± 0.1 cm yr−1 in nearby permafrost (Douglas et al., 2020b). The forest floor surrounding BTL is elevated by ~4.25 m above the present‐day lake surface level, based on high‐resolution LiDAR records (Hubbard et al., 2011), composed of intact Pleistocene‐aged yedoma deposits. The seasonal active layer depth from in‐situ probing shows the upland surrounding Big Trail Lake today is 0.6–1.55 m, indicating that prior to lake formation, seasonal thaw extended through the organic‐rich forest soils into the underlying silt (Unit 1; Figure 1), potentially thawing the soils and subjecting them to mineralization.

Stage 2: Ground ice melt and the forest floor collapse (i.e., thermokarst) lowered the ground surface elevation by 1.9 ± 0.72 m, bringing the ground surface closer to, but still above the lake water table. In this thermokarst terrace environment, peat‐accumulating wetland vegetation (e.g., sedges and herbs) dominated. At the lake–terrace interface, sedge transitions into floating mat vegetation, dominated by Typhus spp. These wetland plants rapidly accumulated peat on top of the eroded forest floor, constituting Lithological Unit 2. This unit of peat is also evident as a distinct drop in magnetic susceptibility (Figure S3).

Stage 3: Further ground ice melt beneath the terrace depresses the ground surface further, leading to flooding of the peat layer. In this shallow, near‐shore environment, silt washes into the lake from thermokarst and rain in the watershed. Accumulation of organic‐rich silt on top of the peat formed Lithological Unit 3.

Stage 4: As ground ice melts, it expands the lake margin laterally and depresses the sediment surface at the point of interest (i.e., sediment core location) further beneath the lake water column, well‐decomposed benthic plants and organisms contribute more heavily to the mineral‐rich gyttja (Lithological Unit 4). Within this depositional context of BTL, the controls on CH4 production were investigated.

Isotopic evidence for the partitioning of CH4 production between shallow sediments and deep taberal sediments

The incubation δ13CCH4 contrasts with the more negative δ13CCH4 of the porewater profiles at the two different sampling sites. This shows that the porewater CH4 is not entirely produced in‐situ. Correspondingly, the δ13CCH4 in the bubbles that migrate from deeper within the talik have more negative δ13CCH4 (−63.8 ± 6.3 ‰, Figure 3, Table S1). The upwards migration of deep‐generated CH4 into the upper meter of sediment appears to be a significant portion of the porewater CH4 that was measured. This overall trend is consistent with observations of more negative δ13CCH4 deeper in sediments, driven by the progressive increase in the prevalence of methanogenesis by CO2 reduction at the expense of acetate fermentation (Galand et al., 2010; Hornibrook et al., 2000; Penning et al., 2005).

To estimate the contribution of two CH4 sources (in‐situ CH4 production and CH4 transported from below) to the porewater profile, two endmember δ13CCH4 are defined: (1) CH4 produced in‐situ and (2) CH4 migrating from deeper in the sediment. Endmember (1) is already measured since it corresponds to the δ13CCH4 from each incubation experiment (Figure 3). An underlying assumption is that the incubation conditions accurately reflect the CH4 production and the resulting isotopic signature that was measured (Hines et al., 2008; Hodgkins et al., 2015). Endmember (2) is taken from the average δ13CCH4 of bubble measurements (−63.8 ± 6.3 ‰) (Figure 3, Table S1). With these endmembers and the porewater δ13CCH4 profile, the relative contribution of deep CH4 in porewater can be estimated. In both sites S and D, the fraction of CH4 produced in‐situ accounts for about 50% of the CH4, below 50 cm and 60 cm, respectively. However, above these depths, the fraction of CH4 produced in‐situ appears to diverge. This is because the difference between the δ13CCH4 of the porewater profile and the incubation becomes greater at site D but lesser at site S. The convergence of the δ13C of porewater CH4 and incubation CH4 in site S is consistent with the observation that in the first 10 cm of the profile, CH4 is not supersaturated.

Measurements of naturally present 14C abundance also show strong support for the δ13CCH4 interpretations given above. The similarity among 14C ages of the CH4 produced in incubation and the TOC of the incubated sediment (Table 1) indicates that during methanogenesis in the sediments, the CH4 comes directly from the breakdown of OM found in‐situ. The CH4 in the sampled surface bubbles reflects older 14C ages, which are in the range of previously reported CH4 seeps in Big Trail Lake (Walter Anthony et al., 2021) and other yedoma thermokarst lakes (Walter Anthony et al., 2016; Walter et al., 2008; Zimov et al., 1997), likely composed of a mixture of modern and Pleistocene‐aged carbon. Since the specific age of the ancient carbon endmember in BTL is unknown and likely heterogeneous, mixing fractions include a large uncertainty. Under the assumption that the ancient carbon source is 35,000 yr BP (Péwé, 1975), and the young CH4 endmember equal to the CH4 produced in the incubation (3300 yr BP), the fraction of ancient CH4 in the bubbles ranges between ~30% and 60%, in both sites.

Since the assumed endmember is of Pleistocene age, this calculation is conservative. Under a younger endmember of 13,650 yr BP, corresponding to the oldest CH4 being emanated from hotspots (Walter Anthony et al., 2021), the fraction of ancient carbon CH4 in the bubbles would be 40–80%. Although the natural settings may be more complex and include OM from various ages along the talik profile, the 14C results clearly show that migration of CH4 from the deeper talik (>1 m) contributes an important portion of CH4 emissions. This finding supports the calculations based on stable isotope mixing estimates, which suggests that 50% of the CH4 in the top meter porewater migrated from deeper layers.

However, the fraction of deep CH4 that equilibrates with the porewater CH4 in the upper meter of talik is unknown. This gap limits the interpretation of the isotopic evidence, since it is probable that some of the deep CH4 does not efficiently equilibrate with porewater CH4 in the upper meter before it escapes to the water column and atmosphere. For example, when bubbles from deeper in the talik are transported quickly through ebullition hotspots. Therefore, the calculations of deep CH4 contribution to the upper meter porewater need to be taken as the minimum value. It constrains that at least 50% of the CH4 emission are from deeper within the talik (>1 m). A greater fraction of lake emissions may come from the deep talik via ebullition. Bubble surveys have shown that in yedoma thermokarst lakes ebullition comprised 89% of total lake CH4 emissions (64% from deep talik and 25% from shallow sediments) while diffusive fluxes were 11% (Sepulveda‐Jauregui et al., 2015; Walter et al., 2006).

The in‐situ CH4 production rates

The alignment of in‐situ CH4 production rates in the top 1 m with two independent approaches provides robust constraint on microbial activity in this system. This is important data because CH4 production rates are reported with three orders of magnitude variations between different studies (de Jong et al., 2018; Heslop et al., 2015). This has critical implications for the biogeochemical cycling in the sediments and the resulting CH4 emissions. Thus, employing two different measurement methods alleviates methodological concerns over CH4 production reported from long‐term incubations (Dannenberg et al., 1997; Ionescu et al., 2015), as the short‐term incubations showed the same rates as the long‐term ones (Figure 4). The use of 14C radiotracer enables an independent measurement of CH4 production rates. It relies on a short‐term (hours to days) incubation rather than CH4 accumulation over many months and measures the “instantaneous” production.

In addition, the nature of the two measurements also removes the possibility that a significant amount of the CH4 produced during the long‐term incubation could be unknowingly oxidized in‐situ. This is because the shorter length of the 14C tracer incubations does not allow for tracer back conversion and thus is not influenced by methane oxidation. Hence, although they are only from one study site, the CH4 production rates presented via the dual‐incubation approach in this study are, thus far, the most reliable estimates of CH4 production in reworked Yedoma‐derived sediment and submerged Holocene boreal forest soils.

The CH4 production rates align with previous studies on thermokarst lake sediments from the Goldstream formation (Martinez‐Cruz et al., 2018), but are lower than a report in which deeper taberal sediments of the talik were included (Heslop et al., 2015). The rates also align with CH4 production of thawed late Pleistocene ice complex deposits (yedoma) and Holocene riverine sediments incubated at 4°C from the Lena River delta (Knoblauch et al., 2018).

Other substrates than acetate and CO2 can also be degraded to CH4 (such as CO, alcohols, or trimethylamines), and may be important CH4 formation pathways in peat soils (Hippe et al., 1979). However, given the similarity of the CH4 production rates between the long‐term and short‐term incubations, it is likely that most of the CH4 production is reliably predicted by considering only acetate fermentation and CO2 reduction.

The pathway of OM degradation is controlled by lithology

The classic reaction stoichiometry of OM degradation in methanogenic environments leads to ⅓ of the CH4 being generated by CO2 reduction, and the rest being driven by acetate fermentation. This is indeed the case for anaerobic digestion (Smith & Mah, 1966) and OM from rice fields (Takai, 1970). Newer studies have demonstrated more complex pathways of carbon transformations that allow deviation from the 1/3:2/3 ratio (Beulig et al., 2018). In BTL, the pathway stoichiometry appears to be respected in Units 2 and 4 where CO2 reduction accounts for about ⅓ of the total CH4 production or less (Figure 5c). In Units 1 and 3, CO2 reduction rate is distinctly higher than ⅓ (Figure 5c).

High Fe2+ concentrations in the sediments of BTL suggest iron reduction in the talik. In this case, reduction of Fe3+ is likely coupled to the oxidation of acetate by microbial metabolisms. This reaction could partially account for the acetate oxidation observed. However, iron reduction rates based on the accumulation of Fe2+ are considerably lower than acetate oxidation rate in BTL, between 0.4 and 2 nmol cm−3 d−1 (Figure S6). Considering the stoichiometric ratio between iron reduction and acetate oxidation (8:1, Eq. S5), iron reduction can only be responsible for a small portion of the observed acetate oxidation (Figure S6). Other electron acceptors such as nitrate and sulfate are only found in low concentrations in BTL (Figure S1, Table S3). Oxidation of acetate by alternative electron acceptors is therefore an unlikely explanation for the measured acetate oxidation rate. The increase in Fe2+ concentrations during the incubations may have been underestimated due to precipitation, sorption of Fe2+ to Fe3+ or peat, but not because of reaction with HS− since it was below detection limit (Table S3).

In several lake sediments, CO2 reduction has been observed to contribute >50% of total CH4 production (Conrad et al., 2011; Ji et al., 2016). The measurements from Unit 1 and 3 are consistent with these observations. However, it is still unclear whether CO2 reduction in the sediments of BTL is due to partial OM degradation that leaves the remaining OM more oxidized, or whether syntrophic acetate oxidation is more important in these units. In the methanogenic zone of the Baltic seafloor, OM has been degrading under anaerobic conditions for thousands of years. There, CH4 production is proceeding almost entirely by CO2 reduction, while the reducing power (H2) needed to drive CO2 reduction appears to be supplied entirely by syntrophic acetate oxidation (Beulig et al., 2018). This is not the case in Units 1 and 3 of BTL. There, the CO2 reduction pathway accounts for a significantly smaller fraction of CH4 production than in the Baltic Sea. Yet, syntrophic acetate oxidation still may be contributing some reducing power to CO2 reduction. If one considers the conventional stoichiometry of the breakdown of OM, the reducing equivalents released by the degradation of cellulose (see introduction) in Units 1 and 3 are not enough to account for the measured CO2 reduction rates. This requires additional reducing equivalents, by incomplete breakdown of recalcitrant types of OM.

In BTL, Units 1 and 3 are Yedoma‐derived sediments of Pleistocene age and consist of remobilized loess (Péwé, 1975). The OM of Units 1 and 3 has undergone a different degradation history than Unit 2; thus, their degradation history may account for the greater fraction of CO2 reduction. For example, in temperate wetland soils, there is a clear relationship between the relative importance of CO2 reduction and older, less reactive sediments (Hornibrook et al., 1997). In sediments of the Goldstream formation, it appears that in Yedoma‐derived sediments of Pleistocene age, CO2 reduction is more prominent than expected in fresh labile OM (Smith & Mah, 1966). Partial degradation of the OM may have occurred in the past thousands of years and affected OM differently depending on depositional history (Heslop et al., 2017). In BTL, climate‐induced variations in active layer depth could have enabled partial degradation of the OM in Unit 1.

Interestingly, these findings offer insight into the δ13CCH4 signatures that were measured. During incubation of Units 1 and 3, the δ13CCH4 is more positive than the δ13C of the CH4 produced deeper within the talik. This suggests that an even greater fraction of the CH4 produced deeper in the talik may be produced by CO2 reduction. Since the pathway of CH4 formation did not influence the total CH4 production in the 14C short‐term incubations and only marginally in the long‐term incubations (Figure 4), little to no direct relationship is established between pathway and overall net CH4 production.

Integrated CH4 production rates are comparable with hotspot emissions from BTL

The top meter of sediment produces up to 35% of the CH4 measured ebullition in BTL. This estimate was obtained by comparing point‐source hotspot seep ebullition of 6.5 mol CH4 m−2 yr−1 on average across the lake (Walter Anthony et al., 2021) with the mean production of Units 3–1, 6.3 ± 2.1 nmol cm−3 d−1 (2.3 ± 0.77 mol m−3 yr−1). Unit 4 was not included because it contains fresh labile OM that is representative of recent lake deposition. The ebullition‐based emission estimate may be high since the value reported by Walter Anthony et al. (2021) was based on bubble samples collected in winter, which appear to have a higher CH4 content than bubble samples collected in summer. On the other hand, estimates based on CH4 ebullition are higher than CH4 emissions because bubbles trapped by lake ice in winter release CH4 to the water column which gets oxidized (Greene et al., 2014). Furthermore, the mean CH4 flux from the lake's eddy covariance tower over 3 years is about 25% lower than estimates based on ebullition (Walter Anthony, unpublished data). The comparison between CH4 ebullition and CH4 production also implies that up to 65% of the CH4 must be from deeper sediments (>1 m). Isotopic evidence discussed earlier in this study constrained the origin of at least 50% of the CH4 as below the top meter and reinforces the assessments. Big Trail Lake has 200 mapped ebullition hotspot seeps that may facilitate CH4 transport from deep within the talik.

It is also possible to compare emissions with CH4 production rates by taking into account the thickness of the thawed talik. Helicopter‐born electromagnetic (resistivity) profiles estimate the depth of permafrost as the depth when resistivity falls below 250 Ohm‐m, which was found to be less than 15 m but most likely over 10 m at both site S and site D (Walter Anthony et al., 2021). Taking these talik depths into consideration and assuming that the CH4 production rates down to the permafrost depth are the same as the average of Units 3–1, the areal CH4 production in BTL would be between 15 and 47 mol CH4 m−2 yr−1. The calculated flux depends on the integration interval, which represents the plausible range of talik depth (10–15 m) (Figure 6).

FIGURE 6

CH4 emissions from BTL (mol m−2 yr−1) are obtained by integrating the CH4 production estimates in the first meter of BTL over the potential range of talik depths (red boxes). The middle of the box corresponds to the mean estimate and the length of the boxes corresponds to 1σ confidence interval. The blue box and whiskers plots are CH4 emission estimates from thermokarst lakes in Alaska (Sepulveda‐Jauregui et al., 2015) and Siberia (Serikova et al., 2019). Blue circles are the outliers of the box and whisker plots. The blue diamond is the CH4 hotspot emission from a pond north of BTL (Elder et al., 2021), and the blue star is the calculated CH4 flux of Big Trail Lake based on winter point‐source ebullition (Walter Anthony et al., 2021)

An emission hotspot located onshore along the margin of a small thermokarst pond just east of BTL (Eastside Pond) has recorded emissions with a maximum of up to 225 mol CH4 m−2 yr−1 (Elder et al., 2021), much larger than the highest production estimates (based on the highest CH4 production and maximum talik depth). This high emission is spatiotemporally limited, resulting from the buildup of CH4 and conduits that channel CH4 produced in a larger volume to a point source of emission (Langenegger et al., 2019). The mean hotspot CH4 flux from this particular hotspot is of 23 mol CH4 m−2 yr−1 (Elder et al., 2021) which fits within the production estimates when assuming a talik depth of 15 m (Figure 6). Thus, one important finding of this study is that it is not necessary to invoke anomalously high CH4 production rates or unique hydrological conditions, to explain high CH4 emissions from deep taberal sediments. Hotspot emissions are the result of the activation of microbial fermentation with rates that resemble those presented in this study, only over large depths of recently thawed permafrost. They make up 0.01% of permafrost land area but account for 4% of the pan Arctic budget, 1.1 Tg CH4 yr−1 (Elder et al., 2021).

The CH4 emission estimates are however less than the CH4 production rates measured in this study when they are combined with the expected talik depth of 10–15 m (Walter Anthony et al., 2021). Indeed, the talik depth needed to produce the average emission of BTL is about 4 m (Figure 6). There are three scenarios that could explain this discrepancy. (1) CH4 production rates are lower in the taberal sediments than those measured in Units 1–3, (2) CH4 oxidation may consume a significant amount of the CH4 production, or (3) the incubations yield higher rates than in‐situ rates.

In scenario (1), the high CH4 production rates observed in the first meter of BTL sediments may be caused by small amounts of fresh OM deposited from the water column and Unit 4 leaching into Units 3–1, thus driving higher CH4 production rates in the topmost meter than in lower taberal sediments. A power‐law decline in OM mineralization rate (Middelburg, 1989) would require very little fresh OM input to produce relatively high rates of CH4 production. The explanation would be consistent with previous observations in Vault Lake showing that high CH4 production rates in uppermost “organic‐rich mud” facies continued after a significant drop in TOC, whereas deep taberal sediment had markedly lower CH4 production rates (Heslop et al., 2015). There is significant evidence that deep taberal sediments have lower CH4 production rates. Ebullition hotspots, which channel deep CH4 migration to the surface show a stronger 13C‐depletion than the CH4 produced in the topmost meter (Figure 3, Table S1), likely because of slower metabolic rates (Valentine et al., 2004). Furthermore, it has been recently shown that in thermokarst settings, older CH4 is produced under slower metabolic rates (Douglas et al., 2020a). Since carbon age is often correlated with depth, this finding reinforces the possibility that CH4 production rates are lower in the deep taberal sediments.

Scenario (1) could also be explained by a slow onset (years) of CH4 production after thaw (Knoblauch et al., 2018). Fe3+ and SO4 2− first need to be consumed for methanogens to thrive (Patzner et al., 2020). Then, the resident microbial community capable of producing CH4 needs to be established. However, the metabolic capacity to produce CH4 may be highly dependent on the depositional history and environmental conditions pre‐thaw (Heslop et al., 2019; Holm et al., 2020). In marine sediments, it has been shown that the taxa exploiting geochemical gradients are not controlled by migration of favorable taxa to the gradients but rather, by the persistence of the initial microbial community deposited at the time of sedimentation (Starnawski et al., 2017). This could significantly inhibit CH4 production in the post thaw taberal sediments. In incubation experiments with thawed permafrost soils from the Lena Delta, the delay between thaw and the onset of CH4 production is sometimes 8 years and in some cases may never occur (Knoblauch et al., 2018). Future studies should carefully study the controls on microbial colonization, growth, and CH4 production in the deep talik.

Scenario (2) considers that CH4 production rates in the deep taberal sediments are similar to the rates measured in the upper sediments. In this case, CH4 oxidation is needed to explain the lower emissions. Methane oxidation was documented within the sediments (Winkel et al., 2019) and the water column of thermokarst lakes (Thalasso et al., 2020). Yet, there is considerable uncertainty about its magnitude. For example, in Vault lake in Alaska, reconciling production estimates with emissions estimates suggest that the vast majority of the CH4 produced may be oxidized aerobically or/and anaerobically before reaching the atmosphere (98.6%) (Heslop et al., 2015). However, in the same lake through an isotopic approach, CH4 oxidation appears to consume 41–83% of production (Winkel et al., 2019). In Goldstream lake located near BTL, around 50% of the CH4 is oxidized under the ice during winter before it has time to escape (Elder et al., 2019; Greene et al., 2014). In BTL, if CH4 production estimates from this study are continuous down to the permafrost, balancing the budget of production with the average emissions requires that ~40–80% of CH4 production is consumed via CH4 oxidation before it reaches the atmosphere (Figure 6). These estimates eliminate the possibility that the vast majority of CH4 production is oxidized before it reaches the surface.

Scenario (3) takes into consideration that despite the considerable effort expended in maintaining the environmental conditions during incubations, it is not possible to maintain conditions exactly like in the sediment. Any type of sediment handling seems to increase rates of organic carbon degradation (Kan et al., 2021; Michaud et al., 2020). In addition, the energetic yield of CH4 production is influenced by the concentration of products and reactants in the sediment (e.g., Beer & Blodau, 2007). As the process occur near the thermodynamic limits, it is that thermodynamics might affect the rates. In the incubations, DIC did not vary significantly from the in‐situ concentration and acetate was low throughout. However, CH4 concentrations in‐situ were at saturation level (Figure S2) while in the incubations, they started near zero. This difference could, in turn, affect methanogenesis, competition between different pathways of methanogenesis or competition for organic substrates by different microbial metabolisms. Future experiments could investigate the effect of the loss of CH4 during sampling on the CH4 production rates. In the broader goal of understanding CH4 emissions from quickly degrading permafrost, the short‐term incubation and the long‐term incubations are an inexpensive, scalable, and reliable way to constrain CH4 cycling in thermokarst lake sediments. The estimates fit within the overall budgets of CH4 emissions obtained by incubation‐independent methods (Figure 6), and provide an independent way to assess CH4 emissions from the ground up.

CONCLUSIONS

This study shows that short‐ and long‐term incubations measure very similar CH4 production rates, providing robust constraint on rates. According to these two independent methods, measured CH4 production in the top meter of the sediments of a young thermokarst lake is ~10 nmol cm−3 d−1. The redeposited silt/loess sediment is dominated by CO2 reduction while the peat appears dominated by acetate fermentation, consistent with patterns previously observed through field measurements of ebullition in Siberian and Alaskan Yedoma thermokarst lakes (Walter et al., 2008). Yet, incubations in the peat layers did not produce CH4 with a different isotopic composition than the silt/loess layers, which were more dominated by CO2 reduction. However, it was clear that the CH4 dissolved in the porewater at the time of coring was more 13C depleted than the CH4 produced during incubations, likely because part of the CH4 dissolved in porewater is produced in different conditions than those encountered in the top meter.

The extrapolated CH4 production aligned well with ebullition estimates from BTL (Figure 6). Thus, CH4 production rates of the magnitude measured in this study can comfortably explain high CH4 emissions in deep taberal sediments without needing to invoke extreme or unique situations. BTL is simply a deep, thawed OM reservoir. Future studies should assess CH4 production throughout the entire depth of thaw because inherently, they could reveal important information on the mechanisms of CH4 production and the root of future CH4 emissions as the permafrost thaws and the talik expands. Finally, quantifying CH4 oxidation is critical to constraining sedimentary CH4 cycling and CH4 production rates in deep taberal sediments.

Supplementary Material

Click here for additional data file.

Supplementary Material

Click here for additional data file.