Deeper snow increases the net soil organic carbon accrual rate in moist acidic tussock tundra: 210Pb evidence from Arctic Alaska

ABSTRACT The net change in the carbon inventory of arctic tundra remains uncertain as global warming leads to shifts in arctic water and carbon cycles. To better understand the response of arctic tundra carbon to changes in winter precipitation amount, we investigated soil depth profiles of carbon concentration and radionuclide activities (7Be, 137Cs, 210Pb, and 241Am) in the active layer of a twenty-two-year winter snow depth manipulation experiment in moist acidic tussock tundra at Toolik Lake, Alaska. Depth correlations of cumulative carbon dry mass (g cm−2) vs. unsupported 210Pb activity (mBq g−1) were examined using a modified constant rate of supply (CRS) model. Results were best fit by two-slope CRS models indicating an apparent step temporal increase in the accumulation rate of soil organic carbon. Most of the best-fit model chronologies indicated that the increase in carbon accumulation rate apparently began and persisted after snow fence construction in 1994. The inhomogeneous nature of permafrost soils and their relatively low net carbon accumulation rates make it challenging to establish robust chronologic records. Nonetheless, the data obtained in this study support a decadal-scale increase in net soil organic carbon accumulation rate in the active layer of arctic moist acidic tussock tundra under conditions of increased winter precipitation.


