Soil carbon stocks in forest-tundra ecotones along a 500 km latitudinal gradient in northern Norway

As shrubs and trees are advancing into tundra ecosystems due to climate warming, litter input and microclimatic conditions affecting litter decomposition are likely to change. To assess how the upward shift of high-latitude treeline ecotones might affect soil organic carbon stocks (SOC), we sampled SOC stocks in the surface layers of 14 mountain birch forest-tundra ecotones along a 500 km latitudinal transect in northern Norway. Our objectives were to examine: (1) how SOC stocks differ between forest and tundra soils, and (2) the relative role of topography, vegetation and climate in explaining variability in SOC stock sizes. Overall, forest soils had higher SOC stocks (median: 2.01 kg m−2) than tundra soils (median: 1.33 kg m−2). However, SOC storage varied greatly within and between study sites. Two study sites had higher SOC stocks in the tundra than in the nearby forest, five sites had higher SOC stocks in the forest, and seven sites did not show differences in SOC stocks between forest and tundra soils. Thus, our results suggest that an upwards forest expansion does not necessarily lead to a change in SOC storage at all sites. Further, a partial least-squares regression (PLSR) model indicated that elevation, temperature, and slope may be promising indicators for SOC stock size at high-latitude treelines. Precipitation and vegetation were in comparison only of minor importance.

that commonly limit plant productivity in boreal ecosystems. Enhanced plant productivity could in turn increase C sequestration, and potentially compensate for C losses via respiration 9,22,23 .
To be able to predict how treeline SOC dynamics may respond to climate change, it is crucial to quantify the size of current SOC stocks and identify the key drivers behind their accumulation. Only very few studies have examined SOC dynamics above the treeline, and these suggest that organic surface soil horizons in the tundra contain larger SOC pools compared to those of the adjacent forest 12,[24][25][26][27] . Conflicting results have been reported for deeper situated mineral soil horizons [24][25][26][27] . The high spatial variability that characterizes treeline ecotones also suggests that SOC stocks may vary substantially across very short distances. At a local scale, for example, topographic features may affect the spatial distribution of SOC stocks by controlling the flow of water and organic matter depositions. Lower slopes and concave curvatures tend to have higher SOC pools because they receive surface runoff and sediments from the surroundings. Convex curvatures and steep slopes, on the other hand, tend to lose organic matter as the topsoil is constantly eroded 28 . Slope aspect, defined as the direction a slope faces, might affect the received solar radiation intensity and create a microclimate that differs considerably from regional climatic conditions 29 . Vegetation community structure and composition is likely to be another key driver behind the accumulation of SOC through differences in litter production and its decomposability 4,30,31 , and SOC stocks are known to vary under the different treeline-forming species found in the Scandes mountains 32 .
Treeline ecotones in the northern Scandes mountains are characterized by sharp borders between mountain birch forests and tundra areas 5 , which in turn encompass marked changes in vegetation, microclimate, and topography across fine spatial scales. Here we have quantified SOC stocks in 14 forest-tundra ecotones along a 500 km latitudinal gradient in northern Norway, and we hereby present results from the currently largest field study on SOC storage in forest-tundra ecotones. The wide latitudinal extent enables examinations of possible regional differences in SOC as a result of precipitation-and temperature variability, and local variations in vegetation and microtopography. Our objectives were to examine: (1) how SOC stocks differ between forest and tundra soils, and (2) the relative role of topography, vegetation, and climate in explaining variability in SOC stock sizes. We focused on SOC in the surface soil (O and A horizons), as these horizons are rich in organic matter and have a high potential for rapid changes in C sequestration under climate change. In line with previous studies that quantified SOC stocks in the organic soil horizons at high-latitude treelines, we expected SOC stocks to be higher in tundra soils above the treeline than in adjacent forest soils.

