Open access

Modelling growing season carbon fluxes at a low-center polygon ecosystem in the Mackenzie River Delta

Publication: Arctic Science
19 June 2023

Abstract

A temporal upscaling study was conducted to estimate net ecosystem exchange (NEE) of carbon dioxide and net methane exchange (NME) for a low-center polygon (LCP) ecosystem in the Mackenzie River Delta, for each of the 11 growing seasons (2009–2019). We used regression models to create a time series of flux drivers from in situ weather observations (2009–2019) combined with ERA5 reanalysis and satellite data. We then used neural networks that were trained and validated on a single growing season (2017) of eddy covariance data to model NEE and NME over each growing season. The study indicates growing season NEE was negative (net uptake) and NME was positive (net emission) in this LCP ecosystem. Cumulative carbon (C) uptake was estimated to be −46.7 g C m−2 (CI95% ± 45.3) per growing season, with methane emissions offsetting an average 5.6% of carbon dioxide uptake (in g C m−2) per growing season. High air temperatures (>15 °C) reduced daily CO2 uptake and cumulative NEE was positively correlated with mean air growing season temperatures. Cumulative NME was positively correlated with the length of the growing season. Our analysis suggests warmer climate conditions may reduce carbon uptake in this LCP ecosystem.

1. Introduction

The Arctic is warming at three to four times the global rate due to Arctic amplification (Flato et al. 2019; Chylek et al. 2022; Rantanen et al. 2022). Warming is causing significant impacts across the Arctic, including permafrost degradation, longer growing seasons, tundra greening, and shifting wetland distribution (Derksen et al. 2019; Frost et al. 2020; Myers-Smith et al. 2020; Turetsky et al. 2020; Kreplin et al. 2021). These impacts can alter the carbon (C) balance of the Arctic and lead to feedbacks that influence the rate of warming (Schuur et al. 2015). Of particular concern are tundra ecosystems in ice-rich permafrost lowlands because they contain high belowground C-stocks, are vulnerable to thermokarst, and are hotspots for CH4 emissions (Olefeldt et al. 2013; Schuur et al. 2015; Kohnert et al. 2018). These ecosystems are also characterized by substantial spatial heterogeneity that is not well resolved by the models used to project the ecological response to climate change (Lara et al. 2020).
Among terrestrial Arctic ecosystems, low-center polygons (LCPs) have considerable spatial heterogeneity in Arctic ecosystems. They develop in ice-rich permafrost lowlands when syngenetic ice wedges deform the land surface, creating an irregular network of elevated polygon rims bounding low polygon centers (Mackay 2000; Liljedahl et al. 2016; French 2017). Carbon dioxide (CO2) and methane (CH4) fluxes vary over small spatial scales (e.g., <20 m) in LCP ecosystems because microtopography influences water table heights and vegetation cover (Olivas et al. 2011; Olefeldt et al. 2013). Polygon centers can emit up to 10 times as much CH4 as adjacent polygon rims (Kutzbach et al. 2004; Sachs et al. 2010). Interannual climate variability has a significant impact on C fluxes from year to year but LCP ecosystems tend to be net growing season CO2 sinks and CH4 sources (Holl et al. 2019; Dengel et al. 2021). Site-specific factors, including vegetation type, active layer thickness, and the relative cover of polygon centers/rims, have a strong influence on the relative magnitude of fluxes between sites (Kutzbach and Wille 2007; van der Molen et al. 2007; Sachs et al. 2008; Wille et al. 2008; Zona et al. 2009, 2010; Runkle et al. 2013).
The Canadian Arctic is underrepresented in the existing network of eddy covariance (EC) sites used to measure ecosystem C fluxes (Delwiche et al. 2021; Pallandt et al. 2021). This means we lack a baseline understanding of C exchange in tundra ecosystems across Canada, and it is difficult to anticipate how ecosystems in this region may respond to climate warming. Our previous study of EC fluxes in an LCP ecosystem in the Mackenzie River Delta was the first such study in Canada (Skeeter et al. 2022). We analyzed half-hourly CO2 and CH4 fluxes measured by EC in an LCP ecosystem on Fish Island in the Mackenzie River Delta over the 2017 growing season (Fig. 1). Mean CO2 and CH4 fluxes were −0.62 (CI95% ± 0.04) g C-CO2 m−2 d−1 and 28.7 (CI95% ± 0.4) mg C-CH4 m−2 d−1, respectively, where negative values indicate net C uptake and positive values indicate emissions from the ecosystem. Unfortunately, we were only able to conduct EC flux measurements at the site for one field season. The mean fluxes observed over the 2017 growing season cannot be considered representative of net ecosystem exchange (NEE) and net methane exchange (NME) exchange over multiseason timescales.
Fig. 1.
Fig. 1. A 30 m resolution ecosystem classification map for Fish Island created following methods outlined in Skeeter (2022). The ecosystem classification was created using a multiyear composite (2015–2022) of LANDSAT8 data obtained from Google Earth Engine and a random forest model trained on areas representative of each ecosystem type. The training area for each ecosystem type was 6.5 ha. The inset in the top right shows a larger scale map of the automated weather station (AWS) and the eddy covariance (EC) station in 2017 along with the boundaries of the 90% cumulative flux footprint for the EC station. Both the main map and inset in the top right are in NAD 1983 UTM Zone 8N. The inset in the top left shows the site's location in the broader Arctic region, using a shapefile of the circumpolar land areas obtained from the Circumpolar Arctic Vegetation Map (Walker et al. 2005). This inset is in the WGS 1984 North Pole Lambert Azimuthal Equal Area coordinate system.
The goal of this study was to build on our previous analysis by investigating how interannual climate variability and warming may influence growing season C exchange in the LCP ecosystem on Fish Island. Using 11 years of data from an automated weather station (AWS) located next to the EC site on Fish Island combined with reanalysis and satellite data, we conducted a temporal upscaling study to model daily growing season NEE and NME from 2009 to 2019. In this paper, we present analysis of growing season C fluxes in the LCP ecosystem, and provide a critical assessment of the upscaling methods.

2. Study site and data

Fish Island (19 km2) is located 10 km south of the Arctic Ocean, in the north-eastern Mackenzie River Delta, Northwest Territories, Canada. Permafrost in this area is continuous and up to 500 m thick (Burn 2017). The island is subject to infrequent episodic flooding during the spring freshet but is protected from the storm surge events that impact exposed coastal locations to the north and west (Morse and Burn 2013; Lantz et al. 2015; Skeeter 2022). LCPs cover approximately half of the land surface of Fish Island (Fig. 1). The LCP ecosystem consists of peat-rich polygon centers separated by low (<20 cm) rims; peat deposits at Fish Island are at least 25 cm thick (Morse et al. 2012). Vegetation is mainly peat moss (Sphagnum spp.), sedge (Carex spp.), and horsetail (Equisetum spp.), along with dwarf willows (Salix spp.) on polygon rims. Air temperatures at Fish Island remain below freezing for much of the year and most biologic activity is restricted to a brief growing season between June and September (Fig. 2). For the purpose of this study, we will define the growing season as the period each year when the ecosystem is snow-free and mean soil temperatures between 5 and 15 cm remain above freezing.
Fig. 2.
Fig. 2. Boxplots of mean daily air temperatures at the Fish Island AWS. The orange lines represent the median, the boxes indicate the interquartile range (IQR) (Q3–Q1), the whiskers indicate Q1 − (1.5*IQR) and Q3 + (1.5*IQR), and the circles represent outliers extending beyond the whiskers.
Our previous study (Skeeter et al. 2022) measured growing season CO2 and CH4 fluxes in the LCP ecosystem at an EC station (69°22′20.20″N, 134°52′51.92″W; WGS 84), 47 m northeast of the Fish Island AWS (Fig. 1). In that study, we developed neural network (NN) models to map the relationships driving fluxes in the LCP ecosystem. One NN modelled NEE, mapping the seasonally variable responses of gross primary productivity (GPP) and ecosystem respiration (ER) to light levels and temperature conditions. Another modelled NME, mapping methanogenesis (CH4 production) and methanotrophy (CH4 consumption), and the transport dynamics that influence the balance between the two. A description of the NN models is given in Section 1 of Appendix A. Figures A2 and A3 show visualizations of the functional relationships driving growing season NEE and NME, respectively.
In this study, we created a flux driver time series, consisting of estimated hourly values of the weather and soil conditions driving NEE and NME that we previously identified. Data from the Fish Island AWS were the basis of our upscaling study. These included hourly air temperature, net shortwave and longwave radiation, rainfall, and wind speed and direction. Since the AWS data were limited in scope (e.g., no subsurface observations), we also incorporated hourly ERA5-Land reanalysis data and daily 500 m resolution MODIS reflectivity data. The ERA5 data were used to estimate weather and soil conditions not directly observed at the AWS. The MODIS data were used to calculate normalized difference snow and vegetation indices (NDSI and NDVI, respectively), to assess daily snow cover and vegetation greenness. Details on the AWS, ERA5, and MODIS data are given in Section 2 of Appendix A.