Introduction
Understanding the arctic carbon cycle today and its likely behavior in the future is essential to recognizing shortand long-term responses and feedbacks to the climate system. The impact of the changing Arctic on the climate system can be substantial because arctic ecosystems harbor as much as 50 percent of the total terrestrial soil carbon pool (Vitharana et al. 2017). Changes in soil organic carbon inventories in arctic regions are expected as a direct response to warmer summers (Chapin and Shaver 1985;Elmendorf et al. 2012;Gallego-Sala et al. 2018) and greater moisture availability (Zhang et al. 2013), including increased winter snow precipitation (Welker, Fahnestock, and Jones 2000;Blanc-Betes et al. 2016;Leffler et al. 2016;Jespersen et al. 2018). However, precipitation predictions are uncertain across the Arctic (some regions are predicted to have more and some are predicted to have less), and this variation is reflected in the annual fluctuations in snow cover and snow water equivalents reported for the past ten years by the annual National Oceanic and Atmospheric Administration Arctic Report Cards (see NOAA 2019).
Increases and decreases in snowfall and winter snowpack can lead to different ecosystem impacts (Bintanja and Selten 2014;Bintanja and Andry 2017). The effects of these changes on soil organic carbon inventories depend on soil carbon cycling processes (Schaeffer et al. 2013;Ricketts et al. 2016;Vitharana et al. 2017) and soil active layer depth (Pattison and Welker 2014;Jespersen et al. 2018). Improved understanding of the mechanisms driving the arctic carbon cycle continues to be at the forefront of research on the role of the "New Arctic" in the global climate system, especially processes that drive winter carbon emissions Lupascu et al. 2018;Natali et al. 2019).
Release of stored arctic carbon into the atmosphere, rivers, and ocean is becoming increasingly apparent (Czimczik and Welker 2010;Mishra et al. 2017;Nowinski et al. 2010;Lupascu et al. 2014Lupascu et al. , 2018Csank et al. 2019). It is being driven by deepening of the active layer as permafrost thaws and by coastal and riverine erosion causing shoreline retreat and export of carbon-rich sediment and water into the Arctic Ocean (Martens et al. 2018). In contrast, arctic tundra is also becoming more productive as seen by increases in the Normalized Difference Vegetation Index values, due in part to greater leaf area associated with shrub density increases (Tape et al. 2012;Leffler et al. 2016;Mekonnen, Riley, and Grant 2018). The net effect of these accelerated carbon losses and productivity gains on the arctic carbon cycle is difficult to quantify because it requires a full picture of the annual CO 2 and CH 4 emissions (Blanc-Betes et al. 2016) and the export of dissolved organic and inorganic carbon by surface water flow (Lehn et al. 2017;Csank et al. 2019), as well as gains in soil organic carbon and biomass due to increases in net primary productivity (Welker, Fahnestock, and Jones 2000;Matamala et al. 2017). It is difficult to depict changes in carbon inventories in permafrost soils because we lack a baseline for net soil carbon accrual rates and how they may be affected by winter precipitation scenarios.
In this study, we explored a 210 Pb-based chronometric approach to obtain an integrative measure of the temporal change in net soil carbon accumulation rate in the active layer organic horizon of moist acidic tussock tundra in response to the construction of a snow fence near Toolik Lake, Alaska, at 68°38′ N, 149°38′ W; elevation 760 m ( Figure 1). We measured and compared profiles of soil organic carbon and 210 Pb (half-life = 22.3 a) along with 7 Be (half-life = 53.3 d), 137 Cs (half-life = 30.2 a), and 241 Am (half-life = 432.2 a) in a set of active layer cores from within and outside the area influenced by a twenty-two-year snow fence experiment. This was a reconnaissance study conducted with a relatively small set of cores to assess whether this approach could give meaningful results over a decadal timescale, given the known complications of 210 Pb-based soil chronology Münnich 1989, 1991;Kaste, Heimsath, and Bostick 2007;Kaste et al. 2011;Landis, Renshaw, and Kaste 2016;Olid et al. 2016;Gonzalez-Meler et al. 2018). Our objective was to ascertain the impact of deeper winter snow accumulation (caused by the snow fence) on the net rate of soil organic carbon accumulation. We assumed that all soil organic carbon in the active layer profile was produced from plant and microbial inputs and that 210 Pb deposited from the atmosphere could be used as a proxy to estimate the net carbon accumulation rate in the organic-rich horizon of the active layer. Key assumptions in this approach were (1) atmospheric Pb deposited onto the tundra surface was fully retained within the active layer, (2) 210 Pb was strongly sorbed to organic matter where it has limited mobility (Landis, Renshaw, and Kaste 2014), and (3) soil organic matter decomposition rate was low in the cold arctic environment. Thus, the profile of 210 Pb deposited onto the organic carbon could allow a quantitative model of the net accumulation rate of soil organic carbon on the timescale of 210 Pb decay. We used 7 Be as an indicator of near-surface 210 Pb mobility (Landis, Renshaw, and Kaste 2016) and made simple modifications to our CRS models to accommodate this effect and also to account for changes in net soil organic carbon accumulation rate and 210 Pb depositional flux. We compared our 210 Pbbased chronologies with the depth profiles of bomb-era fallout radionuclides 137 Cs and 241 Am.
Applications of fallout radionuclides in studies of sediment, soil, and peat accumulation ( 7 Be, 137 Cs, 210 Pb, 241 Am) have been explored by a number of authors (Robbins and Edgington 1975;Münnich 1989, 1991;Urban et al. 1990;Kaste, Heimsath, and Bostick 2007;Kaste et al. 2011;Abril and Gharbi 2012;Kaste 2014, 2016;Shotyk et al. 2015;Olid et al. 2016;Davies et al. 2018;Gonzalez-Meler et al. 2018). The prior studies most relevant to the current study are those dealing with accumulation of soil and peat. There is a growing consensus reflected in the most recent studies that postdepositional downward mobility of 210 Pb into the accumulating organic horizon of soil and peat generally results in a nonexponential profile of 210 Pb activity with depth (e.g., Landis, Renshaw, and Kaste 2016;Olid et al. 2016). The downward motion of 210 Pb is attributed to advective transport of colloidal organic particles with adsorbed 210 Pb, and this generally results in maximum 210 Pb activity at a depth of several centimeters to tens of centimeters below the surface vegetation and litter layer (Landis, Renshaw, and Kaste 2016). The depth profile of 7 Be is a good indicator of the depth of 210 Pb transport because these two radionuclides behave similarly in such systems (Landis, Renshaw, and Kaste 2016). Various numerical and analytical models have been proposed to accommodate near-surface 210 Pb mobility so that chronometric information can be obtained regarding mass accumulation rates and calendar dates as a function of mass depth. Most such models assume steady-state behavior with diffusion-like transport of 210 Pb in the near-surface (typically over a centimeter to decimeter depth range), and derived accumulation rates and calendar dates may be verified with the occurrence of distinct time horizons such as the peak activities of the bomb-pulse radionuclides 137 Cs and 241 Am that occurred in 1963-1964 (Kaste, Heimsath, and Bostick 2007;Landis, Renshaw, and Kaste 2016;Olid et al. 2016). Concordant chronologic results may be obtained between 210 Pb and other fallout radionuclides below the depth of maximum 210 Pb activity in soil profiles (Landis, Renshaw, and Kaste 2016).
Prior to this study, we anticipated that the consequences of deeper snow in winter could lead to lower soil carbon accumulation, due in part to higher rates of winter CO 2 emissions and accelerated summer decomposition (Schimel, Bilbrough, and Welker 2004), and that any changes in vegetation would be insufficient to counter these losses with soil carbon inputs via litter deposition, root growth and exudate inputs, and microbial growth. In contrast, simple models of our 210 Pb data from soil cores beneath artificially deepened winter snow are consistent with a substantial increase in the net rate of soil carbon accumulation relative to cores from ambient snow depth conditions outside the influence area of the snow fence. Our data also indicated a rapid response of soil organic carbon accumulation to increased snow depth, because the increase in net soil organic carbon accumulation rate in cores for which it could be reasonably quantified coincided with the calendar date of snow fence construction.