Materials and methods
Study sites. This study was carried out along a 500 km long latitudinal transect across Norway that was part of a 1100 km long latitudinal transect originally established by Thieme et al. 33 . The study sites (~ 50 m × 200 m) along the transect are forest-tundra ecotones selected for their proximity to sample plots used in Norway's National Forest Inventory. For this study, we revisited the 14 northernmost sites (Fig. 1, Table 1). At all 14 sites, mountain birch (Betula pubescens spp. czerepanovii) is the dominating tree species. The field-and ground layer vegetation vary considerably within and between sites due to highly variable topography and climatic conditions, although it typically consists of a mix of alpine and boreal forest species, among which Betula nana, Salix spp., Vaccinium spp., Empetrum nigrum, Rubus chamaemorus, Calluna vulgaris, Carex spp., and Eriophorum spp. are common in the field layer, and Hylocomium splendens, Pleurozium schreberi, Sphagnum spp., Cladonia spp., and Flavocetraria spp. are common in the ground layer. Detailed information about the vegetation at our study sites is given by Mienna et al. 34 .
Soil sampling and analysis. At each of the 14 sites, we set up two soil sample transects in the forest and one soil sample transect in the tundra. All sample transects had a north-south orientation and consisted of 10 sample points. In the forest, we established each of the two sample transects centered around a mountain birch tree. The two mountain birch trees were located at least 10 m apart. Sample points were selected 0.5, 1.2, 1.6, 2.4 and 4 m northward and southward from the tree trunks. In the tundra, we set up one similar sample transect without a center point tree, but rather centered around a point located on the western border of each study site, 10 m southward from the top. The distance between the forest and tundra sample transects was ~ 200 m.
Soil core samples including the entire organic surface soil (i.e., O and A horizons) were collected at each sample point using a cylindrical soil sampler (Ø = 6.35 cm) during August and September 2020. As there was typically no clear border between the O-and A horizons, we did not attempt to separate these soil layers. Therefore, the two horizons were bulked and only the thickness of the entire OA complex was measured. The samples were frozen at − 20 °C until further analysis. To get precise coordinates of the soil sample points, we used a Topcon HiPer SR receiver in real-time kinetic mode, receiving differential corrections of both the Global Navigation Satellite System and the Global Positioning System.
In the laboratory, the samples were defrosted and fresh plant litter such as leaves and twigs that had not undergone observable decomposition, was removed. The soil samples were dried to constant mass at 40 °C in a drying cabinet. Bulk densities were determined based on the dry matter mass and sample volume. The dried samples were milled and homogenized. Total SOC concentrations were determined by dry combustion using a vario MICRO cube element analyser (Elementar, Hanau, Germany). Volume-based SOC stocks were calculated as: Environmental variables. Data on topographic, vegetation and climatic variables that potentially control the storage of SOC were collected in the field or extracted from existing databases. Topographic attributes (eleva- www.nature.com/scientificreports/ tion, slope, and aspect) were derived from digital elevation models with a resolution of 1 m provided by the Norwegian Mapping Authority using QGIS version 3.16.8. The aspect values were cosine and sine transformed to obtain continuous variables that represent north-south and east-west orientation. To avoid 0 values, we added 2 to both variables, thus obtaining values ranging between 1 (south) and 3 (north) for cos(aspect) and values between 1 (west) and 3 (east) for sin(aspect). In the field, each sample point was assigned a vegetation class (vascular plants, Sphagnum mosses, other mosses, or lichens) based on the dominant species within a surrounding circular plot (Ø = 60 cm) around the sample point. Climatic variables were calculated on the study site-level based on daily weather estimates interpolated from Norway's official weather stations by the Norwegian Meteorological Institute 35 . We calculated mean daily temperature, precipitation, and solar radiation for the period between 1970 and 2020. We also calculated mean summer temperature, precipitation, and solar radiation as the mean of the daily observations within the period from June 1st to September 30th. In addition, we calculated the length of the thermal growing season. Since 5 °C is a commonly accepted thermal threshold for vegetation growth in the Nordic region 36 , we defined the thermal growing season as the period with a mean daily temperature equal to or above 5 °C.