3. Methods

To create the flux driver time series, regression models were trained to estimate each respective flux driver (11 in total) using a combination of observed (AWS/MODIS) and estimated (ERA5) parameters. The models were trained on data from the 2017 field season, using K-fold (K = 30) cross validation to prevent overfitting. Flux driver estimates were then temporally upscaled from May 1st to October 31st each year from 2009 to 2019. The six-month May–October interval was chosen to be fully inclusive of each year's growing season. Section 2 of Appendix A gives a detailed overview of the data and methods used to estimate each flux driver. Cross validation statistics show that the estimated drivers closely matched observed drivers from 2017 (Table 1).
Table 1.
Table 1. Validation metrics including the coefficient of determination (r2) and root mean squared error (RMSE) for the hourly estimates of each flux driver along with validation metrics for hourly and daily estimates of net ecosystem exchange and net methane exchange.
The drivers that could be calculated from AWS data alone are highlighted in blue in Table 1. With the exception of friction velocity, these were estimated using simple linear regression models to scale the AWS observations to their equivalent EC observations. The AWS data closely matched the EC data, but scaling helped account for the small bias between stations (see Table A1). For friction velocity we used a multivariate model to account for the influence of daytime heating on turbulence. Most other drivers were approximated using a combination of AWS and ERA5 reanalysis data; these are highlighted in orange in Table 1. Soil temperatures and water table depth were estimated using ordinary least squares (OLS) regression models with a large (n = 12) set of inputs (Fig. A5). A similar approach was used to estimate the vapor pressure deficit (VPD); first, dewpoint temperature was estimated with OLS regression and then VPD values were calculated using the Clausius–Clapeyron equation (Fig. A4).
Thaw depth mapped the seasonality of both NEE and NME at Fish Island; this parameter was more difficult to estimate. The aim was to generate a physically meaningful estimate informed by the seasonality of spring thaw and winter freeze-up. Thaw depth was estimated as a function of cumulative thaw degree days (mean air temperature > 0 °C) and snow-free days (NDSI < 0) each year (Fig. A6). The estimates were compared with data from two Circumpolar Active Layer Monitoring (CALM; 2009–2017) network sites to ensure they were plausible (CALM 2022). Estimated maximum thaw depths at Fish Island were generally in good agreement with observed maximum thaw depth at both CALM sites: Lousey Point, an upland tundra site 22 km to the east (Spearman ρ = 0.79; p = 0.01), and Reindeer Station, a forested lowland site 86 km to the south (Spearman ρ = 0.93; p < 0.01). Finally, the growing seasons were defined for each year by identifying the first to last snow-free day (NDSI < 0) each year with mean daily soil temperatures above 0 °C at 5 and 15 cm. This definition is somewhat arbitrary and may not fully encompass each growing season, but it was necessary for this study. Factors governing NEE and NME during and immediately following snowmelt are quite dissimilar from those during most of growing season and our model was not trained on flux observations from times when the soil surface was still frozen.
The upscaled flux drivers were input into the NN models developed for our previous study to calculate hourly NEE and NME; 95% confidence intervals (CI95%) were calculated following Skeeter et al. (2022). The flux estimates and their CI95% were then summed daily. To check that our NEE and NME estimates were plausible, we compared modelled hourly fluxes and daily sums to observations from the EC station in 2017 (Table 1; Fig. S8). Daily NEE values were also compared with satellite-derived 9 km resolution SMAP L4 global daily NEE estimates (Table A4; Fig. S9; Kimball et al. 2021). Overall, the comparisons indicated the daily NEE and NME estimates presented in this study are reasonable during times when conditions were similar to those observed during the 2017 field season. The net ecosystem C balance (net ecosystem carbon balance, NECB) was then calculated as the sum of NEE and NME; where NEE and NME are in g CO2-C and g CH4-C, respectively. This summation ignores the dissolved C but, it is a close approximation of the NECB over shorter timescales (e.g., a growing season) at small spatial scales (Baldocchi 2003; Chapin et al. 2006), so we consider it sufficient for this study.

4. Results

Soil temperatures in polygon centers typically remained above 0 °C from mid-June until mid-September (Fig. 3a). The growing season started between 6 and 18 days after snowmelt (NDSI < 0). The median growing season lasted 102 days. The 2018 season was anomalously short due to a late snowmelt and early freeze-up. The 2012 season was considerably warmer than any other, 1.4 °C above the next warmest growing season (Table 2). Rainfall data are missing for 2009, but the 10 years available show significant variability between seasons (Table 2). The NDVI data show rapid green-up occurred in June, peak greenness was reached by early August, and then gradual senescence occurred until snow cover and freezing conditions returned in September/October (Fig. 3b). Maximum peak season NDVI ranged from 0.62 in 2009 to 0.51 in 2019.
Fig. 3.
Fig. 3. Mean daily conditions from 2009 to 2019 in the low-center polygon ecosystem at Fish Island between day of year (DOY) 155 and 285: (a) polygon center soil temperatures at 5 and 15 cm, (b) daily NDSI and NDVI, (c) net ecosystem exchange (NEE), and (d) net methane exchange (NME). Shaded regions in panels (a) and (b) show the maximum and minimum mean daily values; in (c) and (d) they show daily 95% confidence intervals. This DOY range is fully inclusive of year's each growing season. Plots (c) and (d) include flux estimates for some nongrowing season days for reference purposes only; these days were not included in when calculating average/cumulative growing season flux estimates.
Table 2.
Table 2. Conditions by growing season in the Fish Island LCP ecosystem.
NEE was typically negative (CO2 uptake) from around the summer solstice (day of year (DOY) 172) until the end of August (DOY 243) (Fig. 3c). NME was consistently positive (CH4 emission) over the full study period (Fig. 3d). Mean growing season NEE and NME were −0.49 g C-CO2 m−2 d−1 (CI95% ± 0.44) and 27.3 mg C-CH4 m−2 d−1 (CI95% ± 5.6), respectively. Cumulative growing season NEE and NME are shown by year in Table 2. Outside of the growing season, the NN models have reduced confidence and produce less realistic estimates. These days were excluded from our analysis but were included in Fig. 3 for context.
On a daily basis, growing season NEE had a weak negative correlation with mean air temperatures, while NME had a weak positive correlation with air temperature (Figs. 4a and 4b). High air temperatures had a limiting effect on NEE; daily CO2 uptake was greatest when mean air temperatures were below 15 °C. Daily NEE also had a moderate negative correlation with NDVI; most net CO2 uptake occurred on days with NDVI > 0.4 (Fig. 4c). Daily NME was not correlated with NDVI (Fig. 4d). Accounting for light levels and phenology (NDVI), daily air temperatures had a net positive influence over NEE; warmer temperatures lead to higher ER, offsetting a larger proportion of daily GPP (Table 3).
Fig. 4.
Fig. 4. Scatterplots comparing daily NEE (a and c) and NME (b and d) to mean daily air temperatures (a and b) and daily NDVI (c and d) along with their respective coefficients of determination (r2). The r2 values listed are significant to p < 0.05 for each except plot (d).
Table 3.
Table 3. Coefficients for the OLS regression equation: NEE = x1*PPFD + x2*Ta + x3*NDVI + x4, where NEE is daily net ecosystem exchange in C-CO2 m−2 d−1, PPFD is daily total photon flux density in mol m−2 d−1, Ta is mean air temperature in °C, and NDVI is the daily normalized difference vegetation index.
On a per growing season basis, cumulative NEE was positively correlated with mean air temperature, the number of days with mean air temperatures above 15 °C, and the number of days with NDVI above 0.4 (Table 4). Cumulative NME was most strongly correlated with the length of the growing season each year. Since each growing season contained a different number of days, mean fluxes were not directly comparable to cumulative fluxes. Mean growing season NEE had statistically significant positive (p < 0.05) correlations with each of the four parameters listed in Table 4; mean NME was not significantly correlated with any of those parameters.
Table 4.
Table 4. Spearman ρ correlation coefficients and Theil–Sen slope between cumulative growing season fluxes and selected parameters.
The NECB of the LCP ecosystem was estimated to be −47.0 g C m−2 (CI95% ± 45.0) per growing season; CH4 emissions offset between 4.4% and 11.6% of CO2 uptake per season (g C m−2). Overall, longer/warmer growing seasons resulted in reduced CO2 uptake, increased CH4 emissions, and reduced net C uptake (Fig. 5). Cumulative C uptake reached a seasonal minimum between late August and early September, but net emissions during the senescent period in September offset a significant portion of CO2 uptake in some years. The 2012 season was long and warm; cumulative uptake peaked on August 25th, but the site remained snow-free with soil temperatures above zero until October 6th and net C emissions from late August onward offset 47% of the net C uptake that had occurred that season.
Fig. 5.
Fig. 5. The net ecosystem carbon balance (NECB) by growing season, where the NECB is the sum of net ecosystem exchange and net methane exchange. Growing seasons are colour-coded by mean temperatures (dark blue is coldest, dark red is warmest), which correspond to the values listed in Table 2. The year 2017 is shown in orange with an upside-down triangle; this is the year when EC observations were collected by our previous study.