Field site
As part of the ITEX (International Tundra Experiment; Henry and Molau 1997; M. D. Walker et al. 1999), a 2.8 m-high × 60 m-long snow fence was built at the study location in 1994, perpendicular to the prevailing winter southerly winds, to artificially increase leeward snow accumulation to simulate the anticipated increase in winter precipitation with winter climate warming. Within the snow fence experimental plot, leeward snow drift accumulation zones were defined relative to ambient annual snow depth of ~0.5 m as deep (2-3 m), intermediate (0.5-2 m), and low (<0.5 m) (Figure 2; see also Supplementary Figure S1). The entire experimental plot has heterogeneous microtopography and an ~5° eastward slope. The active layer thickness of moist acidic tussock tundra around Toolik Lake averages ~30 cm during the three-to four-month growing season and reaches a maximum thaw depth of 45 to 50 cm by late August. The moist acidic tussock tundra is an ombotrophic bog, and the active layer consists of a 5-to 20-cm-thick organic layer overlying a mineral horizon derived from weathered glacial drift. Total annual precipitation is near 300 mm (~200 mm summer and ~100 mm winter), with snow cover lasting from September to May.
Soil surface temperatures in the snow fence experimental plot reach 10°C to 15°C during summer months and decline to about −7°C in late winter under the deepest part of the snow drift, whereas temperatures as low as −30°C occur under ambient late winter conditions (M. D. Walker et al. 1999;Pattison and Welker 2014). During the early winter, before substantial snow accumulation, small amounts of leaves and windblown organic debris are deposited leeward within 25 m of the snow fence (deep snow zone; Figure 2; M. D. Walker et al. 1999). There is typically a two-to three-week delay in the timing of completion of the spring snowmelt in the deeper drift area leeward of the snow fence (M. D. Walker et al. 1999). Soils are acidic (pH < 5), coarse-loamy, and poorly drained. Vegetation is dominated by tussockforming cottongrass sedge (Eriophorum vaginatum) and intertussock moss (Sphagnum sp.), with wet sedge (Carex bigelowii) and occasional shrubs including dwarf birch (Betula nana) and tealeaf willow (Salix pulchra), in a landscape that is typical of moist acidic tussock tundra across the Alaskan North Slope (M. D. Walker et al. 1999;Leffler et al. 2016). Vegetation cover across the snow fence experimental plot was generally homogeneous, except in the deep snow accumulation zone where the abundance of sphagnum moss and wet sedge was much higher and the abundance of cotton-grass sedge and dwarf birch was much lower (Miller 2014). Standing water is frequently present in the deep snow accumulation zone because it has a slight topographic depression.

Sample collection
Cores were collected in the intertussock areas dominated by sphagnum moss growth at the surface. Soil cores (4.8 cm diameter × 30 cm depth) were taken between 10 a.m. and 3 p.m. in early August 2016, including five cores from within the snow fence plot (one core from deep snow depth zone ~10 m leeward of the snow fence, three cores from intermediate snow depth zone ~35 m leeward of the snow fence, and one core from low snow depth zone ~65 m leeward of the snow fence) and three cores from outside the snow fence plot (intended as controls representing ambient snow conditions; Figure 3; see Supplementary Figure S2 for photographs of typical core sampling sites). Cores were immediately wrapped in Al foil, sealed in plastic bags, and kept frozen until prepared for analysis.

Sample preparation
Frozen cores were sectioned by using a band saw in a −20°C cold room into 1-cm depth increments that were weighed before drying and then air-dried at 60°C to a constant weight for the direct determination of soil water content and bulk dry density (g cm −3 ) at the University of Illinois at Chicago. Core segments were homogenized after drying by using a high-speed rotary blender.
Gamma spectrometer efficiencies were calibrated using measurements of certified radionuclide standards, including DL-1A (CANMET 2017) for 210 Pb and other U-and Th-series nuclides and NIST SRM-4357 (NIST 1997) for 137 Cs. We assumed that the efficiency of the 241 Am peak at 59.5 keV was equal to that determined for the 234 Th peak at 63.3 keV, and we estimated the detector efficiency at the 7 Be energy by linear interpolation between the nearest large peaks of the U-and Th-series nuclides in the DL-1A reference spectra. Specific activities and one-sigma counting errors were calculated using standard statistical techniques (Mook and DeVries 2001) and are reported in millibecquerels per gram dry mass (mBq g −1 ). Reported activities were corrected for detector background and decay from the time of sampling to the start time of the gamma spectrometry measurement. Limits of quantification for radionuclides were defined as three times the standard deviation of the background under the peak used for the activity quantification, and values less than the limit of quantification are not reported.