Statistical analyses.
Statistical analyses were carried out using R (R Core Team, 2021). First, we analyzed how SOC stocks differ between forest and tundra soils. All continuous variables were tested for normality using Shapiro-Wilk tests. Due to the non-normal distributed nature of our dataset, non-parametric Mann-Whitney U tests were used to assess differences in soil properties (soil depth, SOC concentration, and SOC stock) between forest and tundra soils. Differences in SOC stocks between vegetation classes were analyzed with a Kruskal-Wallis test followed by a Dunn test with Bonferroni correction of p-values. www.nature.com/scientificreports/ To evaluate the importance of environmental variables as determinants of SOC stocks, we started by calculating Spearman's rank correlation coefficients among all continuous predictor variables and SOC stocks. Due to strong intercorrelations between several variables, they could not be considered as independent determining factors. To cope with this multicollinearity problem, we used partial-least squares regression (PLSR). PLSR is a multivariate analysis method used to explain variation in one or several response variables by constructing linear combinations of predictor variables. Hence, PLSR is a particularly suitable method to analyze datasets with numerous, highly correlated predictor variables. A detailed mathematical description of PLSR can be found in Abdi 37 and Wold et al. 38 .
We built a PLSR model using the R package pls 39 . Before setting up the model, all continuous variables were standardized (mean = 0, SD = 1) to give them the same prior importance. The dependent variable (SOC stock) was log-transformed to get a normal distribution. To avoid overfitting, the number of components included in the final model was determined systematically using the one-sigma heuristic method with leave-one-out (LOO) cross validation. LOO cross validation consists of excluding each sample observation once, constructing a model without this observation, and predicting the value of the dependent variable. The residuals (i.e., the difference between observed and predicted values) are used to calculate the cross-validated root mean square error of cross validation (RMSEP CV ). The one-sigma heuristic method selects the smallest number of components still capable of producing a RMSEP CV that falls within one standard error of the absolute minimum RMSEP CV 40 . Two important statistics that describe the performance of the final model are the explained variation in the response (R 2 ) and the predicted variation in the response (Q 2 , i.e., goodness of prediction, or cross-validated R 2 ). A large difference between R 2 and Q 2 indicates that the model suffers from overfitting. In PLSR modelling, the relative influence of each predictor variable can be given by the variable importance in projection (VIP), which is the sum of the variable's influence over all model dimensions divided by the total variation explained by the model. Variables with VIP scores ≥ 1 are generally considered to be the most relevant for explaining variability in the dependent variable. The regression coefficients (RC) reveal the direction of the influence of each variable 41,42 . Ethical statement. Experimental research and field studies on plants (either cultivated or wild), including the collection of plant material, complies with relevant institutional, national, and international guidelines and legislation.

Results
SOC stocks in high-latitude treeline surface soils. Soil characteristics varied substantially across the treeline ecotone. Overall, the sampled soil horizons were significantly thicker (p < 0.001) under mountain birch canopies in the forest (median: 3.50 cm) than in the tundra (median: 1.50 cm) (Fig. 2a,b). Forest soils also had greater (p < 0.001) SOC concentrations (median: 44.0%) than tundra soils (median: 41.3%). Forest SOC stocks ranged from 0.04 up to 26.8 kg m −2 with a mean of 3.65 kg m −2 . The median value was only 2.01 kg m −2 , indicating a highly skewed distribution of forest SOC stocks (Fig. 2c). Tundra SOC stocks were significantly smaller than forest SOC stocks (p < 0.001) and covered a narrower range (0.05-14.4 kg m −2 ). Tundra SOC stocks also showed a right-skewed distribution, with most soils being characterized by small SOC stocks (mean: 2.32 kg m −2 , median: 1.33 m −2 ) (Fig. 2d). Moreover, soil characteristics varied greatly within and between study sites. An overview of soil depths and SOC concentrations of forest and tundra soils at each of the study sites can be found as Supplementary Table S1. Out of the 14 study sites, five had greater SOC stocks in the forest than in the nearby www.nature.com/scientificreports/ tundra (p < 0.01). Two sites had greater SOC stocks in the tundra than in the forest (p < 0.05) (Fig. 3)

Topographic, climatic and vegetation controls of treeline SOC stocks.
To identify which variables control the accumulation of SOC stocks in treeline soils, we computed a correlation matrix including different topographic and climatic variables ( Table 2). SOC stocks showed highly significant correlations (p < 0.01) with elevation, mean temperature, mean summer temperature, slope, length of the thermal growing season, and east-west slope aspect. However, many of these variables were inter-correlated and could therefore not be treated as independent determining factors of SOC variability. Thus, we built a PLSR model to extract and rank the most important factors controlling SOC stocks. Our modelling approach resulted in a final PLSR model with three components ( Table 3). The model had a moderate explanatory power accounting for 38.6% of the variability in SOC stocks. The corresponding crossvalidated R 2 cum (Q 2 cum) of 34.4% indicates that the model does not suffer from overfitting. Figure 4 illustrates the VIP values and regression coefficients. Elevation was the most important determinant of SOC stock variability, followed by north-south slope aspect, mean length of the thermal growing season, mean summer temperature, mean solar radiation, slope, mean temperature, and mean summer solar radiation. All other variables had a VIP value < 1 and were thereby of minor importance in explaining variability in SOC stocks. As indicated by the regression coefficients, SOC stocks decreased with increasing elevation and slope, while they www.nature.com/scientificreports/ increased with increasing temperature and solar radiation. SOC stocks were higher on north-oriented slopes than on south-oriented slopes.