5. Discussion

In the absence of long-term observations of carbon exchange in the Mackenzie River Delta region, we investigated the influence of climate variability by extrapolating flux observations from a single field season (2017) to 11 growing seasons. Weather data from the Fish Island AWS along with ERA5 reanalysis and MODIS data were used to create a time series of flux drivers. NNs were used to estimate daily NEE and NME and show that there was considerable interannual variability in growing season C exchange in the LCP ecosystem at Fish Island. Cumulative growing season NEE ranged from −26 g C-CO2 m−2 (CI95% ± 54) in 2012 to −58 g C-CO2 m−2 (CI95% ± 47) in 2019, and NME ranged from 2.2 g C-CH4 m−2 (CI95% ± 0.5) in 2018 to 3.1 g C-CH4 m−2 (CI95% ± 0.7) in 2012. On average, the ecosystem was a CO2 sink and CH4 source; the mean net growing season C uptake (NECB) was −47 g C m−2 d−1 (CI95% ± 45). The 95% confidence intervals for NECB spanned zero in 4 of the 11 years, so the ecosystem may not have been a net C sink every growing season.

5.1. Uncertainty

NNs can make projections beyond conditions they were trained on, but error increases as conditions diverge from the training data (Khosravi et al. 2011). Outside of the growing season, our model projected significant CO2 emissions, particularly during the spring snowmelt period (Fig. 3). These cold season fluxes are likely overestimated, and were not considered in the analysis because we lack cold season flux observations from Fish Island. By limiting our analysis to the growing season (snow-free days with soil temperatures above 0 °C), we restricted our estimates to conditions that were not drastically different than the 2017 field season.
Within the growing season, the largest source of uncertainty in the study was our estimate of thaw depth. The relative influence of thaw depth in the NN models was only 7% and 12% for NEE and NME, respectively (Figs. A2 and A3), so the effects of poorly resolved thaw depth were probably small-to-moderate. Given that maximum thaw depths were well correlated with the two CALM sites, the estimates are likely reasonable enough to guide seasonality in the model. There were no CALM data to help validate the thaw depth estimates for 2018 or 2019. However, the comparison between our upscaled NEE estimates and 9 km resolution satellite derived SMAP L4 estimated NEE (Kimball et al. 2021) suggests model performance in 2018 and 2019 was comparable to earlier years. The SMAP NEE data cannot validate our NEE estimates as they are also modelled estimates, and they are averaged across multiple ecosystems; but they do broadly support that our estimates of NEE are plausible.
We are generally confident in our results as validation statistics indicate both the estimated hourly flux drivers and daily fluxes are in close agreement with observations from 2017. The wide confidence intervals we presented represent the potential ranges of growing season NEE, NME, and NECB this LCP ecosystem and adequately convey the inherent uncertainty associated with this sort of upscaling study. The summation NECB = NEE + NME ignores the dissolved C fluxes which can be significant in tundra ecosystems during snowmelt (Kling et al. 1991; Finlay et al. 2006). However, it is a close approximation of the NECB over shorter timescales (e.g., a growing season) at small spatial scales in flat tundra ecosystems (Baldocchi 2003; Chapin et al. 2006), so it is sufficient for this study.
The LCP ecosystem covers approximately 9.5 km2 of the 19 km2 Fish Island (Fig. 1). It is likely the growing season NECB across this LCP ecosystem falls somewhere in the range calculated in the model. The 95% confidence interval for the NECB per growing season was −92.0 to −2.0 g C m−2, which would equate to −874 to −19 t C uptake per growing season across the LCP area at Fish Island. However, we caution that this is a “back of the envelope” calculation; the study was not intended to upscale flux estimates in both space and time. Location bias (Schmid and Lloyd 1999) caused mean CH4 fluxes observed by the EC station in 2017 to be approximately 1 nmol m−1 s−1 higher than NME in the LCP ecosystem (Skeeter et al. 2022). To correct for this, NME was projected to the mean distribution of polygon rims and centers at 2 m resolution, from the 0.7 km2 microtopography classification map centered on the EC tower (Skeeter et al. 2022). We cannot be certain the microtopographic distribution across the broader LCP ecosystem at Fish Island matches that within the vicinity of the EC station. It is encouraging however, that aircraft measurements of CH4 fluxes over Fish Island (17–35 nmol m−2 s−1) collected in July 2012 and 2013 overlap with our estimates (31.9 to 40.0 nmol m−2 s−1) for mean July NME in 2012 and 2013 (Kohnert et al. 2017).

5.2. Net ecosystem carbon balance

Our findings are broadly in alignment with two longer term EC studies from LCP ecosystems in Alaska (Dengel et al. 2021) and Siberia (Holl et al. 2019) and shorter term analysis for two other Siberian sites (van der Molen et al. 2007; Sachs et al. 2008). Dengel et al. (2021) present a 7-year record (2013–2019) of June–September NEE and mean CH4 fluxes (not gap-filled) for a LCP site near Utqiagvik, AK. This site is colder than Fish Island and experienced greater CO2 uptake than Fish Island; mean daily NEE was −1.31 g C-CO2 m−2 d−1 (1 standard error ± 0.42) at Utqiagvik, compared with −0.49 g C-CO2 m−2 d−1 (CI95% ± 0.44) at Fish Island (Dengel et al. 2021). ER at Fish Island was likely higher due to the warmer temperatures which would explain much of the difference in NEE. They did not present gap-filled estimates of CH4 fluxes, so they cannot be directly compared, but our results suggest the LCP ecosystem at Fish Island is a stronger CH4 source than the one at Utqiagvik (Dengel et al. 2019).
Holl et al. (2019) present cumulative CO2-C uptake over the peak growing seasons (DOY 200 to 235) from 2002 to 2017 for an LCP ecosystem at Samoylov Island, in Siberia's Lena River Delta. They did not present error bounds, but yearly values ranged from −25 to −9 g CO2-C m−2 with a median of −15 g CO2-C m−2. Over the same days (2009–2019), we estimate cumulative CO2 uptake in the Fish Island LCP ecosystem was −26.7 (CI95% ± 12.7) g CO2-C m−2. This suggests the LCP ecosystem at Fish Island is a stronger CO2 sink during the peak of the growing season. Sachs et al. (2008) measured CH4 fluxes at Samoylov Island in 2006 (June 9th–September 19th) and estimated NME was ∼14.0 mg C-CH4 m−2 d−1 (no error bounds provided). This is well below the estimated growing season NME at Fish Island (CI95% (21.7, 32.9) mg C-CH4 m−2 d−1). At Samoylov Island, air and soil temperatures and maximum thaw depth in 2006 were within ranges observed/estimated for the LCP ecosystem at Fish Island. The difference is likely because the two LCP ecosystems had different levels of microtopographic relief and distributions of polygon rims/centers. Polygon rims are elevated relative to the water table and provide aerobic conditions that promote methanotrophy which reduces CH4 emissions (Kutzbach et al. 2004; Sachs et al. 2010). At Fish Island, polygon rims were 10–15 cm high and cover 29% of the land surface (polygon centers covered 66% and the other 5% was ponded troughs). At Samoylov, polygon rims were 50 cm high and covered about 60% of the land area and polygon centers covered the other 40% (Sachs et al. 2008; Skeeter et al. 2022).
Some EC and modelling studies suggest earlier snowmelt and longer/warmer growing seasons will increase CO2 uptake in some low Arctic tundra sites (Lafleur and Humphreys 2008; Luus et al. 2013; Tao et al. 2021), while others have found longer/warmer growing seasons and delayed freeze-up in autumn can reduce or completely cancel out growing season CO2 uptake (Zona et al. 2016, 2022; Commane et al. 2017; Holl et al. 2019; Natali et al. 2019; Bao et al. 2021). Our analysis supports the latter; longer/warmer growing seasons result in decreased net C uptake. While photosynthesis can increase and persist for longer in a warmer growing season, the relative increases in soil respiration (ER + NME) with increasing temperatures more than offset any net increase in GPP. Further, the onset of senescence is strongly influenced by decreasing light levels in September; this limits the potential for warming to increase autumn GPP (Livensperger et al. 2019; Myers-Smith et al. 2019).

