Seasonality buffers carbon budget variability across heterogeneous landscapes in Alaskan Arctic tundra

Arctic tundra exhibits large landscape heterogeneity in microtopography, hydrology, and active layer depth. While many carbon flux measurements and experiments are done at or below the mesoscale (⩽1 km), modern ecosystem carbon modeling is often done at scales of 0.25°–1.0° latitude, creating a mismatch between processes, process input data, and verification data. Here we arrange the naturally complex terrain into mesoscale landscape types of varying microtopography and moisture status to evaluate how landscape types differ in terms of CO2 and CH4 balances and their combined warming potential, expressed as CO2 equivalents (CO2-eq). Using a continuous 4 year dataset of CO2 and CH4 fluxes obtained from three eddy covariance (EC) towers, we investigate the integrated dynamics of landscape type, vegetation community, moisture regime, and season on net CO2 and CH4 fluxes. EC towers were situated across a moisture gradient including a moist upland tundra, a heterogeneous polygon tundra, and an inundated drained lake basin. We show that seasonal shifts in carbon emissions buffer annual carbon budget differences caused by site variability. Of note, high growing season gross primary productivity leads to higher fall zero-curtain CO2 emissions, reducing both variability in annual budgets and carbon sink strength of more productive sites. Alternatively, fall zero-curtain CH4 emissions are equal across landscape types, indicating site variation has little effect on CH4 emissions during the fall despite large differences during the growing season. We find that the polygon site has the largest mean warming potential (107 ± 8.63 g C–CO2-eq m−2 yr−1) followed by the drained lake basin site (82.12 ± 9.85 g C–CO2-eq m−2 yr−1) and the upland site (77.19 ± 21.8 g C–CO2-eq m−2 yr−1), albeit differences were not significant. The highest temperature sensitivities are also at the polygon site with mixed results between CO2 and CH4 at the other sites. Results show a similar mean annual net warming effect across dominant landscape types but that these landscape types vary significantly in the amounts and timing of CO2 and CH4 fluxes.


Introduction
Arctic tundra is characterized by long non-growing periods punctuated with short growing seasons and a high degree of landscape heterogeneity caused by freeze-thaw cycles. Because these ecosystems are largely inundated or frozen for much of year, soil decomposition is slow, resulting in one of the largest terrestrial reservoirs of labile carbon (C) (Hugelius et al 2014). Comprising only 15% of the global land surface, Arctic tundra contain close to one third of the Earth's terrestrial soil C (∼1500 Pg-C) (Zimov et al 2006, Tarnocai et al 2009, Kirschke et al 2013. As this region is undergoing accelerated warming (Serreze and Francis 2006), carbon dioxide (CO 2 ) emissions are increasing, causing a sink-source transition (Schuur et al 2013, Oechel et al 2014, Natali et al 2015, 2019, Commane et al 2017. In addition to CO 2 , methane (CH 4 ) is an important greenhouse gas (GHG) in permafrost and wetland regions as methanogens thrive in areas with large amounts of anaerobic soil (Garcia et al 2000). Arctic wetlands are responsible for ∼15% of global wetland CH 4 emissions and ∼4% of all global CH 4 emissions (Kirschke 2013). As methanogenesis is a temperature sensitive process (Dunfield et al 1993), it is likely that biogenic CH 4 emissions will increase (Tian et al 2012, Lawrence et al 2015. For this reason, attention to natural sources of CH 4 efflux has increased in Arctic regions. Landscape heterogeneity can contribute to the wide range of estimates in CH 4 emission, ecosystem respiration (ER), and gross primary productivity (GPP) (McGuire et al 2012, Treat et al 2018b. This heterogeneity is characterized by variability in soil moisture regime and develops from freeze-thaw cycles into a patchwork of polygonized tundra, thaw lakes, drained lake basins and moist upland tundra (Webber 1978, Zulueta et al 2011, Liljedahl et al 2016. Polygonized tundra occurs from the common development of ice wedges in the soil column. As ice wedges degrade, inundated C rich low-centered polygons become drained high-centered polygons and drier upland tundra (Liljedahl et al 2016). This shift can cause a decrease in CH 4 emissions and an increase in CO 2 emissions (Martin et al 2018). Degrading ice wedges can also form trough ponds that facilitate water movement, change the microbial community and, by extension, change GHG fluxes (Koch et al 2014, Liljedahl et al 2016. Similarly, thaw lakes form from permafrost thawing and subsequent land subsidence (Jorgensen andShur 2007, Huissteden et al 2011). Thaw lakes drain, forming vegetated drained lake basins (Jorgenson and Shur 2007). This indicates that over short distances, the effects of climate change and controls of emissions can be highly variable due to soil moisture content and the resulting plant and microbial communities that develop.
Recent studies have shown emissions of CO 2 and CH 4 occur well into the non-growing season , Oechel et al 2014, Treat et al 2018a, Arndt et al 2019a. Cold period emissions of CH 4 can account for nearly 50% of the yearly budget, and largely occur during fall shoulder periods when air and surface soil temperatures are below freezing, and subsurface temperatures are around 0 • C (20%-30%) , Commane et al 2017, Taylor et al 2018. This period is referred to as the 'zero-curtain' and is associated with the presence of an unfrozen portion of the active layer during freezing while phase transition occurs (Outcalt et al 1990). This highlights the importance of seasonality, yet there is still a paucity of data reflecting how interseasonal C dynamics vary in terms of landscape type. Quantifying ongoing changes to the pan-Arctic carbon budget is important, but cannot be achieved without understanding how variability in landscape scale climate responses affect emissions. By partitioning Arctic tundra ecosystems into sub-landscapes, the variability in timing and magnitude of C fluxes can be better understood. Using a continuous 4 year dataset (2014-2017) of CO 2 and CH 4 fluxes obtained from three eddy covariance (EC) towers, each in a distinct landscape type, this study aims to quantify the integrated dynamics involved in CH 4 , CO 2 , and CO 2 + CH 4 (expressed as CO 2 equivalent (CO 2 -eq), hereafter referred to as combined C) budgetary contributions due to landscape type, vegetation community composition, and seasonality.