Discussion
SOC stock sizes across the treeline ecotone. As boreal forests are currently advancing into tundra regions at many locations 18 , it is important to understand how SOC dynamics differ between forest and tundra soils. In this study, we sampled 14 treeline ecotones along a 500 km latitudinal gradient in northern Norway. Contrary to our expectations, our results indicate that SOC stocks in the surface soil (O and A horizons) are generally greater in mountain birch forests than in the tundra above the treeline. However, this general pattern does not hold true for all the study sites. Two study sites had higher SOC stocks in the tundra than in the nearby forest, while seven sites did not show differences in SOC storage between forest and tundra soils. At present, there is no consensus on how the vegetation transition across the treeline ecotone affects SOC storage. Overall, there are three key predictions 11,43 that may explain why we did not observe the same pattern at all study sites. First, carbon cycle models predict that the decrease in aboveground biomass and litter fall along the transition from forest into generally low-productive tundra results in smaller tundra SOC pools 44 . Second, previous field studies in the Scandes mountains suggest that SOC stocks in the organic horizon are smaller observations. Study sites are arranged from southernmost site (i.e., Humpen) to northernmost site (i.e., Suohpavuopmi). Study sites with statistically significant (Mann-Whitney U tests, p < 0.05) differences between SOC stocks in forest and tundra soils are indicated with an asterisk. Table 2. Correlation matrix of SOC stocks with climatic and topographic variables. SOC soil organic carbon stock (kg m −2 ), MT mean temperature (°C), MST mean summer temperature (°C), MP mean precipitation (mm), MSP mean summer precipitation (mm), MR mean solar radiation (Wm 2 s −1 ), MSR mean summer solar radiation (Wm 2 s −1 ), GS mean length of the thermal growing season (days). *Level of significance: p < 0.05. **Level of significance: p < 0.01. www.nature.com/scientificreports/ beneath mountain birch forest soils than beneath tundra soils due to higher soil respiration rates 12,25-27 . One commonly proposed explanation is that shrubs and trees trap more snow, which insulates the soil and thereby increases microbial activity 45,46 . Parker et al. 31 , however, found the effect of snow cover to be small, and suggested that variation in decomposition rates can rather be explained by litter quality and activity of the soil microbial community during summer. Third, forest and tundra soils may have similar SOC stocks if higher litter inputs in the forest are counterbalanced by accellerated decomposition rates. This hypothesis has been suggested for a forest-tundra ecotone in the Ural mountains, where total SOC stocks did not vary across the treeline despite differences in litter decomposibility 24 . Similar results have been found in a long-term warming experiment, where tussock tundra changing into shrub tundra increased the decomposer activity in the mineral soil horizon without altering SOC pools 47 . Our large-scale dataset demonstrates that SOC stocks do not always decrease from the forest to the tundra as has previously been reported for surface soils in mountain birch forest-tundra ecotones 6,[25][26][27] . The contrasting results between different study sites strongly suggest that changes in SOC stock sizes across the treeline depend on the extent to which local environmental conditions alter the balance between litter input and soil respiration. It should be noted that we reported SOC stocks for the entire OA complex, while earlier work gave estimates for individual soil layers 12,[25][26][27] . Thus, direct comparisons of SOC stock sizes with previous studies should be done with caution because of methodological differences. To get a full picture of SOC dynamics at the forest-tundra ecotone, soil layers should ideally be sampled down to the bedrock. However, as the organic soil horizon contains the most labile and temperature-sensitive SOC 7 , our study is a valuable contribution to the scarce data on SOC stock sizes across high-latitude treeline ecotones.  www.nature.com/scientificreports/ Within-and between site variability in SOC stock sizes. Since the quantification of SOC stock sizes through field studies at remote mountainous areas is cost-and time-intensive, the identification of variables controlling SOC storage and their incorporation into prediction models could be a promising alternative. Using PLSR, we ranked the relative importance of a few variables in explaining the high within-and between site variability in SOC stocks. Note that we included variables measured at two spatial scales. Climatic variables were measured at the study site level, while information on topography and vegetation was derived for each of the individual soil sample points. Among the considered variables, elevation was the most important determinant of SOC stock size. In previous studies, elevation has often be found to be a good predictor of SOC as it integrates the effects of temperature and precipitation, and reflects erosional and depositional processes 48 . Our model, however, did not assign substantial explanatory power to precipitation (i.e., mean precipitation and mean summer precipitation). Temperaturerelated variables (i.e., south-north slope aspect, length of the thermal growing season, mean summer temperature, mean solar radiation, mean temperature, mean summer solar radiation), on the other hand, were considered to be of major importance in explaining SOC stock size variability. Overall, high elevations were negatively correlated with temperature and SOC stocks. Although low temperatures reduce microbial activity, and consequently limit decomposition, the increase of freezing episodes and duration of snow cover also limit plant productivity 49 . Thus, the low SOC stocks found at high elevations suggest that plant productivity may be the main driver of SOC accumulation in treeline soils. Moreover, the decrease in SOC with elevation may be associated with soil erosion along the hillslope of each study site. Within each site, sample points at high elevations were typically located at shoulderslope and backslope positions, while sample points at lower elevations were situated at footslope positions. Lower slope positions receive water runoff and sediments from upper slope positions, and may thus accumulate more SOC 28,29 . Easily derived from digital elevation models and meteorological datasets, elevation, temperature and slope seem to be promising indicators for SOC stock size at treeline ecotones.
Treeline ecotones form a mosaic of vegetation communities associated with nutrient availability, snow depth, wind-, and drainage conditions 50 . In this study, we sampled a wide variety of vegetation types, ranging from the edge of poorly-drained grass bogs to extremely dry and wind-exposed lichen patches. We therefore considered vegetation as a potential driver of SOC stock size variability. Although our PLSR model suggested that vegetation was only of minor importance in predicting SOC stocks compared with climate and topography, a Kruskal-Wallis test did reveal differences in SOC storage among vegetation types. Sample points dominated by Sphagnum spp. had particularly high SOC stocks, which may be explained by their high water retention capacity, low litter quality, and anoxic soil environments 51 .
Implications of a shifting treeline for SOC stocks. Patterns in SOC stocks along the current vegetation transition from the forest to the tundra may hold clues regarding the consequences of future upwards treeline shifts 25 . Assuming that our comparison between adjacent forest and tundra soils represents a plausible 'spacefor-time' approach, our data suggest that changes in SOC storage in response to treeline shifts will be highly spatially variable due to local site characteristics. The future storage of SOC in treeline soils, however, is impacted by a large number of interacting biotic and abiotic factors, several of which we did not take into account. For instance, SOC stocks will be impacted by direct effects of warming on decomposition 12 . Also, the response of SOC dynamics to treeline shifts is likely to change as the forest matures 26 . However, since we sampled 14 study sites, of which half had similar SOC stock sizes in forest and tundra soils, our study provides strong evidence that a shifting treeline does not necessarily lead to a change in SOC storage at all sites.

Conclusion and outlook
Climate change is driving upward shifts of boreal forests into currentlytreeless tundra, and there is considerable concern that this vegetation shift may provide a positive feedback (i.e. further strengthen warming) on the global climate by reducing carbon storage. In our study covering 14 forest-tundra ecotones, only two sites had higher SOC stocks in the tundra than in the forest. Thus, we predict that the upwards advance of mountain birch forest into tundra communities does not necessarily lead to a decrease in SOC storage at all sites. Further, our results suggest that elevation, temperature, and slope may be promising indicators for SOC stock size. We strongly recommend future studies to continue identifying remotely measurable biotic and abiotic indicators of SOC stock size in forest-tundra ecotones, and we are convinced that remote sensing techniques represent highly promising alternatives to the regular monitoring of SOC stocks through time consuming and resource demanding field-and laboratory work.

Data availability
The dataset generated during the current study is available from the corresponding author on reasonable request.