6. Conclusions

We are confident this LCP ecosystem was a net growing season C sink between 2009 and 2019. However, we cannot claim the site was a C sink on an annual basis because we do not include cool season emissions that could have offset growing season uptake in some years. Given the relationships between reduced net C uptake and warmer temperatures, we anticipate that climate warming will reduce growing season net C uptake in this LCP ecosystem. Further observations are needed to better understand and project how climate change will affect the carbon balance of this region. Continuous annual measurements are not necessarily practical in a remote Arctic location like Fish Island, but any future study should at least make efforts to include the weeks before snowmelt through to the end of early winter freeze-up to provide better bounds on annual fluxes.
Warming at Fish Island may halt ice wedge growth at Fish Island and lead to thermokarst. When polygon rims collapse in this ecosystem, they form ponded troughs that may result in increased CH4 emissions. If Fish Island were revisited, a nested study with multiple towers at different heights could be used to better assess the role of microtopographic features. By including different proportions of polygon rims, centers, and degraded troughs, space-for-time substation could be used to anticipate how future changes to thaw depth and hydrology in the LCP ecosystem may influence the carbon balance.

Acknowledgements

We would like to thank Shawne Kokelj and Casey Beel with the Government of the Northwest Territories Department of Environment and Natural Resources for providing the AWS data for Fish Island. Thanks to Oliver Sonnentag for reviewing a draft of the manuscript and suggesting data sources. Thanks to the staff at the Aurora Research Institute in Inuvik for providing logistical support. We especially thank Edwin Amos for the local expertise and invaluable assistance he provided in the field, and Rick Ketler for helping to design and install the EC system used in the 2017 study.

References

