Minor contribution of small thaw ponds to the pools of carbon and methane in the inland waters of the permafrost-affected part of the Western Siberian Lowland

The WSL is comprised of taiga, forest-tundra, and tundra biomes. Permafrost is abundant north of 60°N, changing from sporadic to discontinuous and continuous northward (figure 1). The climate is humid semi-continental with mean annual temperatures ranging from −2.8 °C in the south of the cryolithozone (Surgut region) to −10.4 °C in the north (Gyda Peninsula). The annual precipitation ranges from 600 mm in the south to 350 mm in the north. Peat thickness ranges from 1 to 3 m and is not directly linked to the latitude (Kremenetski et al 2003). The bogs of the region were subjected to freezing during the Subboreal period (~ 4500 y.a). After the increase in temperature and precipitation during the Subatlantic period (2500 y.a.), the thermokarst started (Ponomareva et al 2012, Pastukhov et al 2016). It produced thaw ponds and thermokarst lakes that were formed due to progressive thawing of frozen peat rich in ice (Pokrovsky et al 2011, Grosse et al 2016). On flat watershed divides, which dominate the territory (palsa peat bog and polygonal tundra), the lakes are subjected to a natural cycle of growth and draining which is well established in the WSL (Kirpotin et al 2011) but is also evidenced in other parts of the circum-Arctic zone (Hopkins 1949, French 1996). Because of the homogeneous nature of the surrounding substrate (Quaternary clays and sands of fluvio-glacial and lacustro-glacial origin underlying the surface peat deposits), the thermokarst lakes are generally round in shape. As such, all the lakes and thaw ponds inventoried in this work are considered as having a thermokarst origin. While satellite imagery cannot resolve the origin of small ponds, our terrestrial observations in addition to the very high number of small ponds suggest these ponds were formed after thermokarst subsidence rather than remaining in the basin following large lake drainage.

The range of thermokarst lake and thaw pond sizes in the permafrost-affected part of the WSL is from a few m² to hundreds of km² (Manasypov et al 2014). For convenience, we will distinguish between large lakes (>10 000 m² or 1 ha), small lakes (500–10 000 m²), and very small thaw ponds (< 500 m²). This conventional distinction is based on hydrochemical and gas regimes of lakes and ponds; there is an abrupt rise of DOC, CO₂, and CH₄ in water bodies smaller than 500 m² (Pokrovsky et al 2011, 2013, Shirokova et al 2013, Manasypov et al 2014, Polishchuk et al 2015b). This classification is also consistent with that of Holgerson and Raymond (2016) who considered lakes < 1000 m² as very small water bodies. Muster et al (2017) followed the distinction of Rautio et al (2011) and defined ponds as bodies of standing water with a surface area smaller than 10 000 m².

Water samples were collected from a PVC boat for large lakes or directly from the lake center for small (< 50 m diameter) water bodies (thaw ponds) from the middle of the water column. The sampling was performed during July and August from 2010–2014 (Pokrovsky et al 2013, 2014, 2016, Shirokova et al 2013, Manasypov et al 2014, 2015, Pavlova et al 2016, Loiko et al 2017) and in August 2015 and 2016 (this study). The typical depth of sampling was between 0.5 and 0.25 m depending on the lake depth (see section 3.3 below). A full list of collected lakes and relevant physico-geographical parameters together with the measured GHG and DOC is given in table S1 of the supplementary information available at stacks.iop.org/ERL/13/045002/mmedia. For GHG analyses, unfiltered water was sampled in 25 m L serum bottles that were closed without air bubbles using vinyl stoppers and aluminum caps and immediately poisoned by adding 0.2 ml of saturated HgCl₂ using a two-way needle system. In the laboratory, a headspace was created by displacing approximately 40% of water with N₂ (99.999%). Two 0.5 m L replicates of the equilibrated headspace were analyzed for their concentrations of CH₄ and CO₂, using a Bruker GC-456 gas chromatograph equipped with flame ionization and thermal conductivity detectors. After every ten samples, a calibration of the detectors was performed using air liquid gas standards (i.e. CH₄ = 145 ppmv and CO₂ = 3000 ppmv). Duplicate injections of the samples showed that results were reproducible within ±5%. The specific gas solubility for CH₄ and CO₂ (Weiss 1974) were used in the calculation of the total CH₄ and CO₂ content in the vials, which was then recalculated to µmol/L of the initial waters.