Carbon concentrations
Total carbon contents (weight percentage, dry basis) were determined by using a zero-blank elemental analyzer (Costech Analytical, EC4010; Valencia, CA) at the University of Illinois at Chicago. Total carbon concentrations are assumed to equal total soil organic carbon, because inorganic carbon solid phases are unstable at the low pH of the moist acidic tundra active layer. No further acidification was performed during sample preparation. Precision of carbon concentrations was ±1 to 2 weight percentage, based on replicate analyses of samples and standards.

Pb data
Profiles of unsupported 210 Pb activity in soil cores were analyzed with a modification of the constant rate of supply (CRS) model (Appleby and Oldfield 1978;Robbins 1978), to obtain estimates of apparent net carbon accumulation rates (mg cm −2 yr −1 ) from best-fit linear regressions of cumulative carbon dry mass (mg cm −2 ) vs. ln(A 0 /A Z ). Because the 1994 construction of the snow fence caused substantial increases in the annual accumulation of winter precipitation across the experimental site, we considered that 210 Pb deposition rates also could have increased leeward of the snow fence after 1994.
Sediment or soil mass accumulation rates can be calculated using the CRS model, assuming a constant flux of 210 Pb to the sediment, and this model allows changes in sedimentation rate to be identified graphically. In this case, the integrated activity of unsupported 210 Pb in the profile is expressed as where A Z is the cumulative unsupported 210 Pb activity (Bq cm −2 ) below sediment depth Z and A 0 is the cumulative excess 210 Pb (Bq cm −2 ) for the entire profile where λ is the decay constant of 210 Pb (0.03108 yr −1 ), C Z is the unsupported 210 Pb activity (mBq g −1 ) in the core section at depth z, and ρz is the dry bulk density (g cm −3 ) of the core section at depth z. Unsupported 210 Pb activity was obtained by subtracting the mean value of supported 210 Pb as determined from deeper core sections exhibiting approximately constant total 210 Pb vs. depth. Supported 210 Pb activities in seven of the cores were calculated from four to ten deeper sections and ranged from 25 ± 4 to 37 ± 8 mBq g −1 . An anomalously high value of 54 ± 7 mBq g −1 was determined for core INT-3 Open circles indicate cores reflecting ambient winter snow conditions; filled circles indicate cores for which snow depth was affected by snow fence. Snow fence is approximately 600 m south of Toolik Lake at ~50 m elevation above lake level. Prevailing wind direction is southerly.
in which only two deeper sections were available for this calculation.
The approximate age t (years) of material at a given depth is given by and the net mass accumulation rate (r) (mg cm −2 yr −1 ) is given by where λ is the decay constant of 210 Pb (0.03108 yr −1 ) and M is the slope of the line of cumulative dry mass (g cm −2 ) vs. ln(A 0 /A z ). For the purpose of this study, we transformed cumulative dry mass to cumulative carbon mass by taking the product of the dry bulk density (g cm −3 ) and the weight percentage carbon for each 1-cm core segment and then integrating downward from top to bottom of the core profile. Simple modifications of the CRS model to account for step changes in carbon accumulation and 210 Pb flux are presented and discussed below with reference to the data obtained in this study.

Calculation of unsupported 210 Pb, 137 Cs, and 7 Be inventories
The inventories of unsupported 210 Pb, 137 Cs, and 7 Be in each core were calculated as follows: where I is the total inventory of each radionuclide (mBq cm −2 ), ρ i is the bulk dry density of core segment i, A i is the specific activity of each radionuclide in core segment i (mBq g −1 ), and Δz i is the thickness of core segment i (cm), and the sum is over all core segments.

Calculation of 210 Pb depositional flux
The depositional flux of 210 Pb for each core location is calculated as follows: where Q is the depositional flux of 210 Pb (mBq cm −2 yr −1 ), λ is the decay constant of 210 Pb (0.03108 yr −1 ), and I is the unsupported 210 Pb inventory (mBq cm −2 ) as calculated by Equation (6).