Aubinet M., Vesala T., Papale D.(Eds.). 2012. Eddy Covariance: A Practical Guide to Measurement and Data Analysis (2012 edition). Dordrecht, NY, Springer.
Baldocchi D.D. 2003. Assessing the eddy covariance technique for evaluating carbon dioxide exchange rates of ecosystems: past, present and future. Global Change Biology, 9(4): 479–492.
Bao T., Xu X., Jia G., Billesbach D.P., Sullivan R.C. 2021. Much stronger tundra methane emissions during autumn freeze than spring thaw. Global Change Biology, 27(2): 376–387.
Burn C.R. 2017. Mackenzie Delta: Canada's principal Arctic Delta. In Landscapes and landforms of Western Canada. Edited by Slaymaker O. Springer International Publishing, Cham. pp. 321–334.
CALM. 2022. Circumpolar Active Layer Monitoring Network - International Permafrost Association.Available from https://ipa.arcticportal.org/products/gtn-p/calm  [accessed 12 March 2022].
Chapin F.S., Woodwell G.M., Randerson J.T., Rastetter E.B., Lovett G.M., Baldocchi D.D., et al. 2006. Reconciling carbon-cycle concepts, terminology, and methods. Ecosystems, 9(7): 1041–1050.
Chylek P., Folland C., Klett J.D., Wang M., Hengartner N., Lesins G., Dubey M.K. 2022. Annual mean Arctic amplification 1970–2020: observed and simulated by CMIP6 climate models. Geophysical Research Letters, 49(13): e2022GL099371.
Commane R., Lindaas J., Benmergui J., Luus K.A., Chang R.Y.-W., Daube B.C., et al. 2017. Carbon dioxide sources from Alaska driven by increasing early winter respiration from Arctic tundra. Proceedings of the National Academy of Sciences, 114(21): 5361–5366.
Delwiche K.B., Knox S.H., Malhotra A., Fluet-Chouinard E., McNicol G., Feron S., et al. 2021. FLUXNET-CH4: a global, multi-ecosystem dataset and analysis of methane seasonality from freshwater wetlands. Earth System Science Data, 13(7): 3607–3689.
Dengel S., Zona D., Sachs T., Aurela M., Jammet M., Parmentier F.J.W., et al. 2013. Testing the applicability of neural networks as a gap-filling method using CH4 flux data from high latitude wetlands. Biogeosciences, 10(12): 8185–8200.
Dengel S., Billesbach D., Torn M.S. 2021. Influence of tundra polygon type and climate variability on CO2 and CH4 fluxes near Utqiagvik, Alaska. Journal of Geophysical Research: Biogeosciences, 126(12): e2021JG006262.
Derksen C., Burgess D., Duguay C., Howell S., Mudryk L., Smith S., et al. 2019. Changes in snow, ice, and permafrost across Canada. In Canada's changing climate report. Edited by Bush E., Lemmen D.S. Government of Canada, Ottawa, Ont., Canada. pp. 194–260.
Finlay J., Neff J., Zimov S., Davydova A., Davydov S. 2006. Snowmelt dominance of dissolved organic carbon in high-latitude watersheds: implications for characterization and flux of river DOC. Geophysical Research Letters, 33(10).
Flato G., Gillett N., Arora V., Cannon A., Anstey J. 2019. Modelling Future Climate Change; Chapter 3 in Canada's Changing Climate Report. pp. 74–111.
French H.M. 2017. Ground freezing, permafrost and the active layer. In The periglacial environment. 4th ed. John Wiley & Sons, Ltd. pp. 63–85.
Frost G.V., Bhatt U.S., Epstein H.E., Myers-Smith I., Phoenix G.K., Berner L.T., et al. 2020. Arctic Report Card 2020: Tundra Greenness.
Gevrey M., Dimopoulos I., Lek S. 2003. Review and comparison of methods to study the contribution of variables in artificial neural network models. Ecological Modelling, 160(3): 249–264.
Heskes T. 1997. Practical confidence and prediction intervals. In Advances in neural information processing systems. MIT Press, Cambridge. pp. 176–182.
Holl D., Wille C., Sachs T., Schreiber P., Runkle B.R.K., Beckebanze L., et al. 2019. A long-term (2002 to 2017) record of closed-path and open-path eddy covariance CO2 net ecosystem exchange fluxes from the Siberian Arctic. Earth System Science Data, 11(1): 221–240.
Hornik K. 1991. Approximation capabilities of multilayer feedforward networks. Neural Networks, 4(2): 251–257.
Irvin J., Zhou S., McNicol G., Lu F., Liu V., Fluet-Chouinard E., et al. 2021. Gap-filling eddy covariance methane fluxes: comparison of machine learning model predictions and uncertainties at FLUXNET-CH4 wetlands. Agricultural and Forest Meteorology, 308–309: 108528.
Joabsson A., Christensen T.R., Wallén B. 1999. Vascular plant controls on methane emissions from northern peatforming wetlands. Trends in Ecology & Evolution, 14(10): 385–388.
Khosravi A., Nahavandi S., Creighton D., Atiya A.F. 2011. Comprehensive review of neural network-based prediction intervals and new advances. IEEE Transactions on Neural Networks, 22(9): 1341–1356.
Kimball J., Jones L.A., Arthur E., Tobias K., Rolf R. 2021. SMAP L4 Global Daily 9 km EASE-Grid Carbon Net Ecosystem Exchange, Version 6 [data set]. NASA National Snow and Ice Data Center DAAC.
Kling G.W., Kipphut G.W., Miller M.C. 1991. Arctic lakes and streams as gas conduits to the atmosphere: implications for tundra carbon budgets. Science, 251(4991): 298–301.
Knox S.H., Jackson R.B., Poulter B., McNicol G., Fluet-Chouinard E., Zhang Z., et al. 2019. FLUXNET-CH4 synthesis activity: objectives, observations, and future directions. Bulletin of the American Meteorological Society, 100(12): 2607–2632.
Kohnert K., Serafimovich A., Metzger S., Hartmann J., Sachs T. 2017. Strong geologic methane emissions from discontinuous terrestrial permafrost in the Mackenzie Delta, Canada. Scientific Reports, 7(1): 1–6.
Kohnert K., Juhls B., Muster S., Antonova S., Serafimovich A., Metzger S., et al. 2018. Toward understanding the contribution of waterbodies to the methane emissions of a permafrost landscape on a regional scale: a case study from the Mackenzie Delta, Canada. Global Change Biology, 24(9): 3976–3989.
Kreplin H.N., Santos Ferreira C.S., Destouni G., Keesstra S.D., Salvati L., Kalantari Z. 2021. Arctic wetland system dynamics under climate warming. WIREs Water, 8(4): e1526.
Kutzbach L., Wille C. 2007. The exchange of carbon dioxide between wet arctic tundra and the atmosphere at the Lena River Delta, Northern Siberia. Biogeosciences, 4: 869–890.
Kutzbach L., Wagner D., Pfeiffer E.-M. 2004. Effect of microrelief and vegetation on methane emission from wet polygonal tundra, Lena Delta, Northern Siberia. Biogeochemistry, 69(3): 341–362.
Lafleur P.M., Humphreys E.R. 2008. Spring warming and carbon dioxide exchange over low Arctic tundra in central Canada.Global Change Biology, 14(4): 740–756.
Lai D.Y.F. 2009. Methane dynamics in northern peatlands: a review. Pedosphere, 19(4): 409–421.
Lantz T., Kokelj S., Fraser R. 2015. Ecological recovery in an Arctic delta following widespread saline incursion. Ecological Applications, 25: 140915094202006.
Lara M.J., McGuire A.D., Euskirchen E.S., Genet H., Yi S., Rutter R., et al. 2020. Local-scale Arctic tundra heterogeneity affects regional-scale carbon dynamics. Nature Communications, 11(1): 4925.
Lee S.-C., Christen A., Black T.A., Jassal R.S., Briegel F., Nesic Z. 2021. Combining flux variance similarity partitioning with artificial neural networks to gap-fill measurements of net ecosystem production of a Pacific Northwest Douglas-fir stand. Agricultural and Forest Meteorology, 303: 108382.
Liljedahl A.K., Boike J., Daanen R.P., Fedorov A.N., Frost G.V., Grosse G., et al. 2016. Pan-arctic ice-wedge degradation in warming permafrost and its influence on tundra hydrology. Nature Geoscience, 9(4): 312–318.
Livensperger C., Steltzer H., Darrouzet-Nardi A., Sullivan P.F., Wallenstein M., Weintraub M.N. 2019. Experimentally warmer and drier conditions in an Arctic plant community reveal microclimatic controls on senescence. Ecosphere, 10(4): e02677.
Lucht W., Schaaf C.B., Strahler A.H. 2000. An algorithm for the retrieval of albedo from space using semiempirical BRDF models. IEEE Transactions on Geoscience and Remote Sensing, 38(2): 977–998.
Luus K.A., Kelly R.E.J., Lin J.C., Humphreys E.R., Lafleur P.M., Oechel W.C. 2013. Modeling the influence of snow cover on low Arctic net ecosystem exchange. Environmental Research Letters, 8(3): 035045.
Mackay J.R. 2000. Thermally induced movements in ice-wedge polygons, western Arctic coast: a long-term study. Géographie Physique et Quaternaire, 54(1): 41–68.
Moffat A.M., Papale D., Reichstein M., Hollinger D.Y., Richardson A.D., Barr A.G., et al. 2007. Comprehensive comparison of gap-filling techniques for eddy covariance net carbon fluxes. Agricultural and Forest Meteorology, 147(3–4): 209–232.
Morse P.D., Burn C.R. 2013. Field observations of syngenetic ice wedge polygons, outer Mackenzie Delta, western Arctic coast, Canada. Journal of Geophysical Research: Earth Surface, 118(3): 1320–1332.
Morse P.D., Burn C.R., Kokelj S.V. 2012. Influence of snow on near-surface ground temperatures in upland and alluvial environments of the outer Mackenzie Delta, Northwest Territories. Canadian Journal of Earth Sciences, 49(8): 895–913.
Muñoz S.J. 2019. ERA5-Land hourly data from 1981 to present [data set]. Copernicus Climate Change Service (C3S) Climate Data Store (CDS).
Myers-Smith I.H., Grabowski M.M., Thomas H.J.D., Angers-Blondin S., Daskalova G.N., Bjorkman A.D., et al. 2019. Eighteen years of ecological monitoring reveals multiple lines of evidence for tundra vegetation change. Ecological Monographs, 89(2): e01351.
Myers-Smith I.H., Kerby J.T., Phoenix G.K., Bjerke J.W., Epstein H.E., Assmann J.J., et al. 2020. Complexity revealed in the greening of the Arctic. Nature Climate Change, 10(2): 106–117.
Natali S.M., Watts J.D., Rogers B.M., Potter S., Ludwig S.M., Selbmann A.-K., et al. 2019. Large loss of CO2 in winter observed across the northern permafrost region. Nature Climate Change, 9(11): 852–857.
Olefeldt D., Turetsky M.R., Crill P.M., McGuire A.D. 2013. Environmental and physical controls on northern terrestrial methane emissions across permafrost zones. Global Change Biology, 19(2): 589–603.
Olivas P.C., Oberbauer S.F., Tweedie C., Oechel W.C., Lin D., Kuchy A. 2011. Effects of fine-scale topography on CO2 flux components of Alaskan coastal plain tundra: response to contrasting growing seasons. Arctic, Antarctic, and Alpine Research, 43(2): 256–266.
Pallandt M., Kumar J., Mauritz M., Schuur E., Virkkala A.-M., Celis G., et al. 2021. Representativeness assessment of the pan-Arctic eddy-covariance site network, and optimized future enhancements. Biogeosciences Discussions, 1–42.
Papale D., Valentini R. 2003. A new assessment of European forests carbon exchanges by eddy fluxes and artificial neural network spatialization. Global Change Biology, 9(4): 525–535.
Ramirez J.A., Baird A.J., Coulthard T.J., Waddington J.M. 2015. Ebullition of methane from peatlands: does peat act as a signal shredder? Geophysical Research Letters, 42(9): 3371–3379.
Rantanen M., Karpechko A.Y., Lipponen A., Nordling K., Hyvärinen O., Ruosteenoja K., et al. 2022. The Arctic has warmed nearly four times faster than the globe since 1979. Communications Earth & Environment, 3(1): 1–10.
Runkle B.R.K., Sachs T., Wille C., Pfeiffer E.-M., Kutzbach L. 2013. Bulk partitioning the growing season net ecosystem exchange of CO2 in Siberian tundra reveals the seasonality of its carbon sequestration strength. Biogeosciences, 10(3): 1337–1349.
Sachs T., Wille C., Boike J., Kutzbach L. 2008. Environmental controls on ecosystem-scale CH4 emission from polygonal tundra in the Lena River Delta, Siberia. Journal of Geophysical Research: Biogeosciences, 113: G00A03.
Sachs T., Giebels M., Boike J., Kutzbach L. 2010. Environmental controls on CH4 emission from polygonal tundra on the microsite scale in the Lena river delta, Siberia. Global Change Biology, 16(11): 3096–3110.
Sarle W.S. 1995. Stopped training and other remedies for overfitting. In Proceedings of the 27th Symposium on the Interface of Computing Science and Statistics. pp. 352–360.
Schaaf C. 2021. MODIS User Guide V006 and V006.1. Available fromhttps://www.umb.edu/spectralmass/terra_aqua_modis/v006/introduction [accessed 12 October 2021].
Schaaf C.B., Gao F., Strahler A., Lucht W., Li X., Tsang T., et al. 2002. First operational BRDF, albedo nadir reflectance products from MODIS. Remote Sensing of Environment, 83(1–2): 135–148.
Schäfer K.V.R., Duman T., Tomasicchio K., Tripathee R., Sturtevant C. 2019. Carbon dioxide fluxes of temperate urban wetlands with different restoration history. Agricultural and Forest Meteorology, 275: 223–232.
Schmid H.P., Lloyd C.R. 1999. Spatial representativeness and the location bias of flux footprints over inhomogeneous areas. Agricultural and Forest Meteorology, 93: 195–209.
Schuur E.A.G., McGuire A.D., Schädel C., Grosse G., Harden J.W., Hayes D.J., et al. 2015. Climate change and the permafrost carbon feedback. Nature, 520(7546): 171–179.
Skeeter J. 2022. Using machine learning to identify and map controls of growing-season carbon dioxide and methane fluxes in the Mackenzie Delta region. Ph.D. thesis, University of British Columbia.
Skeeter J., Christen A., Henry G.H.R. 2022. Controls on carbon dioxide and methane fluxes from a low-center polygonal peatland in the Mackenzie River Delta, Northwest Territories. Arctic Science, 8(2): 471–497.
Tao J., Zhu Q., Riley W.J., Neumann R.B. 2021. Warm-season net CO2 uptake outweighs cold-season emissions over Alaskan North Slope tundra under current and RCP8.5 climate. Environmental Research Letters, 16(5): 055012.
Turetsky M.R., Abbott B.W., Jones M.C., Anthony K.W., Olefeldt D., Schuur E.A.G., et al. 2020. Carbon release through abrupt permafrost thaw. Nature Geoscience, 13(2): 138–143.
van der Molen M.K., Huissteden J., Parmentier F.J.W., Petrescu A.M.R., Dolman A.J., Maximov T.C., et al. 2007. The growing season greenhouse gas balance of a continental tundra site in the Indigirka lowlands, NE Siberia. Biogeosciences, 4(6): 985–1003.
Walker D.A., Raynolds M.K., Daniëls F.J.A., Einarsson E., Elvebakk A., Gould W.A., et al. 2005. The circumpolar Arctic vegetation map. Journal of Vegetation Science, 16(3): 267–282.
Wanner W., Strahler A.H., Hu B., Lewis P., Muller J.-P., Li X., et al. 1997. Global retrieval of bidirectional reflectance and albedo over land from EOS MODIS and MISR data: theory and algorithm. Journal of Geophysical Research: Atmospheres, 102(D14): 17143–17161.
Wille C., Kutzbach L., Sachs T., Wagner D., Pfeiffer E.-M. 2008. Methane emission from Siberian arctic polygonal tundra: eddy covariance measurements and modeling. Global Change Biology, 14(6): 1395–1408.
Zona D., Oechel W.C., Kochendorfer J., Paw U.K.T., Salyuk A.N., Olivas P.C., et al. 2009. Methane fluxes during the initiation of a large-scale water table manipulation experiment in the Alaskan Arctic tundra. Global Biogeochemical Cycles, 23: GB2013.
Zona D., Oechel W.C., Peterson K.M., Clements R.J., Paw U.K.T., Ustin S.L. 2010. Characterization of the carbon fluxes of a vegetated drained lake basin chronosequence on the Alaskan Arctic Coastal Plain. Global Change Biology, 16(6): 1870–1882.
Zona D., Gioli B., Commane R., Lindaas J., Wofsy S.C., Miller C.E., et al. 2016. Cold season emissions dominate the Arctic tundra methane budget. Proceedings of the National Academy of Sciences of the United States of America, 113(1): 40–45.
Zona D., Lafleur P.M., Hufkens K., Bailey B., Gioli B., Burba G., et al. 2022. Earlier snowmelt may lead to late season declines in plant productivity and carbon sequestration in Arctic tundra ecosystems. Scientific Reports, 12(1): 3986.