For DOC analyses, collected waters were filtered onsite in pre-washed 30 mL polypropylene Nalgene® bottles through single-use Minisart filter units (0.45 µm pore size, Sartorius, acetate cellulose filter). The first 20 ml of filtrate were discarded. Blanks were used to control the level of contamination induced by sampling and filtration. DOC blanks of filtrate never exceeded 0.1 mg L⁻¹ which is quite low for the organic-rich pore waters sampled in this study (i.e. 10–100 mg L⁻¹ DOC). DOC was analyzed using a Carbon Total Analyzer (Shimadzu TOC VSCN) with an uncertainty better than 3% (see Prokushkin et al 2011 for methodology).

Satellite imagery from the Landsat-8 Operational Land Imager (30 m MR) [available at the USGS Global Visualization Viewer: http://glovis.usgs.gov] together with HR (2 m)) Kanopus-V scenes (available at the Roscosmos website: http://eng.ntsomz.ru/) were used to map the lake distribution over the full territory of permafrost-affected WSL (figure 1). We used MR Landsat-8 images collected at the end of July–beginning of August 2013–2014. HR Kanopus-V images were collected during summers 2013–2015. This free-ice period corresponds to the minimal coverage of the territory by lakes and minimal seasonal variation of lake water levels. A combination of images from the two years was necessary because a single year's observation could not provide full coverage of the territory.

For the entire permafrost-affected WSL territory of 1.05 million km², we used a superposition of 134 MR images. Images were treated using standard tools of ArcGIS 10.3 software (Kennedy et al 2013). For automatic identification of lakes we used the Fmask algorithm (Zhu et al 2015). It allows resolving the lakes even with some cloud coverage and uses several spectral ranges of Landsat-8 using a set of empirical spectral indexes. As input data we used the intensity of the upper atmosphere reflectance and the brightness temperature presented as nominal units of pixel brightness. First, for the mosaic of Landsat-8 imagery, cloud masks were defined for individual images. Then the cloud masks were removed from the images and replaced by cloud-free fragments taken from the adjacent images during the same period of observation. The minimal pixel on the MR image was 30 × 30 m², thus limiting the minimal lake size of 0.5 ha or 5000 m². This corresponds to six pixels of 30 m × 30 m (900 m²), which is sufficient for reliably distinguishing lakes from background noise.

A two-stage procedure was developed to separate thermokarst lakes and thaw ponds from other water objects. First, we masked coastal zones and rivers, based on vector layers of OpenStreetMap and State Water Register (http://textual.ru/gvr/). Second, we superposed this mask on a raster layer of recognized water objects forming a layer of lakes. The uncertainty of lake surface area measurement using Landsat-8 for lakes with sizes larger than 10⁵ m² is 3.5% and for lakes larger than 2 × 10⁴ m² is 7% (Bryksina and Polishchuk 2013, Kornienko 2014).

It is important to note that thermokarst lakes and thaw ponds of western Siberia are different from glacial lakes of other boreal and subarctic regions. The latter are often developed on moraine and crystalline rocks and located within non-isometric elongated basins that were formed by glacier movement. In contrast, thermokarst lakes in palsa bogs of western Siberia occurring within frozen 1–3 m peat deposits are rarely oval and are usually round and isolated (Kirpotin et al 2009, 2011, 2014, Shirokova et al 2013, Pokrovsky et al 2014). According to our field observations and detailed lake mapping over 950 km latitudinal gradient of western Siberia, the fraction of lakes with an irregular shape is less than 10% (Pokrovsky et al 2011, Manasypov et al 2014, 2015, Polishchuk et al 2017). For this reason, the influence of the irregularity of the lake border on the determination of the total lake area can be neglected. This is also confirmed by a study of the tortuosity of thermokarst lake border lines demonstrating that the deviation from circumference in lakes of sporadic, discontinuous and continuous permafrost zones produces less than 5% uncertainty on total lake area evaluation by Landsat-8 (Polishchuk and Polishchuk 2012).

HR Kanopus-V images were used for the quantification of small lakes and very small thaw ponds. For the reliable identification of lakes, at least five pixels are necessary. The minimal pond size inferred from the Kanopus-V images was 20 m² as the spatial resolution of 2 m yielded a pixel size of 4 m² allowing the required five pixels to be distinguished on the 20 m² area. The scenes were processed using ArcGIS 10.3 within 78 test sites (TSs), which were evenly distributed over sporadic, discontinuous, and continuous permafrost zones of the WSL (figure 1). The characteristics of these TSs are given in table 1. The number of lakes and ponds on each TS ranged from 97–2583. All TSs had similar sizes (30 km²) and covered approximately 2000 km² which represents 0.002% of the overall permafrost-affected territory of the WSL.

We used a binary classification algorithm in ArcGIS 10.3 analogous to the 'density slice' method used by Muster et al (2017). This algorithm is based on the visual detection of certain threshold values of spectral brightness that allow us to distinguish two object classes: 'water' and 'not water'. The threshold value had to be selected individually for each TS due to differences in illumination, acquisition geometry, and different sensor spectroradiometry. The panchromatic imagery exhibited a strong contrast between water surface and the surrounding vegetation which allowed an adequate selection of threshold values needed for classification. The submerged moss and macrophytes were not detected and non-submerged moss and plants along the lake coasts had no influence on lake water surface area.

To build the diagrams of lake distribution in a wide range of lake sizes, we used 21 partial ranges of lake size taken within a logarithmic scale, i.e. 20–50 m², 50–100 m², 100–200 m², 200–500 m² and so on until 200 000 000 m². To combine MR and HR images we developed a three-step procedure. First, we constructed diagrams for the number and area distribution of all thermokarst lakes detectable on the MR images. Then, we built similar diagrams for small thermokarst lakes and very small thaw ponds using HR images of the 78 TSs across the permafrost-affected part of WSL. At this stage, the total number and area of lakes located within the TSs were extrapolated to the full territory of the permafrost-affected part of WSL, as described below.

Considering the values of lake numbers and areas in the lake size intervals suitable for Landsat-8 mapping as the most reliable, we required that within the lake size intervals where Landsat-8 and Kanopus-V images overlap, the lake number and area determined from HR scenes on TSs were to be equal to values extrapolated to the full studied territory. We believe that such extrapolation is valid because the TSs are evenly distributed across the permafrost-affected part of the WSL. The extrapolated (calculated) values for the total number (_Ncj) and total lake area (_Scj) in each jth interval of HR images-based histogram were obtained from:

$$\begin{equation} {N_{cj}} = {N_{ej}} \times K \end{equation} \tag{ 1 }$$

$$\begin{equation} {S_{cj}} = {S_{ej}} \times K \end{equation} \tag{ 2 }$$

where N_ej and S_ej are the total number and area of lakes, respectively, that are experimentally determined based on HR scenes of all TSs, j is the number of intervals of the histogram (figure 2), and K is the empirical extrapolation coefficient which was averaged over all lake sizes and therefore independent in terms of size range. Finally, we combined primary MR and HR histograms into integral histograms of lake numbers and overall areas. These histograms included lakes and small and very small thaw ponds spread over seven orders of magnitude in size.

Zoom In Zoom Out Reset image size

The uncertainty of extrapolating data from the 78 TSs to the whole territory of the permafrost-affected part of the WSL was evaluated via a comparison between data of HR and MR in the lake size interval where the HR and MR scenes overlap (figure 2). Using the number of lakes in the seven intervals of superposition, we found the relative difference in lake numbers using equation 3:

$$\begin{equation} \beta_{i} = (m_{i}-n_{i})/n_{i} \end{equation} \tag{ 3 }$$

where m_i and n_i represent the number of lakes in the ith overlapped interval, determined using HR and MR images, respectively. In this equation, n_i is chosen as a reference number given that HR images were treated for the full territory of the permafrost-affected part of WSL. The average β_i value of the seven overlapped intervals was equal to 15% which is acceptable for lake size distribution estimation.

Zoom In Zoom Out Reset image size

The MR images allowed resolving lakes of the minimal size of 5000 m² corresponding to six pixels which is necessary to distinguish the lakes from the surrounding landscape. The total number of lakes ranging from 5000 m² to 2 × 10⁸ m² was 727 700 with a total water area of 60000 km². The histograms of lake numbers and lake area distribution over discrete size ranges are given in the blue columns in figures 2(a) and (b), respectively. The HR (Kanopus-V) scenes were used to quantify the number and areas for small thaw ponds. The empirical extrapolation coefficient (equations 1 and 2) was calculated as an average of the ratio of lake number and area determined by MR scenes for the entire permafrost-affected part of the WSL territory to the total lake number and area from all HR images (at all TSs) in each lake size interval where the HR and MR histograms overlap. The histograms have a seven-interval wide overlap range from 5000 to 10⁶ m² (figure 2). The extrapolation coefficient values calculated from the data on the number and area of lakes were equal to 85 and 86 respectively. Thus, K = 85.5 was defined as an arithmetic mean of these values.

The calculated numbers and total areas of lakes in each interval of lake size using equations 1 and 2 yielded the diagrams of small lakes and very small thaw ponds distribution on the territory, shown as yellow columns in figures 2(a) and (b) for the lake size range of 20–10⁶ m². The integration (stitching) of HR and MR-based lake diagrams can be performed at the lake size of 20 000 m² with a typical uncertainty on lake determination at 7% (Bryksina and Polishchuk 2013). This point is shown in figure 2 as a red vertical line. This last step of integrated histograms yielded a lake size and area distribution between 20 and 2 × 10⁸ m² so that lakes smaller than 2 × 10⁴ m² were processed using HR scenes, and the lakes >20 000 m² were treated from MR images.

The total number of lakes was equal to 8.7 million and their total area was 6.4 × 10¹⁰ m². Very small thaw ponds (< 500 m²) provided 44.2% of total lake numbers, but their area contributed to less than 0.1% of total area coverage. Big lakes (>10 000 m²) which accounted for only 9.1% of lake numbers provided 91.1% of all lake areas.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The number of lakes and lake area–lake size distribution diagrams obtained in this work for the permafrost-affected part of WSL are in fairly good agreement with the inventory of all lakes of the circum-Arctic permafrost region (Muster et al 2017); deviation between the two estimations does not exceed several percent (figure 3). A comparison with the world lake database of Cael and Seekell (2016) is shown in figure 4 where empirical distribution lakes in Sweden and throughout the world are shown. It can be seen that for lakes smaller than 10⁴–10⁵ m², at both the planetary and regional (WSL, Sweden) level, there is a sizeable difference between the results of direct measurement and calculation from the power law of lake distribution. Specifically, the power law extrapolation strongly overestimates the small lake number compared to empirical data. As a result, the use of power law can also lead to the overestimation of water and solute amounts in small lakes and very small thaw ponds. Indeed, a plot of the relative contribution of lakes of different size to the total area of lakes demonstrates sizeable differences in the fraction of small lakes calculated by Holgerson and Raymond (2016) for the whole planet and our direct measurements of very small thaw ponds in the permafrost-affected part of the WSL territory (figure 5). The two dependences have a generally similar shape and lakes having an area of 10⁵ m² provide maximal contribution to the total areal coverage. However, significant differences are pronounced for small water bodies (< 1000 m²); the relative contribution of lakes evaluated by Holgerson and Raymond (2016) for the world is several times higher than that directly measured in western Siberia. For very small ponds of 100 m², this difference is a factor of 6.5. As such, the estimation of total concentrations and emissions of CH₄ and CO₂ from small lakes on a worldwide scale, performed based on Monte Carlo numerical modeling using the Pareto power law (Holgerson and Raymond 2016), cannot be applied to very small thaw ponds of western Siberia. In the WSL direct remote sensing using HR images from 78 TSs allows for a more rigorous lake and pond inventory and thus provides a solid basis for the evaluation of GHG and carbon pools.

Zoom In Zoom Out Reset image size

It is important to note that the use of Muster et al's (2017) Pond and Lake circum-Arctic database for the WSL territory is not possible. In the work by Polishchuk et al (2017), only two sites (Central Yamal and the Surgut region) were inventoried where lake coverage (15%–20%) was much higher than that of the total WSL permafrost-affected zone (5%–10%). In the southern part of the WSL cryolithozone, within the peatland 'Surgutskoe Polesye' of the northern taiga biome (61.5–63°N), the proportion of open water ecosystems is around 19% based on HR (1–3 m) Google Earth inventory (Terentieva et al 2016). However, the average limnicity of the permafrost-affected part of the WSL assessed in this study is consistent with values reported in a similar context for two peat plateau areas of northeast European Russia (7%–13.6%, Sjöberg et al 2013) and the recent evaluation of lake area density based on geo-statistical modeling (Messager et al 2016).

To estimate the amount of GHG and dissolved carbon in lakes of the permafrost-affected part of the WSL, published and new collected concentrations of DOC, CO₂, and CH₄ in 146 thermokarst lakes north of the permafrost boundary were used (Pokrovsky et al 2013, 2014, 2016, Shirokova et al 2013, Manasypov et al 2014, 2015, 2017, Pavlova et al 2016, Loiko et al 2017, table S1 of the supplementary information). The CO₂ and CH₄ concentrations in small (< 500 m²) thaw ponds range from 29–778 µmol CO₂ L⁻¹ and 0.02–47 µmol CH₄ L⁻¹ with average values equal to 172 and 9 µmol L⁻¹, respectively. This is significantly (p < 0.05) higher than the range of CO₂ and CH₄ in thermokarst lakes >500 m²; 29–332 with an average of 92 µmol L⁻¹ and 0.02–8.6 with an average of 0.98 µmol L⁻¹, respectively, as illustrated in figures 6(a) and (b).

Zoom In Zoom Out Reset image size

The DOC concentration in ponds is a factor of 1.5–2.0 higher than that in lakes (figure 6 (c)). Note that the ratio of CO₂ to CH₄ decreases with the decrease in water body area (figure 7) and that this is consistent with the global world trend (Holgerson and Raymond 2016). The concentrations of CO₂, CH₄, and DOC were not found to be significantly (p < 0.05) dependent on latitude across the WSL permafrost-affected zone surface waters, soil solutions, and suprapermafrost waters (Manasypov et al 2014, Raudina et al 2017). For this reason, as a first approximation, we will neglect any latitudinal trend in C concentrations in ponds and lakes and consider only lake area as the main controlling parameter of GHG and DOC concentration. It is important to note that the thermokarst lakes and thaw ponds of the WSL are extremely shallow (large lakes are 0.5–1.5 m deep and small thaw ponds are < 0.5 m deep) and are not stratified in temperature. We sampled the surface and bottom layer and have not found any measurable chemical stratification in DOC and CO₂ or O_2. Lakes were fully oxygenated over 0.5–1.5 m depth (Shirokova et al 2013, Manasypov et al 2015). The latter allowed us to assume the constant profile of CH₄ concentration in the water column.

The depth of WSL thaw ponds and thermokarst lakes depends on the size of the water body (S) and can be approximated by the polynomial equation 4:

$$\begin{equation} \begin{array}{l@{}l@{}l@{}} {\rm Depth} = 0.00007 & & \, \times \,S^{3}-0.0049 \times S^{2} + 0.1141\\ & & \times \, S - 0.0124\, (R^{2} = 0.552) \end{array} \end{equation} \tag{ 4 }$$

which is based on field measurements of thaw ponds and lakes of 0.2 to 2 × 10⁷ m² (Polishchuk et al 2017). This information provides a solid background for the estimation of GHG and DOC reservoirs in the lakes of western Siberia.

The masses of CH₄, CO₂, and DOC in WSL lakes and ponds for seven lake size intervals are listed in table 2. The relative proportions of GHG and DOC in lakes of different size are given in table 3. The maximal contribution to the CH₄ pool in the permafrost-affected part of the WSL is provided by lakes between 10⁵ and 10⁶ m² surface area (figure 8). Lakes larger than 10⁴ m² provide >90% of the total CH₄ amount whereas the share from small thaw ponds (< 1000 m²) is only 2%. According to these data, small thaw ponds contribute less to overall CO₂ and DOC levels than to CH₄ level.

Zoom In Zoom Out Reset image size

Remote sensing studies of the permafrost zone of western Siberia demonstrated that the number of newly formed small thermokarst lakes (5000–20 000 m²) over the past three decades exceeds, by a factor of 20, the number of large lakes which tend to disappear during the same period (Bryksina and Polishchuk 2015). This increase in the number of small and very small thaw ponds over the past few decades was demonstrated through an examination of HR (2 m) Kanopus scenes. Assuming the conservative behavior of carbon in WSL thermokarst lakes and ponds and ignoring the short and long-term dynamics of carbon cycling, the draining of large thermokarst lakes and the appearance of small thaw ponds suggest a sizeable decrease in CO₂, CH₄, and DOC pools in thermokarst water bodies of the WSL. Thus, a 10%–20% decrease in the area of large (>20,000 m²) lakes suggested in earlier studies (Smith et al 2005, Riordan et al 2006, Kravtzova and Bystrova 2009, Shiklomanov et al 2013) will lead to a corresponding decrease in CH₄, CO₂, and DOC concentrations in ponds and lakes of the permafrost-affected part of the WSL. At the same time, even a two to threefold increase in the number of small (< 10 000 m²) lakes and a tenfold increase in very small (< 500 m²) thaw pond number and area will not produce any sizeable (i.e. >2%) increase in GHG and DOC pools in WSL ponds and lakes. Taken together, the permafrost thaw and climate warming can only decrease the pools of GHG and organic carbon in thermokarst water bodies of western Siberia. Assuming a solely diffusive flux, which is governed by the GHG concentration gradient, this will decrease the CO₂ and CH₄ emission potential to the atmosphere. In addition, the lateral export of DOC from the lakes and ponds to the hydrological network will also decrease. In this regard, our study corroborates the previous results of van Huissteden et al (2011) and Kessler et al (2012) that thermokarst lake drainage leads to lower GHG emissions.

There are other indications that climate warming and the permafrost thaw may not positively impact the GHG emissions from lentic aquatic ecosystems. For example, the landscape-scale model that included the entire life cycle of thaw lakes in northern Siberia (initiation, expansion, drainage, and eventual re-initiation) demonstrated that methane emissions from thaw lakes in eastern Siberia (Yakutia) are an order of magnitude less alarming than previously suggested (van Huissteden et al 2011). Moreover, the use of the Integrated Global System Model demonstrated that, on the whole the Arctic and boreal permafrost-scale increases in atmospheric CH₄ and its radiative forcing (which result from the fact that thawed, inundated emission sources under saturated (anaerobic) soils) are small, particularly when weighed against human emissions (Gao et al 2013). Considering the results presented in this work on aerobic (open lake and pond) areas, we believe the overall potential threat of permafrost thaw on GHG emission from Arctic and high latitude boreal permafrost-bearing regions may be overestimated. We do note, however, that detailed mapping of western Siberian wetland complexes (i.e. Peregon et al 2008, 2009, Terentieva et al 2016) is necessary to integrate available DOC and GHG concentrations and fluxes into a global ecosystem-specific model for this territory. Here, mapping of fine-scale heterogeneity of WSL's wetland landscapes, which is currently possible only for CH₄ emission (i.e. Bohn et al 2007, Glagolev et al 2011), will constitute the main challenge for both the DOC and CO₂ inventories. In this regard, quantifying the amount of water and carbon in waterlogged hollows—especially those abundant in the southern part of the permafrost-affected WSL territory (i.e. Terentieva et al 2016)—that presumably contain high concentrations of GHG and DOC, may represent a bottleneck of C reservoir inventory in inland waters of the Siberian wetlands. Furthermore, seasonal and daily variability of GHG fluxes in the WSL (Kim et al 2011, Sabrekov et al 2014) may impede quantitative upscaling of available carbon and GHG concentrations on a year-round basis. For this, landscape-scale modeling of process-based lake biogeochemistry (i.e. Tan et al 2017) is needed.

It is further important to note that, despite a net lake area loss due to drainage of large lakes, thermokarst expansion can offset the decrease in GHG emissions and in some regions (i.e. Yedoma Ice Complex) small lake area increase can increase GHG emissions despite a net loss of large lake area (Walter Anthony et al 2016). Moreover, the ultimate response of GHG to the thermokarst lake cycle may be a function of the permafrost type (continuous versus discontinuous), the permafrost carbon stock in terrestrial landscape surrounding lakes (Yedoma, peatland or non-Yedoma mineral soils), the rate of thermokarst expansion, and/or the balance of gross lake area loss versus gain.

Another uncertainty involves the relative role of GHG ebullition versus diffusion which cannot be assessed solely from CO₂ and CH₄ concentrations. Although ebullition is known to dominate over emission in many boreal and Arctic ponds and lakes including those in Siberia (i.e. Walter et al 2006, Anthony et al 2010) there is no relationship between diffusive emissions and ebullition (Wik et al 2016, Sepulveda-Jauregui et al 2015). In lakes of permafrost-free regions of the WSL, ebullition strongly dominated over diffusion in the middle taiga zone but was comparable to diffusion in the southern taiga zone (Sabrekov et al 2017). As such, the response of GHG emission to frozen peat thaw and permafrost boundary change in the WSL may be more complex than predicted exclusively from GHG total concentrations in lentic waters.