Results and discussion
The analytical results of this study are compiled in Supplementary Information Table S1. Depth profiles for specific activities (mBq g −1 ) of unsupported 210 Pb, 7 Be, 137 Cs, and 241 Am are shown in Figure 4. Soil organic C concentrations in cores from the study site ranged from a maximum of 40 to 47 weight percentage in the organic horizon (uppermost 4-15 cm of core). Vegetation (mainly sphagnum moss) and a minor amount of plant litter comprise the top 2 to 4 cm of most cores (Supplementary Figure S3), with the remaining organic layer consisting of variably compacted, partially humified plant material with minimal mineral content. Profiles of dry bulk density and weight percentage carbon vs. depth show that the organic layer has dry bulk density in the range of 0.05 to 0.2 g cm −3 and >30 weight percentage carbon (Supplementary Figures S4 and S5). Organic carbon decreases sharply below the organic layer in a 2-to 4-cm-thick transitional layer to about 2 to 5 weight percentage in the underlying fine-grained mineral horizon that has dry bulk density in the range from 0.8 to 1.4 g cm −3 .
The activity of 7 Be was measured only in five of the cores (those counted within a few months of collection). and in these cores it was found to have declining activity with depth in the upper 2 to 4 cm ( Figure 4). The disappearance of 7 Be with depth generally corresponded to the observed maxima of unsupported 210 Pb that ranged from 399 to 537 mBq g −1 in the control (CTL) cores, 488 mBq g −1 in the reduced snow depth (LOW) core, 476 to 681 mBq g −1 in the intermediate snow depth (INT) cores, and 1,305 mBq g −1 in the deep snow depth (DEEP) core. Unsupported 210 Pb exceeded supported 210 Pb in all cores by a factor of ten or more at the depth of maximum 210 Pb activity. The activities of 137 Cs were measured in all cores and each core had a welldefined maximum activity. Maximum activities occurred over a range of core segment depths from 0 to 1 cm up to 9 to 10 cm. Activities of 241 Am were detectable in all but one core profile.