Appendix A

A.1. Neural network models

This section gives a brief summary of the neural network (NN) models developed in our previous study to map the functional relationships driving net ecosystem exchange (NEE) (CO2 exchange) and net methane exchange (NME) (CH4 exchange), respectively. For an overview of the data collected during the 2017 study and a comprehensive discussion of the NN models, we direct the reader to Skeeter et al. (2022). NNs are universal approximators; they make no prior assumptions about the distribution of a dataset and with enough hidden nodes they are capable of mapping any continuous function to an arbitrary degree of accuracy (Hornik 1991; Khosravi et al. 2011). Given this, an NN will fit any pattern, real or artificial, so care must be taken to ensure model estimates are physically plausible. Our previous study took a number of steps to prevent overfitting and ensure our models were plausible.
Both models (one for NEE and one for NME) were ensembles of 30 randomly initialized feed forward single layer NN, with each NN trained on a bootstrapped sample of the 2017 eddy covariance (EC) data. Early stopping was used to prevent overfitting (Sarle 1995; Heskes 1997). Ensemble means and standard deviations were used to calculate hourly NEE and NME estimates with 95% confidence intervals (Khosravi et al. 2011). First, overparametrized models were trained on a large set (n = 21) of potential inputs. Then, weights method (Gevrey et al. 2003) was used to prune unnecessary inputs and the models were retrained using only the most important inputs (drivers). This approach had not been applied to EC data before, distinguishing our work from other EC studies that have treated NN as “black box” models (e.g., Papale and Valentini 2003; Moffat et al. 2007; Dengel et al. 2013; Knox et al. 2019; Schäfer et al. 2019; Lee et al. 2021).
Figure A1a shows the fit between observed CO2 fluxes and modelled NEE, Fig. A2 shows observed modelled relationships driving NEE. The NN model was not set up to explicitly partition NEE into its component fluxes, ecosystem respiration (ER) and gross primary productivity (GPP), where NEE = ER − GPP (Aubinet et al. 2012). However, modelled NEE provided a strong fit with the EC observations (r2 = 0.95, RMSE = 0.33 µmol m−2 s−1) and the NN correctly identified the relevant drivers of both ER and GPP. Photosynthetic photon flux density (PPFD) was the primary driver of GPP, which was modulated by vapor pressure deficit (VPD). The NN mapped the exponential relationship between soil temperature and ER; the response function was similar to the commonly used Q10 approach but was weighted to account for the differential respiration rates of polygon rims and centers. Thaw depth mapped seasonality in the NN, which was a conflation of seasonal phenology changes and the increasing volume of substrate available to respiration. This model also included a one-hot-coded “daytime” variable denoting whether the sun was above/below the horizon to explicitly distinguish between day and night.
Fig. A1.
Fig. A1. Comparison between observed half-hourly fluxes of CO2 and CH4 and modelled half-hourly NEE (a) and NME (b), respectively.
Fig. A2.
Fig. A2. The modelled partial first derivatives of the drivers of net ecosystem exchange (NEE). Each subplot shows the partial first derivative of NEE with respect to the flux driver listed in the title. The derivates are unitless because the flux data were z-score normalized before NN training. They represent relative changes per unit change in the given flux driver. The ensemble mean is shown in red and the 95% confidence interval in blue. Plots are ordered by the relative influence (RI) of each driver, which is listed in the title; the RI is calculated from the sum of squared partial derivatives (SSD) of each driver normalized by the sum of SSD values across all drivers. The grey bars show histograms of the distribution of training samples for each driver, representing the range of conditions the model was trained on.
Fig. A1b shows the fit between observed CH4 fluxes and modelled NME; Fig. A3 shows the modelled relationships driving NME. Methane fluxes are complex; the literature does not suggest any broadly applicable, well-defined functional relationships governing ecosystem scale NME (Irvin et al. 2021). The balance of methanogenesis and methanotrophy dictates NME and both components are driven by soil temperature, moisture, and substrate availability (Lai 2009). However, substrate availability can change over short timescales and the rate of CH4 transport varies by emission pathway (diffusion, ebullition, or vascular transport), limiting methanotrophy (Joabsson et al. 1999; Ramirez et al. 2015). Modelled NME provided a moderate fit with the EC observations (r2 = 0.79, RMSE = 5.89 nmol m−2 s−1). The NN model showed that at subdaily timescales, NME was driven by the dependence of methanogenesis on GPP, the dependence of vascular transport on net radiation (Rn), ventilation of peat under turbulent conditions, and the relative contributions of polygon centers/rims to flux observations. At longer timescales, soil temperature, water table depth, and thaw depth were the important drivers of NME.
Fig. A3.
Fig. A3. The modelled partial first derivatives of the drivers of net methane exchange (NME). Each subplot shows the partial first derivative of NME with respect to the flux driver listed in the title. The derivates are unitless because the flux data were z-score normalized before NN training. They represent relative changes per unit change in the given flux driver. The ensemble mean is shown in red and the 95% confidence interval in blue. Plots are ordered by the relative influence (RI) of each driver, which is listed in the title; the RI is calculated from the sum of squared partial derivatives (SSD) of each driver normalized by the sum of SSD values across all drivers. The grey bars show histograms of the distribution of training samples for each driver, representing the range of conditions the model was trained on.
The relative areal fractions of polygon centers (66%) and rims (29%) at the ecosystem scale were slightly different than the median distribution of polygon centers (63%) and rims (23%) in the flux footprint climatology; the flux footprint climatology does not sum to 100% because the footprint function is technically infinite. We estimated that location bias (placement of the EC station) caused our NME estimates to be approximately 1.0 nmol m−2 s−1 above background levels in the low-center polygon (LCP) ecosystem at Fish Island. We did not find location bias impacted our NEE estimates because flux observations were not as sensitive to changes in the flux source area. We did find CO2 fluxes varied between polygon centers and rims at the plot scale, but this variability averaged out at the footprint scale because the magnitude of the differences was less drastic than for CH4.

A.2. Driver estimation

A.2.1. Fish Island AWS