Study area
EC study sites are located on continuous permafrost tundra on the North Slope of Alaska, near Utqiaġvik ( figure 1(a)). The sites include the Barrow Environmental Observatory (US-Beo) ( figure 1(b)), Biocomplexity Experiment South (US-Bes) (figure 1(c)) and a site near the NOAA Earth System Research Laboratory (US-Brw) (figure 1(d)). These sites were chosen for long-term continuous data acquisition as they capture dominant landscape variability of the region. US-Bes is in a drained lake basin containing the wettest soils (table 1), with the water table above the surface for most of the growing season and is dominated by wet sedges and sphagnum moss (Davidson et al 2016a). US-Brw is a moist upland tundra containing the driest soils (table 1) and is dominated by graminoids and lichens (Kwon et al 2006). US-Beo is characterized by ice wedge polygon formations that arise from the freeze thaw cycle (Webber 1978). Due to these polygon formations, US-Beo is a mixed landscape, exhibiting both inundated and drained areas and consists of wet sedge/sphagnum moss dominated areas as well as drier graminoid/lichen dominated areas (Davidson et al 2016a).

Eddy covariance and meteorological data
CH 4 and CO 2 fluxes were estimated at half-hourly intervals from year-round data collected at 10 Hz following the procedures outlined from LI-COR® EddyPro®. A double rotation was applied to the axis rotations of three-dimensional wind speeds according to Wilczak et al (2001) and a block averaging interval was used to define turbulent fluctuations. An in situ/analytic correction, according to Ibrom et al (2007), was applied to the gas analyzer data as the greenhouse gas analyzer (GGA) has a closed path. Quality flags were output within datasets according to Mauder and Foken (2011) and data that did not pass the quality requirements were removed. An internal chamber pressure of ⩾20.67 kPa (155 torr) in the GGA indicates line blockage or instrument failure, and these data were eliminated. Additionally, a turbulence threshold was applied, identifying conditions with insufficient turbulence (indicated by low friction velocity (u * ⩽ 0.1 m s −1 )), and those data were removed in accordance with Reichstein et al (2005). A moving window of 2 weeks was applied and fluxes that were three standard deviations away from the mean were removed as outliers for CH 4 and CO 2 fluxes. EC tower site and instrumentation information can be found in table 1 (Goodrich et al 2016, Arndt et al 2020. Meteorological data were obtained at 30 s intervals and averaged into half-hourly means. Each of the three sites had independent meteorological instrumentation and measurements used in defining the peak thaw and zero-curtain period. Meteorological instrumentation used in analysis included soil water content (SWC) (Campbell Scientific® CS616 Water Content Reflectometer), soil temperature (Omega Engineering™, type T thermocouples), air temperature and humidity (Vaisala, HMP 45), and photosynthetically active radiation (PAR; LI-COR® LI-190R quantum sensor). At each site, soil temperature was measured at 0, −5, −15, and −30 cm from the soil surface and soil moisture probes were inserted from the soil surface to a depth of 20 cm, providing the average soil moisture in the top 20 cm of the soil column. Measurements of soil moisture were taken at one location at US-Brw and US-Bes, and two locations at US-Beo (a high center polygon and polygon trough) for better relief representation then averaged over the growing seasons within the study period. Data from both EC and meteorological instrumentation were collected using datalogger/multiplexer arrays from Campbell Sci-entific® (i.e. CR-3000, CR-23X and AM-1632). Thaw depth was measured weekly during the growing season by probing the land surface to the extent of the active layer along a transect with a small diameter metal rod.

Data/statistical analyses
EC tower footprints were estimated with the analytical footprint model of Korman and Meixner (KM) (2001) using the R package 'FReddyPro' v1.0 (Xenakis 2016). The KM model calculates the density function of the footprint contribution for a two-dimensional area surrounding the EC tower and was used to estimate the landscape area in which 80% of fluxes originated by averaging half-hourly single flux footprints during 2016. WorldView-3 (Maxar Technologies) imagery (1.24 m multispectral resolution) of Utqiaġvik, AK was used to show variability in the WorldView normalized difference water index (NDWI) at each of the three sites. As SWC is measured at one to two locations within each EC tower footprint and site hydrology can be variable, NDWI is used in conjunction with SWC to characterize differences in site moisture regime. The WorldView NDWI is calculated as the normalized difference between the coastal band (R c , 400-450 nm) reflectance and the second near infrared band (R NIR2 , 860-1040 nm) reflectance (equation (1)). This is because it has been shown that the R c and R NIR2 bands show a better soilwater separation than the typical green (510-580 nm) and the first near infrared band (770-895 nm) combinations and is indicative of surface water moisture levels (Maglione et al 2014). The NDWI results in values between −1 and 1 where more positive values represent wetter landscapes. Pixels outside of the footprint of the EC towers as well as those representing structures within the tower footprints were masked for statistical analysis. Before calculating statistics, pixels were aggregated into 3 × 3 pixel grids using the nearest neighbor approach to avoid bias by over sampling and to reduce high variability. A pairwise Wilcoxon rank sum tests was used to compare NDWI values among sites. Additional imagery for NDWI analysis across study period can be found in supplementary information (figure S1 (available online at stacks.iop.org/ERL/16/035008/mmedia)) Daily average fluxes were calculated in R V 3.6.2 (R Core Team 2019) and R Studio software using the 'data.table' package and were calculated with a minimum of 30 half-hourly samples per day to ensure proper representation of the diurnal patterns of CH 4 and CO 2 fluxes. Daily averages of CH 4 and CO 2 fluxes for purposes of the examination of site variability and C budgets are beneficial as they show an accurate representation of the systems while qualifying data complexity and size. Data gaps are unavoidable in the harsh conditions characteristic of Arctic environments. These gaps can be a result of power or network outages as well as instrument failure. Total data coverage is between 61% and 71% depending on the site with the best coverage during the summer and fall periods. Detailed data coverage information by season and annual totals can be found in the supporting information (table S1). Data gaps were filled using random forest machine learning (R package, 'missForest') utilizing a 300-decision tree design. Model validation can be found in supporting information ( figure S2). Comparisons of model validation between the default half-hourly data output and daily averages show models perform much better when using a daily average. This method reduces the 'noise' and is therefore better equipped to inform machine learning processes. The beginning of the growing season was defined as the period where the top five cm of the soil are above 0 • C, ending at the onset of the zero-curtain. The zero-curtain was defined as the period during the fall shoulder beginning when soil temperature of the top 5 cm of the soil are less than 0 • C for 3 or more days, ending when the temperature at −15 cm (roughly the middle of the active layer, table 1) dropped below −0.75 • C for 3 or more days. Non-growing season, as defined here, includes both winter and spring, beginning at the end of the zerocurtain period and ending at the beginning of the growing season.
ER and GPP were partitioned from net ecosystem exchange (NEE) according to Lasslop et al (2010) using the 'REddyproc' package in R (Wutzler et al 2018), as nighttime data is unobtainable during Arctic summer. Temperature response relationships were calculated using a weekly mean to reduce noise and to better represent annual trends. Models showing temperature response curves for CH 4 were calculated using soil temperature, rather than air temperature, as this has been shown to act as a better predictor for CH 4 fluxes (Arndt et al 2019a) (figure S3). All temperature response curves were linearized by log transformation and compared using analyses of covariance (ANCOVA) to test for homogeneity among the regressions with the 'car' R package (Fox and Weisberg 2019) and Q 10 values were calculated from temperature response regressions with the 'respirometry' R package (Birk 2020).

Site moisture regime
NDWI as well as the KM model footprint of each EC tower were used to establish differences in surface water content within each EC tower footprint (figure 2). Using a pairwise Wilcoxon rank sum tests, NDWI was found to be significantly different (p < 0.001) among each of the three sites. The results of the NDWI, showing levels of surface moisture, agreed with soil moisture data at the sites (table 1) with US-Bes showing the wettest conditions (NDWI (mean ± standard error) = −0.079 ± 0.005), US-Beo showing intermediate NDWI levels (NDWI = −0.144 ± 0.003) and US-Brw with the lowest NDWI supporting its position as the driest site in the study (NDWI = −0.204 ± 0.003). Further analyses of NDWI show that imagery acquired 24 July 2016 was representative of site differences and while some variability occurs, these positions are maintained (figure S1) as changes in vegetation community or hydrology at the landscape scale happens over longer periods of time (Liljedahl et al 2016, Arndt et al 2019b.

Carbon budget variability
Seasonal fluctuations (growing season to zerocurtain) of CO 2 budgets are most pronounced at the driest site, US-Brw, and least pronounced at wettest site, US-Bes ( figure 4(a)). Due to the higher zero-curtain CO 2 emissions dampening growing season CO 2 uptake, US-Brw is the weakest CO 2 sink on average. This trend holds with drier sites generally emitting larger amounts of CO 2 during the zero-curtain offsetting much of the uptake during the growing season (US-Brw: 75%; US-Beo: 62%; US-Bes: 25%). Growing season CH 4 emissions are highest at US-Beo and US-Bes ( figure 4(b), p < 0.05). However, zero-curtain CH 4 emissions are roughly equal across all landscapes. This indicates that the percent contribution to local CH 4 budget varies by site during the zero-curtain (US-Brw: 45%; US-Beo: 34%; US-Bes: 32%) and total non-growing season (US-Brw: 56%; US-Beo: 48%; US-Bes: 43%). Site variability in combined C budgets are most pronounced during the growing season, yet still occur, driven by CO 2 emission, during the zero-curtain period (figure 4(c)). As US-Beo is characterized by low centered polygons, it exhibits characteristics of both US-Brw and US-Bes. This accounts for the larger values observed in CH 4 emission relative to US-Brw and the larger values in net CO 2 uptake relative to US-Bes during the growing season. The effect is that the mixed landscape, US-Beo, is the largest mean C contributor (table 2). More detailed information regarding C budget by year/season can be found in supplementary information (table S2).

Temperature response rates
ER temperature response relationships show that air temperature increases in predictive strength as site wetness decreases (figures 6(a)-(c)). In the inter-site comparison of linearized ER regressions, US-Beo and US-Brw show significant similarity (p < 0.05) while US-Bes is significantly different than both US-Brw (p = 0.12) and US-Beo (p = 0.076). Q 10 is highest at US-Beo (3.5), followed by US-Brw (2.5) and US-Bes (2.2). Methane temperature response relationships show similar predictive strength across all sites (figures 6(d)-(f)). Q 10 for CH 4 is again highest at US-Beo (4.6), followed by US-Bes (4.2) and US-Brw (3.1). Contrary to ER, in the inter-site comparison of linearized CH 4 regressions, US-Beo and US-Bes show significant similarity (p < 0.05) while US-Brw is significantly different than both US-Beo (p = 0.11) and US-Bes (p = 0.14).

CO 2 seasonality and annual budget
Seasonality appears to be the dominant factor in annual variation in NEE; however, we find the magnitude of that effect is dependent on local variation in site hydrology ( figure 4). US-Brw exhibited    the highest GPP yet was the weakest mean annual CO2 sink in comparison with the other wetter sites (table 2). Conversely, US-Bes exhibited the lowest GPP yet was the strongest mean annual CO2 sink (table 2). This is because sites with higher summer GPP exhibited larger zero-curtain CO2 emissions. The lower SWC at US-Brw likely leads to a substantially larger portion of the soil column under oxic conditions, therefore supporting aerobic respiration, increasing CO 2 emissions. While it is probable that this constitutes much of the budget disparity, the contribution of CH 4 oxidation to ER in primarily methanogenic areas has been found to be up to 35% (Nielsen et al 2019). As US-Brw contained the deepest active layer (table 1) and largest GPP  (table 2), soils may contain higher amounts of photosynthates and labile C. This can increase methanogenesis rates deeper in the soil column (Dorodnikov et al 2011) and fuel CO2 producing methanotrophs closer to the surface that would be more active under the oxic conditions in drier areas relative to waterlogged areas (Megonigal and Schlesinger 2002). Alternatively, as US-Bes is inundated, methanotrophy is likely substantially lower during the zero-curtain period when the highest amount of CO2 loss is observed. Moreover, the lower GPP at US-Bes may result in more recalcitrant C and thus lower ER rates. The US-Beo site exhibited intermediate growing season GPP and zero-curtain CO2 emission relative to the other two sites likely due to US-Beo exhibiting a mixed landscape regarding the prevalence of drained and inundated areas ( figure 4(a)).
ER temperature response relationships show that US-Beo had the strongest temperature dependence (Q 10 = 3.5) relative to US-Brw and US-Bes (Q 10 = 2.5 and Q 10 = 2.2) ( figure 6). This may be linked to the intermediate soil moisture of US-Beo. Higher SWC can limit the temperature sensitivity of soil respiration in wetland regions due to the restriction of oxygen and thus, aerobic respiration (Chen et al 2018). Alternatively, respiration can be limited by lower SWC through a reduction in microbial mobility and substrate diffusion (Grant and Rochette 1994). The temperature sensitivity of belowground respiration can also be dependent on productivity by providing photosynthates as substrates (Hartley et al 2006). This may contribute to US-Bes having the lowest temperature sensitivity, as this site exhibited the lowest productivity (table 2).
Similar studies of annual CO 2 budgets show that the largest annual net CO 2 loss is seen during the nongrowing season, particularly associated with early winter respiration (Oechel et al 2014, Commane et al 2017, Euskirchen et al 2017, agreeing with data presented here. As early winter respiration comprises a large part of the CO 2 budget, further increases in zero-curtain duration will likely result in winter CO 2 emissions that exceed growing season uptake (Arndt et al 2019b). Each of the sites in this study were found to act as a weak sink for CO 2 , yet others have reported relatively strong annual source signals from similar systems (Commane et al 2017, Eurskirchen et al 2017. This highlights the need for the monitoring and greater representation of mesoscale (⩽1 km) processes in climate projections, particularly at subgrid scales.

CH 4 seasonality and annual budget
Growing season CH 4 emissions are lowest at US-Brw ( figure 4(b)). This is likely due to lower soil moisture increasing the volume of soil experiencing aerobic conditions not conducive to methanogenesis (Garcia et al 2000). The polygon tundra site, US-Beo, had similar CH 4 emissions relative to the inundated US-Bes site despite US-Bes being the site with the higher SWC ( figure 4(b)). This may be related to the vegetation community composition. As US-Bes contains a lower percent cover of sedges (Davidson et al 2016b), US-Beo may produce more photosynthates, like acetate, that could leach into surrounding waterlogged soil in polygonal environments, further fueling methanogenesis via acetoclastic methanogenic pathway (King et al 2002, Dorodnikov et al 2011. Further, sedge density is positively correlated to CH 4 emissions as sedges provide a pathway for CH 4 through the vegetation to the atmosphere (Lai 2009, Andresen et al 2017. Interannual variability in CH 4 emissions is low at US-Bes compared with the other sites (figure 5). As US-Bes is consistently inundated, interannual differences in snow melt and rainfall would have a larger impact on SWC and by proxy, oxygen availability at US-Brw, and US-Beo, possibly explaining this variability.
Contrary to zero-curtain CO 2 emission having site dependent variability, zero-curtain CH 4 emissions are roughly equal across all sites. This may be due to the frozen surface soils, creating an ice 'cap' , and limiting oxygen diffusion into the soil column thereby equalizing oxygen availability and by extension, methanogenesis, across sites. This shows that variability in the growing season, rather than the zerocurtain, may have a stronger impact on annual CH 4 variability across different landscapes. However, zerocurtain emissions may increase as the zero-curtain extends longer into the winter with a warming climate (Arndt et al 2019a). Zero-curtain CH 4 contributions were found to be higher than in previous works due to the length of the study period capturing the range of annual variability . Methane temperature response relationships also indicate that temperature dependence was strongest at US-Beo (Q 10 = 4.6), followed by US-Bes (Q 10 = 4.2), and US-Brw (Q 10 = 3.1) (figure 6). As US-Beo, like US-Bes, contains large amounts of anaerobic soil, temperature sensitivity of methane production would be stronger than US-Brw where soil moisture is a limiting factor.

Annual combined carbon budget
The largest mean combined C emissions were from the mixed landscape US-Beo, exhibiting both inundated and drained areas. These polygonized landscapes comprise close to 65% of the Alaskan coastal plain (Lara et al 2018) and contain both anaerobic areas that produce large quantities of CH 4 as well as drained areas where aerobic respiration can readily occur. US-Beo exhibited the strongest temperature response, for both CO 2 and CH 4 . On this basis, it is possible that further climate change may disproportionately increase C emissions from polygonized landscapes as rising temperatures will support increased production and emission of CO 2 and CH 4 . However, rising temperatures will likely coincide with polygon succession and hydrologic transitions (Liljedahl et al 2016). These hydrologic transitions can significantly alter annual carbon budgets (Kittler et al 2017), stressing the importance of monitoring landscape heterogeneity for in these regions for carbon budget estimation.
The interplay of CO 2 and CH 4 dynamics are affected strongly by both seasonality and by mesoscale landscape variability. Though the variability in summer C emissions is significant among the landscapes studied, zero curtain releases of CO 2 and CH 4 tend to offset site differences. This acts as a buffer to variability and leads to similar annual combined C budgets across the sites studied. However, the differences in timing and magnitude of CO 2 and CH 4 fluxes elevate the importance of mesoscale processes for restricting uncertainty in Arctic model projections. Arctic regions make up the largest portion of uncertainty in climate global climate models (IPCC 2014). Pan-Arctic models are typically run at coarse scales that describe landscape heterogeneity by the dominant landscape. It has been stressed that a higher degree of spatial and temporal coverage is needed (Natali et al 2019) and that representation of wet and dry tundra at a finer scale (⩽4 km 2 ) can result in a threefold reduction in model error (Lara et al 2020). Data presented here demonstrate the need for this improvement, particularly for models that can represent mesoscale landscape heterogeneity and subsequent differences in seasonal carbon emission patterns.

Conclusions
Although the northern coastal tundra region in Alaska continues to be a weak CO 2 sink in all observed landscapes, CH 4 emissions push the region to have a net warming effect on the atmosphere. Data show the site with the largest mean GPP experienced the lowest mean annual CO 2 uptake, while the site with the lowest mean GPP experienced the highest mean annual CO 2 uptake. This is primarily due to zero-curtain CO 2 emissions and indicates that zero-curtain CO 2 emissions are positively correlated with growing season GPP. Despite site variability in growing season CH 4 emissions, zero-curtain CH 4 emissions are nearly equal across sites. This implies that the percent contribution of zero-curtain CH 4 emissions to annual CH 4 budget varies by site and can be larger than previously thought, being as high as 45% of the yearly budget from the zero-curtain period alone and over half of the yearly budget from the total non-growing season (including the zero-curtain). Tundra exhibiting both inundated and drained areas are the largest mean annual combined C source and show a stronger ER and CH 4 temperature response than either largely inundated or drained areas. These results show that local variation in site hydrology, seasonality and interannual variability in regional temperature work in tandem to determine carbon balance. This interaction may be indicative of a variable response under further climate change, yet seemingly lacking the strength to cause strong differences in annual C budgets at this time. As both wetting and drying of Arctic tundra has been reported, differential landscape development in response to climate change and subsequent C budget divergence may occur. Without improved representation of landscape heterogeneity, this potential divergence would likely confound long term global model predictions further, as changes would occur at sub-grid scales.

Data availability statement
The data that support the findings of this study are openly available at the following URL/DOI: https:// arcticdata.io/catalog/view/doi:10.18739/A2X34MS1B.