Carbon accumulation rates determined from 210 Pb profiles
To obtain a net rate of soil organic carbon accumulation for each core profile, we evaluated the linear correlations of cumulative organic carbon (g cm −2 ) with ln(A 0 /A Z ). We refer to this type of data analysis henceforth as a CRS plot. The slope of a best-fit linear regression between these two parameters in a typical CRS plot yields the net organic carbon accumulation rate. Changes in slope on a CRS plot may correspond to changes in net carbon accumulation rate. Carbon accumulation rates (mg cm −2 yr −1 ) derived from linear regression slopes are listed in Table 1 with errors, along with the total inventories of unsupported 210 Pb, 137 Cs, and 7 Be and the steady-state equivalent depositional fluxes of 210 Pb.
In order to assist our interpretation of the observed changes in slope observed in CRS regression lines for these cores, we did a simple numerical experiment (described in Supplementary Information Addendum and Table S2) to demonstrate the effects of step changes in arbitrary values of the net carbon accumulation rate and the 210 Pb depositional flux on the chronologic records obtained from CRS plots. We evaluated hypothetical changes in these parameters corresponding with the time of construction of the snow fence in 1994; these are portrayed in a demonstrative plot (  For our measured core profiles, we estimated carbon accumulation rates from the slopes of linear regressions and considered the 210 Pb depositional flux as a fitting parameter that would allow the observed changes in slope to begin at a calendar date of 1995. This resulted in a fairly unique solution in terms of both carbon accumulation rate and the apparent change in 210 Pb flux. In order to account for the expected near-surface downward mobility of 210 Pb (Landis, Renshaw, and Kaste 2016), we excluded from our linear regressions any core segment in which 7 Be was detected. In cores where 7 Be was not measured, we excluded the top two core segments (top 2 cm) from the CRS plot regression analyses on the assumption that these were most likely affected by downward 210 Pb mobility. The chronologies implied by the resulting modified CRS models are considered accurate only below the depth limit of significant downward 210 Pb mobility (Landis, Renshaw, and Kaste 2016). Individual modified CRS plots prepared for each of the measured core profiles are shown in Figure 6. Three "control" cores (CTL prefix in sample ID) were collected. These sites were designated in the original experimental design to represent ambient snow conditions; that is, outside the leeward snowdrift influence of the snow fence. Two of these cores  showed good linear correlations over the entire depth range in which unsupported 210 Pb was present (Figure 6), consistent with a constant 210 Pb depositional flux and constant net organic carbon accumulation rates of 1.5 ± 0.1 to 1.7 ± 0.1 mg cm −2 yr −1 over the past ~100 years (Table 1). In these two cores, 241 Am was detected in the core segment that would be Figure 6. Modified CRS plots for measured core profiles showing cumulative organic carbon (g cm −2 ) vs. ln(A 0 /A Z ). Note differences in ranges of x-and y-axis values. Dashed lines are best-fit linear regressions. Numerals adjacent to data points correspond to calendar years calculated from ln(A 0 /A Z ). Net carbon accumulation rates (mg cm −2 yr −1 ) were obtained from slopes of regression lines. Data points having 7 Be (detected or assumed) are excluded from regressions. expected to contain bomb-pulse radionuclides based on its estimated calendar date, supporting the validity of the chronologic record estimated from the bulk values of each core segment. However, the 137 Cs peaks in these cores were found in segments well above the depth of 241 Am, indicating likely mobility of 137 Cs in the low-pH pore waters of the active layer. In core CTL-4, 241 Am was also detected in segments having calendar dates after the bomb-pulse peak. The third control core CTL-5 was, contrary to expectation, apparently influenced by the snow fence because its modified CRS plot showed a substantial increase in slope after 1994 ( Figure 6). This core profile indicated that the carbon accumulation rate increased substantially after the construction of the snow fence (Table 1). The location of CTL-5 was ~65 m NNW of the western edge of the snow fence. Wind directions during winter at Toolik are generally southerly (katabatic winds), with occasional southeasterly or westerly trends (S. Sturm and Stuefer 2013). Occasional southeasterly winds during winter may thus have caused consistent snowdrift pile-up at CTL-5 following the snow fence construction. The unsupported 210 Pb inventory in CTL-5 is higher than that in both CTL-1 and CTL-2 (Table 1). This further supports our interpretation of CTL-5 as being influenced by enhanced snow accumulation that was caused by the snow fence, because additional snow accumulation would be likely to cause an increase in the 210 Pb depositional flux.
Most of the cores affected by enhanced snow accumulation leeward of the snow fence (i.e., CTL-5, LOW-2, INT-1, and INT-5) had modified CRS plots best described by two-slope models having a distinct change in slope at depth, with concomitant increases in net carbon accumulation rate and 210 Pb depositional flux ( Figure 6). Two of the cores (INT-3 and DEEP-4) were difficult to fit with modified CRS models because they had insufficient numbers of core segments below the depth of measured 7 Be activity with which to define values for two different slopes. For the cores having successful two-slope modified CRS models, the carbon accumulation rates were given by the slopes of the linear regressions, and the 210 Pb fluxes were treated as a fit parameter that was adjusted to obtain a calendar date of 1995 for the change in slope. All of the two-slope models required increases in 210 Pb flux (ranging from +82 to +310 percent) to fit a calendar date of 1995 for the change in slope. This implies that the increase in winter snow accumulation leeward of the snow fence was accompanied by a roughly proportional increase in 210 Pb deposition or possibly that there had been an overall increase in 210 Pb flux deposited by summer precipitation after 1994. A gradual increase in summer 210 Pb flux in the tundra region may be likely as the thaw depth increases to allow more release of 222 Rn from the active layer.
The portions of the successful two-slope modified CRS models corresponding to the period before the change in slope indicated net carbon accumulation rates ranging from 0.7 ± 0.1 to 1.8 ± 0.3 mg cm −2 yr −1 , which bracket those obtained from the CRS models for the two control cores CTL-2 (1.5 ± 0.1 mg cm −2 yr −1 ) and CTL-4 (1.7 ± 0.1 mg cm −2 yr −1 ). The slopes of the two-slope modified CRS models beginning in 1995 (i.e., after snow fence construction) yielded uniformly larger carbon accumulation rates ranging from 2.8 ± 0.2 to 4.6 ± 0.3 mg cm −2 yr −1 (Table 1), roughly two to three times as much as those of the control cores.

Inventories of unsupported 210 Pb, 137 Cs, and 7 Be
Total inventories and annual fluxes of unsupported 210 Pb, 137 Cs, and 7 Be in the cores were calculated as described in the Methods section and are listed in Table 1. The total inventories of unsupported 210 Pb range from 79 ± 2 to 481 ± 5 mBq cm −2 . The values for control cores CTL-2 and CTL-4 (105 ± 3 and 79 ± 2 mBq cm −2 , respectively) are in good agreement with unsupported 210 Pb inventories obtained previously for sediments from Toolik Lake (75 mBq cm −2 ; Cornwell 1985) and for intertussock and ridgetop tundra soils (95 and 104 mBq cm −2 , respectively) from the Imnavait Creek watershed, ~16 km east of Toolik Lake (Cooper et al. 1991). Mean annual 210 Pb depositional fluxes obtained from our Toolik area cores ranged from 2.5 to 15.0 mBq cm −2 yr −1 . Inventories and fluxes of 210 Pb are lowest in control cores CTL-2 and CTL-4 (3.3 ± 0.1 and 2.5 ± 0.1 mBq cm −2 yr −1 , respectively) and moderately higher in most cores affected by enhanced snow accumulation (3.8 ± 0.1 to 4.5 ± 0.1 mBq cm −2 yr −1 ), with an anomalously high value measured in core DEEP-4 (15.0 ± 0.2 mBq cm −2 yr −1 ), which is further discussed below. The general trend of these values is consistent with that expected from increased snow depth, which appears to be accompanied by a higher depositional flux of atmospheric 210 Pb. However, because snow fence construction occurred in 1994, only one 210 Pb half-life (approximately twenty-two years) had elapsed between snow fence construction in 1994 and our core sampling in 2016, so the observed total inventories and fluxes mainly reflect prevailing conditions prior to 1994. Using the modified CRS models that we generated for this study, we verified that the best-fit post-1994 increases in 210 Pb depositional flux reproduced the observed increases in unsupported 210 Pb inventories relative to those of the control cores. Calculated pre-1994 210 Pb depositional fluxes obtained from our best-fit modified CRS model parameters ranged from 2.1 to 2.9 mBq cm −2 a −1 , with a mean value of 2.6 ± 0.4, which is in reasonable agreement with the 210 Pb depositional fluxes obtained from our measurements of unsupported 210 Pb inventories for the control cores (2.5 ± 0.1 and 3.3 ± 0.1 mBq cm −2 a −1 ).
The 210 Pb depositional fluxes that we derived from our measurements (Table 1) are in good agreement with a range of values from other high-latitude studies. A study of fifty-five sediment cores from five small arctic lakes near our Toolik field site yielded fluxes ranging from 2.9 to 7.7 mBq cm −2 yr −1 (Fitzgerald et al. 2005). A global compilation of 210 Pb depositional fluxes showed that the mean value at nine sites in 60° to 70° N latitude was 3.7 ± 3.5 mBq cm −2 yr −1 (Baskaran 2011). The 210 Pb flux measured in Greenland ice cores Summit and Dye 3 ranged from 0.4 to 2.5 mBq cm −2 yr −1 (Dibb 1992).
The inventories of 137 Cs in our cores ranged from 45 ± 1 to 104 ± 1 mBq cm −2 . Values obtained for the control cores CTL-2 and CTL-4 (83 ± 1 and 57 ± 1 mBq cm −2 , respectively) were within the range of those affected by the snow fence (45 ± 1 to 104 ± 1 mBq cm −2 ). These values agree well with 137 Cs inventories measured in tundra samples collected in 1990 from the Imnavait Creek watershed (Cooper et al. 1995) that had a range of 42 to 97 mBq cm −2 (reported range of 77 to 177 mBq cm −2 , corrected for twenty-six years of decay from 1990 to 2016). Cooper et al. (1995) concluded that there was no evidence for widespread downslope 137 Cs transport in the Imnavait Creek watershed, although they noted that the elevated 137 Cs inventories in some of the Toolik Lake sediment cores near lake inlets indicated some apparent mobility. In several of our cores (INT-1,   137 Cs peaks were concordant with 210 Pb chronologies, but most of the other cores showed apparent evidence for 137 Cs mobility generally toward the surface. This may be attributable to the low pH of pore waters (<5) along with possible advective transport of dissolved Cs through pore water flow or surface runoff during the spring snowmelt season and possible uptake by plant roots. Further study of 137 Cs activities in tundra pore waters and streams around Toolik Lake may provide useful constraints on the possible mechanism(s) and extent of 137 Cs mobility in this system.
The inventories of 7 Be in the cores affected by the snow fence ranged from 0.8 ± 0.1 to 16.2 ± 0.2 mBq cm −2 , in remarkably good agreement with the range of 0.6 to 15.5 mBq cm −2 measured in tundra samples of the Imnavait Creek watershed (Cooper et al. 1991). The range of 7 Be inventories in snowpack cores collected from that watershed in 1989 and 1990 was 1.0 to 8.1 mBq cm −2 , similar to that measured in tundra samples. These inventories represent a subannual time period because of the relatively short (53 d) half-life of 7 Be.
The anomalously high 7 Be inventory of 16.2 ± 0.2 mBq cm −2 in core DEEP-4 is correlated with its anomalously high unsupported 210 Pb inventory of 481 ± 5 mBq cm −2 (Table 1). This may indicate that short-term lateral advective water flow carrying particles having adsorbed 7 Be and 210 Pb into the local topographic lows (i.e., incipient thermokarst) that developed in the area of deepest snow accumulation. An additional factor that could help explain these anomalous values in core DEEP-4 is the observation of windblown leaves and surficial organic debris deposited leeward within 25 m of the snow fence in the deep snow zone (Figure 2; M. D. Walker et al. 1999).

Summary and conclusions
The principal objective of this study was to evaluate the effect of deeper winter snow accumulation on the net accumulation rate of soil organic carbon in the active layer of moist acidic tussock tundra, following the construction of a snow fence in 1994 as part of the ITEX (Jones et al. 1998;Welker, Fahnestock, and Jones 2000;Welker et al. 2005;Jespersen et al. 2018). We explored the application of 210 Pb as a chronometer to constrain decadal-scale carbon accumulation rates in the active layer. Based on the values we obtained from our modified CRS models for net soil organic carbon accumulation rates in the control cores and the pre-1994 portions of the other cores, we estimated the average background (i.e., ambient snow conditions) net soil organic carbon accumulation rate in the active layer of the study area at ~14 g C m −2 yr −1 . Increases in winter snow depth over a period of twenty-two years (since construction of the snow fence in 1994) resulted in apparent increases of this rate at this site by ~90 to 310 percent to rates of ~26 to 44 g C m −2 yr −1 . This represents a substantial increase in the soil carbon sink strength, as predicted by a number of climate warming models (e.g., Charman et al. 2013;Jiang et al. 2016;Chaudhary, Miller, and Smith 2017). Similar values of carbon accumulation rates were found in a chronologic study of soil profile TFS2 sampled ~50 m west of our study area (at a location where the dominant vegetation consisted of Sphagnum capillifolium, Aulacomnium turgidum, and Salix reticulata), which showed a substantial increase in carbon accumulation rate from ~14 g C m −2 yr −1 to ~43 g C m −2 yr −1 corresponding to an apparent ecosystem change from minerotrophic to oligotrophic conditions after the Little Ice Age at ~1850 CE (Taylor et al. 2019).
These results suggest that though deeper snow leads to elevated winter CO 2 emissions, the positive effects of deeper snow on shrub growth, aboveground biomass and net CO 2 fixation, and presumably greater root biomass collectively overcome the gaseous losses of C, leading to the positive net effect of increasing active layer soil carbon accumulation rates. This positive and negative feedback contest between ecosystem C losses (gaseous emissions) and gains (combinations of C fixation, biomass accumulation, and soil C formation) has been at the center of some of the earliest studies of tundra responses to changes in conditions with ITEX (Welker et al. 1997). The findings in this study suggest that over decades of deeper snow in winter, soil carbon accumulations of tussock tundra will be appreciably greater due in large part to vegetation changes associated with shrub growth and elevated soil carbon inputs (M. Sturm et al. 2005;Tape et al. 2012;Lynch et al. 2018) and lower organic matter decomposition rates as soil water content increases with snow (Blanc-Betes et al. 2016), even though CO 2 emissions in winter are higher.
However, abiotic and biotic interactions may change the magnitude, extent, redistribution, and even longerterm net effect of this increased soil carbon sink. For instance, thawing of previously frozen soils may initially stimulate rates of microbial degradation of organic carbon, and as accessible carbon pools are depleted over time and soil conditions change, microbial communities may shift toward more effective decomposition and metabolic strategies (Ricketts et al. 2016;Melillo et al. 2017). Competing rates of soil methanogenesis versus respiration, which respectively emit CH 4 and CO 2 to the atmosphere, depend in large part on whether the soils are wet or dry (Blanc-Betes et al. 2016). These changes in soil carbon accumulation will not occur in isolation from other climate change-induced effects. Permafrost thaw may change the ecohydrology of the system (Blanc-Betes et al. 2016;Throckmorton et al. 2016;Jespersen et al. 2018), and an even greater dominance of the landscape by shrubs with more leaf area will continue to increase carbon fixation (Leffler et al. 2016;Lynch et al. 2018), leading to greater carbon allocation to roots. Denser shrub canopies may cause shading and soil cooling, countering greater soil carbon emissions (Mekonnen, Riley, and Grant 2018) and amplifying the biocomplexity of the tussock tundra's response to climate change (D. A. Walker, Epstein, and Welker 2008;Flower et al. 2019).
The multidecadal net effect of the complex ecohydrologic interactions in the arctic tundra are difficult to predict (Welker et al. 1997). The gradual expansion of shrubs along with the effects of warmer winters on the duration and magnitude of winter snowpack may add to the longer-term complexity of carbon source-sink relations in this arctic ecosystem Natali, Watts, and Rogers et al. 2019). Nonetheless, it is likely that the rapid climate change now being observed in the Arctic will continue, and the consequences of changing winter snow depth accompanying climate change will likely be one of the most important factors in the future evolution of the arctic ecosystem (Welker, Fahnestock, and Jones 2000;Natali et al. 2019). The approach demonstrated here may provide a new tool for better understanding of the response of the arctic tundra's soil organic carbon inventory to the changing climate.