In August 2008 an AWS was installed in the LCP ecosystem at Fish Island (69°22′19.7″N, 134°52′55.9″W), operated by the Canadian Department of Indian Affairs and Northern Development, Water Resources Division. The AWS recorded hourly data on a CR1000 data logger equipped with a CFM100-XT data module (Campbell Scientific Inc., Logan, UT, USA; CSI). The AWS measured wind speed (U) and direction at 2.9 m with a 05103AP-10-L Wind Monitor (RM Young, Traverse City, MI, USA), air temperature (Ta) at 2 m with a 44212 Temperature Thermilinear Probe (CSI), incoming and outgoing longwave (LW) and shortwave (SW) radiation at 1.55 m with a CNR2 4-component net radiometer (Kipp & Zonen, Delft, NL). On August 21st 2009, an unheated TE525M Tipping Bucket (Texas Electronics, Dallas, TX, USA) was installed to measure rainfall (rain). Data were nearly complete over the full period of record, except rainfall, which are missing until August 21st 2009 and during periods of snowfall. Data gaps of 2 h or less were filled with linear interpolation and longer gaps were filled using the average of corresponding hourly values for the day before and after the gap. We obtained the AWS data (spanning August 28th 2008–August 11th 2020) through personal communications with the Government of the Northwest Territories Department of Environment and Natural Resources.
Over the 2017 study period, AWS data were in close agreement with hourly averages of equivalent observations from our EC station; our EC station was located 47 m northeast of the AWS in 2017. Regression models were used to scale each AWS variable to its equivalent EC variable to remove bias between the two observations (Table A1). Since friction velocity was not measured at the AWS, we used a slightly different approach to estimate this variable. An exponential regression model was fit to the half-hourly EC data:
where the coefficients were a = 0.082 m s−1, b = 0.96 m s−1, and c = 6.0 × 10−5 m3 J−1. The model provided a reasonable fit (r2 = 0.92; RMSE = 0.039 m s−1) and accounted for the influence of diurnal heating on friction velocity. Then, hourly u* was estimated from the scaled estimates of U and NR measured by the AWS using these coefficients.
Table A1.
Table A1. The slope and intercept of the simple linear regression models used to scale the data from the automated weather station (AWS) at Fish Island to equivalent values from the eddy covariance (EC) station that was 47 m NE of the AWS in 2017.

A.2.2. ERA5 reanalysis

One-hour resolution ERA5-Land reanalysis data produced by the European Center for Medium-Range Weather Forecasting were used to complement the AWS data (Muñoz 2019). Hourly ERA5 data from January 1st 2008 to December 31st 2020 were obtained at 11 km grid resolution from Google Earth Engine. Fish Island was covered by parts of three grid cells, so the spatially weighted mean of each hourly value was calculated over Fish Island. The data included estimates of air temperature and dewpoint at 2 m, total precipitation (liquid equivalent), fractional snow cover, soil temperature and volumetric water content. The soil data were averaged across depth profiles: (1) 0–7 cm, (2) 7–28 cm, and (3) 28–100 cm.
The reanalysis estimates of air temperature were in strong agreement with observations from the EC station in 2017 (r2 = 0.90; RMSE = 1.5 K) and with observations from the AWS over its full period of record (r2 = 0.96; RMSE = 3.0 K). Reanalysis estimated rainfall data were included in this study because a rain gauge was not installed at the AWS until late August 2019 and because the AWS rainfall data did not include snow. On an hourly basis, the observed and estimated rainfall were poorly correlated (r2 = 0.19), but when summed over 48 h intervals the correlation improved (r2 = 0.63). Monthly ERA5 precipitation totals were about 74% higher than AWS observations, but the values were reasonably correlated (r2 = 0.82).
Reanalysis estimated dewpoint temperature was in moderate agreement with the EC station in 2017 (r2 = 0.71; RMSE = 1.7 K). These data were necessary to model the VPD at Fish Island, because humidity was not measured at the AWS. The reanalysis appears to underestimate the maritime influences at Fish Island (Fig. A4a). For example, when winds were off the ocean from the northwest (315°–360°), estimated dewpoint temperature was less accurate (r= 0.33; RMSE = 2.2 K) compared with when winds were from the southeast (135°–180°) with no maritime influence (r= 0.86; RMSE = 1.2 K). To estimate VPD, we estimated the dewpoint from ERA5 data, then calculated the VPD using the Clausius–Clapeyron equation with the estimated dewpoint and observed air temperatures. Because the dewpoint estimate was only accurate for certain wind sectors, we used an ordinary least squares (OLS) regression model with seven inputs to scale the dewpoint estimates and account for the continental bias of the reanalysis data. The seven OLS inputs included ERA5 estimate dewpoint and air temperature, along with five AWS observations of Ta, SW, LW, and the u and v wind components. Inputs were z-score normalized before training and relative magnitude of each coefficient is shown in Fig. A4b. The OLS scaled dewpoint temperature provided a better fit with observations from the EC station in 2017 (r2 = 0.86; RMSE = 1.0 K) (Fig. A4c). Before calculating the VPD, any dewpoint estimates that were above the Ta at the AWS were fixed to Ta. The final VPD estimates provided a reasonable fit with observations from 2017 (r2 = 0.88; RMSE = 120 hPa) (Fig. A4d).
Fig. A4.
Fig. A4. (a) Observed dewpoint temperatures at the EC station versus reanalysis estimates from 2017. (b) Z-score normalized OLS coefficients used to estimate dewpoint temperature at Fish Island. (c) Observed dewpoint temperatures at the EC station versus OLS estimated dewpoint temperature. (d) Observed vapor pressure deficits (VPD) at the EC station versus estimated VPD.
The ERA5 data were also needed because the AWS data did not include any subsurface observations. The ERA5 soil data are integrated over depth levels; no single estimate is directly comparable to any observation. Table A2 shows the correlations between soil conditions observed at EC station in 2017 to reanalysis estimates. The ERA5 soil temperatures for levels one and two were moderately correlated with soil temperatures at 5 cm from polygon rims and polygon center, respectively, but the reanalysis tended to over/under estimate the daily maximum/minimums. Polygon rims are elevated and drier than polygon centers, which explains why they had better associations with reanalysis estimates at levels 1 and 2, respectively. The relationships between reanalysis estimates with temperatures at 15 cm were weak for all levels, but correlations did improve when averaged daily (not shown). Reanalysis estimates of soil moisture at both levels were well correlated with water table depth observations. Each soil parameter was estimated using an OLS regression model, following a similar procedure as was used to estimate the dewpoint temperature. Each model had nine inputs that were z-score normalized before training (Fig. A5). The inputs included the five ERA5 subsurface parameters: ERA5 estimated rainfall (rolling 48 h totals); AWS observations of Ta, SW, LW; and rolling 6-h average air temperatures. Each estimated parameter provided a reasonable fit (r2 ≥ 0.88) with its corresponding observed values in 2017.
Fig. A5.
Fig. A5. Z-score normalized OLS coefficients used to estimate soil temperatures and water table deficits at Fish Island along with their K-fold cross-validation statistics.
Table A2.
Table A2. A correlation matrix (Pearson's r) showing the association between mean hourly values from the EC station over the 2017 field season and ERA5 estimates over the same time period.

A.2.3. MODIS data

Daily 500 m resolution Nadir Bidirectional Reflectance Distribution Function Adjusted Reflectance (NBAR) MODIS data were obtained using Google Earth Engine. The NBAR data product is a daily composite of surface reflectance covering every pixel of the global land surface; it is derived from an algorithm that models the surface reflectivity using multidate, multiangular, cloudfree, atmospherically corrected, surface reflectance observations (Wanner et al. 1997; Lucht et al. 2000; Schaaf 2021). Daily NBAR data are modelled using observations from MODIS instruments on the Terra and Aqua satellites over a 16-day period, with the image date centered on the 9th day, weighted data as a function of quality, coverage, and distance (temporal) from the day of interest (Schaaf et al. 2002; Schaaf 2021). The NBAR data aims to minimize missing data due to cloud cover while providing representative estimates of reflectivity (Schaaf 2021). The NBAR data had good coverage, but there were missing data in September (4%) and October (25%) due to cloudy conditions and low light levels; these data gaps were backfilled linearly. The NBAR data were spatially averaged across the LCP ecosystem at Fish Island (defined by the landscape classification map in Fig. 1) and used to calculate daily normalized difference vegetation index (NDVI) (Bands 1 & 2) and Normalized Difference Snow Index (NDSI) (Bands 4 & 6) for the periods between May 1st and October 31st, each year from 2009 to 2019.
When NDSI > 0, it indicates a pixel is snow covered, so days were classified as snowy if NDSI > 0. Comparing the NDSI data (excluding gaps) with ERA5 fractional snow cover in May/June, the values were moderately correlated (r2 = 0.61) but the reanalysis estimated earlier snowmelt dates than were observed by satellite; in September/October, the two metrics were in better agreement (r2 = 071). The NDSI data were combined with AWS mean daily air temperatures (Ta) and used to calculate a daily thaw degree day index:
were and . Then, thaw depth (TD) was then calculated by fitting a limited number (n = 50) thaw depth measurements from 2017 to an exponential regression model:
Thaw depth was measured manually on 5 days (day of year 174, 184, 213, 233, and 258; 10 replications each day) in 2017. The model was trained with K-fold (K = 10) cross validation; the coefficients a = 0.035 and b = 0.53; the fit with the 2017 thaw depth observations was r2 = 0.92 and RMSE = 0.02 m (Fig. A6a). With the exponential function, thaw depth increases rapidly with snowmelt and continued to increase monotonically over the season until freezing conditions returned in late fall (Fig. A6b). Thaw depth remained fixed at its maximum value until the end of the year; in reality, the active layer refreezes from the top and bottom and a portion of the active layer remains thawed well into early winter (Natali et al. 2019). However, there is no reasonable way to incorporate these dynamics here, and estimating C fluxes during freeze-up is beyond the scope of this analysis.
Fig. A6.
Fig. A6. (a) Thaw depth estimates from fit to observations from 2017, (b) thaw depth each year, (c) maximum annual thaw depth at Fish Island (y axis) plotted against Lousey point (x axis), and (d) maximum annual thaw depth at Fish Island (y axis) plotted against Reindeer Station (x axis).
We obtained maximum annual thaw depth data from Lousey Point and Reindeer Depot, two Circumpolar Active Layer Monitoring network sites (CALM, 2022). Lousey Point is an upland site on Richard's Island, 22 km east of Fish Island and Reindeer Depot is a forested site in the delta 86 km south of Fish Island. The data included values through 2017; both CALM stations had thicker active layers than Fish Island, but these data were needed to provide an independent assessment of our estimates at Fish Island (Figs. A7c and A7d). Both stations Lousey Point (Spearman ρ = 0.79; p = 0.01) and Reindeer Depot (Spearman ρ = 0.93; p < 0.01) provided a reasonable fit with estimates of maximum annual thaw at Fish Island. To check for the plausibility of early season thaw depth estimates, May and June soil temperature estimates were compared with mean daily soil temperatures at 5 and 15 cm with linear regression models. They had reasonable fits (r2 ≥ 0.57) and the 0 °C intercepts with estimated thaw depths were reasonably close to the respective depth levels of the temperature observations (Table A3). Given the simplicity of the thaw depth estimation procedures, these relationships are encouraging.
Fig. A7.
Fig. A7. Comparison between observed hourly fluxes of CO2 and CH4 and modelled NEE (a) and NME (b), respectively, and daily gap-filled hourly fluxes of CO2 and CH4 and modelled NEE (c) and NME (d). Hourly observations were the hourly means of half-hourly flux data originally presented in Skeeter et al. (2022). The gap-filled daily sums of gap-filled half-hourly flux data originally presented in Skeeter et al. (2022).
Table A3.
Table A3. Coefficients of simple regression lines Y = mX + b, where Y was estimated daily thaw depth in May/June and the X parameters were the corresponding estimates of mean daily soil temperatures in May/June.

A.3. Evaluating the upscaled fluxes

Figure A7shows a comparison between estimated and observed fluxes. Net ecosystem exchange closely matched observed values and the low mean bias error (MBE) indicated that there was not a significant bias in estimated NEE. Net methane exchange did not match observations as closely, and the estimates had a small positive MBE. The bias was because the flux driver time series did not account for the influence that variable source area had on hourly fluxes observed by the EC station in 2017; we projected NME to the mean distribution of polygon rims/centers at the site. Daily sums show a similar fit for NEE and an improved fit for NME (Figs. 3c and 3d). Overall, this comparison suggests the estimates of daily NEE and NME presented in this study are reasonable during times when conditions were most similar to the training data used to develop the NN models. The NN models are capable of making projections beyond conditions they were trained on; but error increases as conditions diverge from the training data, and care must be taken to ensure conditions are not fundamentally different than the training data. By limiting our analysis to days with mean soil temperatures at 15 cm > 0 °C, we are restricting our upscaled estimates to conditions that are not drastically different than the 2017 study period.
For an additional check of the plausibility of our daily NEE estimates, we compared them to SMAP L4 global daily 9 km resolution NEE data (Kimball et al. 2021). These satellite estimates of NEE are calculated from 500 m resolution MODIS data, 9 km resolution Soil Moisture Active Passive (SMAP) observations, and 1/4 degree resolution Goddard Earth Observing System Model, Version 5 (GEOS-5) reanalysis data. The SMAP L4 data product is based on a landscape classification map that does not include tundra or wetlands; the LCP ecosystem at Fish Island is classified as shrubland. The 9 km SMAP data are averaged across multiple ecosystems (LCP, shrub tundra, small lakes, etc.), so they are not directly comparable. However, the comparison does help show that our growing season NEE estimates are reasonable (Table A4; Fig. A8). The 95% confidence intervals for our NEE estimates contain overlap with the SMAP L4 estimates of each year. The correlation between the two NEE estimates by year was Spearman ρ = 0.70 and p = 0.18. Overall, daily NEE estimated by our model matched the seasonal patterns of the SMAP data reasonably well (r2 = 0.66) and there was not a significant disparity between the 5 years of overlapping data.
Fig. A8.
Fig. A8. Comparison between NN modelled daily net ecosystem exchange (NEE) in the LCP ecosystem at Fish Island and 9 km resolution SMAP L4 estimates of daily NEE (Kimball et al. 2021). The comparison is limited to June–August because cloudy conditions and limited daylight in September reduced the availability of SMAP data for September.
Table A4.
Table A4. Mean growing season net ecosystem exchange (NEE) in the LCP ecosystem at Fish Island modelled by our NN compared with 9 km resolution SMAP L4 estimates of mean NEE from 2015 to 2019 (Kimball et al. 2021).

Information & Authors

Information

Published In

cover image Arctic Science
Arctic Science
Volume 9Number 3September 2023
Pages: 689 - 709

History

Received: 12 July 2022
Accepted: 12 April 2023
Accepted manuscript online: 19 May 2023
Version of record online: 19 June 2023

Data Availability Statement

All data and code used in this analysis are available on Github: https://github.com/June-Skeeter/FishIsland_Flux_Analysis. High-frequency data are available on request to JS at [email protected].

Key Words

  1. carbon fluxes
  2. polygonal tundra
  3. permafrost
  4. climate change
  5. machine learning

Authors

Affiliations

June Skeeter [email protected]
Department of Geography, The University of British Columbia, Vancouver, BC V6T 1Z2, Canada
Author Contributions: Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Project administration, Visualization, Writing – original draft, and Writing – review & editing.
Andreas Christen
Environmental Meteorology, Faculty of Environment and Natural Resources, Albert-Ludwigs Universität Freiburg, Freiburg, Germany
Author Contributions: Funding acquisition, Supervision, and Writing – review & editing.
Department of Geography, The University of British Columbia, Vancouver, BC V6T 1Z2, Canada
Author Contributions: Funding acquisition, Supervision, and Writing – review & editing.
Greg Henry served as Consulting Editor at the time of manuscript review and acceptance; peer review and editorial decisions regarding this manuscript were handled by Gwyneth MacMillan.

Author Contributions

Conceptualization: JS
Data curation: JS
Formal analysis: JS
Funding acquisition: AC, GH
Investigation: JS
Methodology: JS
Project administration: JS
Supervision: AC, GH
Visualization: JS
Writing – original draft: JS
Writing – review & editing: JS, AC, GH

Competing Interests

The authors declare that there are no competing interests.

Funding Information

Funding for this study was provided by the Canada Foundation for Innovation (CFI, Grant #33600 to AC and GH) and the National Science and Engineering Research Council of Canada (NSERC) through Discovery Grant RGPIN-2017-03958 (AC), Northern Research Supplement Grant RGPNS-503529 (AC), and Accelerator Supplement 507854 (AC).

Metrics & Citations

Metrics

Other Metrics

Citations

Cite As

Export Citations

If you have the appropriate software installed, you can download article citation data to the citation manager of your choice. Simply select your manager software from the list below and click Download.

There are no citations for this item

View Options

View options

PDF

View PDF

Figures

Tables

Media

Share Options

Share

Share the article link

Share on social media

Cookies Notification

We use cookies to improve your website experience. To learn about our use of cookies and how you can manage your cookie settings, please see our Cookie Policy.
×