Environmental and taxonomic controls of carbon and oxygen stable isotope composition in Sphagnum across broad climatic and geographic ranges

Rain-fed peatlands are dominated by peat mosses (Sphagnum sp.), which for their growth depend on nutrients, water and CO2 uptake from the atmosphere. As the isotopic composition of carbon (12,13C) and oxygen (16,18O) of these Sphagnum mosses are affected by environmental conditions, Sphagnum tissue accumulated in peat constitutes a potential long-term archive that can be used for climate reconstruction. However, there is inadequate understanding of how isotope values are influenced by environmental conditions, which restricts their current use as environmental and palaeoenvironmental indicators. Here we tested (i) to what extent C and O isotopic variation in living tissue of Sphagnum is speciesspecific and associated with local hydrological gradients, climatic gradients (evapotranspiration, temperature, precipitation) and elevation; (ii) whether the C isotopic signature can be a proxy for net primary productivity (NPP) of Sphagnum; and (iii) to what extent Sphagnum tissue δ18O tracks the δ18O isotope signature of precipitation. In total, we analysed 337 samples from 93 sites across North America and Eurasia using two important peat-forming Sphagnum species (S. magellanicum, S. fuscum) common to the Holarctic realm. There were differences in δ13C values between species. For S. magellanicum δ13C decreased with increasing height above the water table (HWT, R2 = 17 %) and was positively correlated to productivity (R2 = 7 %). Together these two variables explained 46 % of the between-site variation in δ13C values. For S. fuscum, productivity was the only significant predictor of δ13C but had low explanatory power (total R2 = 6 %). For δ18O values, approximately 90 % of the variation was found between sites. Globally modelled annual δ18O values in precipitation explained 69 % of the between-site variation in tissue δ18O. S. magellanicum showed lower δ18O enrichment than S. fuscum (−0.83 ‰ lower). Elevation and climatic variables were weak predictors of tissue δ18O values after controlling for δ18O values of the precipitation. To summarize, our study provides evidence for (a) good predictability of tissue δ18O values from modelled annual δ18O values in precipitation, and (b) the possibility of relating tissue δ13C values Biogeosciences, 15, 5189–5202, 2018 www.biogeosciences.net/15/5189/2018/ G. Granath et al.: Environmental and taxonomic controls of carbon and oxygen stable isotope composition 5191 to HWT and NPP, but this appears to be species-dependent. These results suggest that isotope composition can be used on a large scale for climatic reconstructions but that such models should be species-specific.

Abstract. Rain-fed peatlands are dominated by peat mosses (Sphagnum sp.), which for their growth depend on nutrients, water and CO 2 uptake from the atmosphere. As the isotopic composition of carbon ( 12,13 C) and oxygen ( 16,18 O) of these Sphagnum mosses are affected by environmental conditions, Sphagnum tissue accumulated in peat constitutes a potential long-term archive that can be used for climate reconstruction. However, there is inadequate understanding of how isotope values are influenced by environmental conditions, which restricts their current use as environmental and palaeoenvironmental indicators. Here we tested (i) to what extent C and O isotopic variation in living tissue of Sphagnum is speciesspecific and associated with local hydrological gradients, climatic gradients (evapotranspiration, temperature, precipitation) and elevation; (ii) whether the C isotopic signature can be a proxy for net primary productivity (NPP) of Sphagnum; and (iii) to what extent Sphagnum tissue δ 18 O tracks the δ 18 O isotope signature of precipitation. In total, we analysed 337 samples from 93 sites across North America and Eurasia us-ing two important peat-forming Sphagnum species (S. magellanicum, S. fuscum) common to the Holarctic realm. There were differences in δ 13 C values between species. For S. magellanicum δ 13 C decreased with increasing height above the water table (HWT, R 2 = 17 %) and was positively correlated to productivity (R 2 = 7 %). Together these two variables explained 46 % of the between-site variation in δ 13 C values. For S. fuscum, productivity was the only significant predictor of δ 13 C but had low explanatory power (total R 2 = 6 %). For δ 18 O values, approximately 90 % of the variation was found between sites. Globally modelled annual δ 18 O values in precipitation explained 69 % of the between-site variation in tissue δ 18 O. S. magellanicum showed lower δ 18 O enrichment than S. fuscum (−0.83 ‰ lower). Elevation and climatic variables were weak predictors of tissue δ 18 O values after controlling for δ 18 O values of the precipitation. To summarize, our study provides evidence for (a) good predictability of tissue δ 18 O values from modelled annual δ 18 O values in precipitation, and (b) the possibility of relating tissue δ 13 C values to HWT and NPP, but this appears to be species-dependent. These results suggest that isotope composition can be used on a large scale for climatic reconstructions but that such models should be species-specific.

Introduction
Peatlands in temperate, boreal and arctic regions form large reservoirs of carbon, which are vulnerable to release under expected changes in global climate and land management (Rydin and Jeglum, 2013;Loisel et al., 2014). Because peat decomposes slowly and gradually accumulates, it preserves information on past peatland ecosystem dynamics and responses to allogenic and autogenic forcings. Palaeoenvironmental studies of peat may, therefore, help to anticipate the future responses of these globally important ecosystems to climate change (Loader et al., 2016). Past climate and local hydrology can be estimated using a variety of biotic and biogeochemical proxies, including the δ 13 C and δ 18 O values of organic material (e.g. van der Knaap, 2011;Royles et al., 2016;Lamentowicz et al., 2015). However, the environmental (e.g. climate) and biotic (e.g. species identity) controls of isotope differentiation in peatland-dwelling plants are still poorly understood and current assumptions regarding these controlling factors are yet to be tested on larger spatial scales.
Sphagnum mosses are the most dominant peat-forming plant group in acidic peatlands. The composition of stable isotopes of carbon and oxygen in their tissues is affected by different environmental conditions, operating through their impact on fractionation processes. When not submerged, carbon isotope signals in bulk tissues or components such as cellulose depend mainly on the concentration and isotopic composition of CO 2 in the chloroplast, which alters isotope discrimination during biochemical fixation of CO 2 (Farquhar et al., 1989;O'Leary, 1988). In mosses, the CO 2 concentration in the chloroplast, [CO 2 ] c , is determined by temperature, light availability, CO 2 partial pressure and, most importantly, plant water status (Finsinger et al., 2013;van der Knaap et al., 2011;Ménot and Burns, 2001;Ménot-Combes et al., 2004;Royles et al., 2014;Skrzypek et al., 2007a;Kaislahti Tillman et al., 2013). When wet, external water films on leaf surfaces impede diffusion and [CO 2 ] c is lowered (Rice and Giles, 1996;Rice, 2000;Williams and Flanagan, 1996); consequently, the proportion of fixed 13 C increases due to internal drawdown of the preferred isotope 12 C. When submerged, assimilation of respired or methane-derived CO 2 can alter [CO 2 ] and also the isotopic composition of C in Sphagnum (Raghoebarsing et al., 2005). Even when not submerged, respiratory carbon can be refixed by Sphagnum (Turetsky and Wieder, 1999;Limpens et al., 2008). Given that respired CO 2 is isotopically lighter than that in the atmosphere, it may also contribute to variation in tissue isotope values. Despite many detailed studies, there remains uncertainty about how the multiple controls on 13 C isotope values combine to determine isotopic composition, and how universal the proposed mechanisms are on a global scale. This uncertainty currently restricts the utility of C isotope signals as palaeoclimatic or palaeoenvironmental indicators in peatlands (Loader et al., 2016).
Oxygen isotope values in moss tissues depend on the isotopic composition of the water sources, enrichment associated with evaporation from the moss surface and biochemical fractionation (Dawson et al., 2002). Once on the plant, 18 O present in water equilibrates with that in atmospheric CO 2 prior to fixation as well as being incorporated directly during hydrolysis reactions, especially during the initial stages of carbon fixation (Gessler et al., 2014;Sternberg et al., 2006). Hence, variation in tissue oxygen isotopes reflect environmental conditions that control source water (rainfall, snowfall, groundwater) as well as fractionation caused by evaporation prior to fixation, which is controlled by micrometeorological conditions (mainly temperature, relative humidity and incident energy) (Daley et al., 2010;Moschen et al., 2009;Royles et al., 2013;Kaislahti Tillman et al., 2010). Oxygen isotope composition has, therefore, been used to reconstruct climatic conditions and to infer the dominant water source in peatlands (Aravena and Warner, 1992;Ellis and Rochefort, 2006;van der Knaap et al., 2011). Ongoing measurements of oxygen isotopes in precipitation across the globe (Bowen, 2010;IAEA/WMO, 2015) have generated models that predict spatial patterns in oxygen isotope composition of precipitation based on temperature, elevation, atmospheric residence time and circulation patterns (e.g. Bowen, 2010). Once isotopic composition of the source water is accounted for, variation in moss tissue isotopic values should be largely determined by fractionation that accompanies evaporation from the surface of plants. How well oxygen isotopes in Sphagnum tissues reflect atmospheric water or plant surface water depends on local weather conditions such as precipitation, air temperature and humidity. For example, Bilali et al. (2013) suggest that oxygen isotopes in Sphagnum mosses from maritime bogs will track variation in precipitation patterns whereas isotopic values in continental habitats will be more dependent on summer temperature, as temperature and humidity are more variable in those regions. On local scales, oxygen isotope values also vary as a function of temperature and humidity. Aravena and Warner (1992) found differences that correspond with changes in microtopography. Elevated microsites (hummocks) were enriched in 18 O, which they ascribed to higher evaporation compared to that of neighbouring wet depressions (hollows). However, as with 13 C, there remains uncertainty in how 18 O signatures relate to environmental factors and species identity and to what extent global 18 O patterns in precipitation dominate over local processes.
Stable isotopes can also serve as indicators of net primary productivity (NPP) (Rice and Giles, 1996;Williams and Flanagan, 1996;Rice, 2000). However, few studies have explored these relationships in the field. In a multispecies 5192 G. Granath et al.: Environmental and taxonomic controls of carbon and oxygen stable isotope composition comparison of peat mosses, Rice (2000) found that plants with higher relative growth rates had lower discrimination against 13 C and therefore were more enriched in 13 C. This was attributed to the local environment, with fast-growing plants of wetter microhabitats having thicker water films that inhibit CO 2 diffusion into the plant and to species-specific differences in maximum rates of photosynthesis. Both factors would reduce internal [CO 2 ] and thereby lower discrimination. In line with this, a warming experiment by Deane-Coe et al. (2015) reported a positive relationship between moss net NPP and δ 13 C values for tundra mosses (Dicranum, Pleurozium, Sphagnum). Clearly, carbon isotope values show promise as indicators of peat moss contemporary growth and potentially as an NPP proxy in palaeoecological studies. This could be particularly valuable to differentiate productivity and decomposition controls in long-term carbon accumulation studies. To date, we are not aware of attempts to explore the robustness of these relationships across large spatial scales.
Together, tissue carbon and oxygen isotope compositions are controlled both by environmental factors at micro-and macro-scales, and by species-specific differences that relate to water balance and carbon dynamics in peat mosses. Palaeoecological studies rely on such environment-isotope relationships for environmental reconstructions (Ellis and Rochefort, 2006;van der Knaap et al., 2011). The underlying mechanisms are, however, rarely fully explored using known environmental gradients (but see Ménot and Burns, 2001 for an example) or only tested across narrow bands of environmental variation, often with sets of correlated environmental factors (Loader et al., 2016). Moreover, interactions with biotic factors such as species identity have received little attention despite the large variations in Sphagnum species dominance commonly observed down peat cores (e.g. Ménot and Burns, 2001). Here we aim to provide a robust, cross-scale evaluation of how environmental factors and species identity influence the C and O isotope compositions of Sphagnum using two common and widely distributed peat-forming species (S. magellanicum and S. fuscum) that are primarily rain-fed. To achieve this, we performed an unprecedented large sampling campaign across the Holarctic realm.
Specifically, we (i) investigated relationships between C and O isotope values and factors known to influence plant water availability (height above the water table -HWT, temperature, evaporation and precipitation) and CO 2 partial pressure (elevation), and tested whether their effects were modified by species identity; (ii) tested the prediction that Sphagnum tissue δ 13 C values are associated with NPP; and (iii) tested whether tissue δ 18 O in rain-fed Sphagna is predicted by the δ 18 O isotope signature in precipitation but modified by negative relationships with precipitation and positive ones with temperature/evaporation. Across these objectives we examined how C and O isotope values varied with scale (within-peatland versus between-peatlands) and to what ex-tent HWT and NPP could explain variation within and between peatlands.

Study species and collection sites
Our study focused on two common peat-forming Sphagnum species, S. fuscum (Schimp.) H. Klinggr. (circumpolar distribution) and S. magellanicum Brid. (cosmopolitan distribution). In general, these species are confined to primarily rain-fed peatlands (bogs) and described as hummock (S. fuscum) and lawn (S. magellanicum) species. However, S. magellanicum is a species with a very broad niche and found in a range of habitats with varying degrees of groundwater influence (Flatberg, 2013). These species are easy to identify, but recent research has shown that the dark European morph of S. fuscum is conspecific to the North American S. beothuk (Kyrkjeeide et al., 2015), and S. magellanicum has been shown to consist of two genetically diverged morphotypes (Kyrkjeeide et al., 2016) that recently were separated at the species level (Hassel et al., 2018). Unpublished genetic data suggest that samples collected in our study consist of both S. magellanicum morphs (approximately 50/50) and possibly one or two samples of S. beothuk (Narjes Yousefi, personal communication, 2018). Hence, we here treat our species as aggregates (i.e. species collectiva), S. fuscum coll. and S. magellanicum coll.
The two species were sampled across the Holarctic region at a total of 93 sites ( Fig. 1; Table S1 in the Supplement) at the end of the growing season. To make comparisons between species and between sites possible, we focused on habitats where both species can be found and have low influence of surrounding groundwater. Thus, we only sampled bogs (including a few poor fens with ombrotrophic character) and open (no tree canopy) habitats. Sampling was conducted mainly during 2013, but a few sites were sampled at a similar time of year in 2014. At each site two patches (minimum 10 m apart) for each species were sampled (except for 11 sites that contained only one patch for one species). At each sampling patch we recorded moss growth, HWT and GPS coordinates, and collected a moss sample (78 cm 2 and 5 cm deep) at the end of the growing season (September to November depending on location and generally coincided with when there was a risk of the first snowfall to occur). Moss samples were dried (24 h at 60-65 • C) within 72 h or immediately frozen and later thawed and dried. The apical part (the capitula, top 1 cm) of the dried plant shoots was used for isotope analysis, while the stem section was used for bulk density estimation to calculate moss NPP.

Isotope determination
Ten capitula from each patch were selected and finely chopped with a single-edge razor by hand and mixed. Ca- At some sites only one of the two Sphagnum species was sampled, indicated by red triangles or black circles. Otherwise sites contained both species (blue crosses). The map is centred on the North Pole and has an orthographic projection. Geographical ranges: latitude 41.6-69.1 N, elevation 2 -1829 m a.s.l. See Table S1 for details.
pitula were chosen as they reflect the most recently fixed organic matter and should relate better to recent growing season conditions. In Sphagnum, δ 13 C from the capitulum is similar to that of branches within the top 15 cm of plants but is approximately 1-2 ‰ less negative than stems (Loader et al., 2007). For δ 18 O, the offset between branches and stems is around 1 ‰ (Moschen et al., 2009). Standard deviations of repeated samples were 0.6 ‰ and 0.7 ‰ for δ 13 C and δ 18 O, respectively. Approximately 0.5 mg of dry sample was packed in tin cups for δ 13 C analyses, and ∼ 0.2 mg in silver cups for δ 18 O analyses. Samples were analysed at Union College (Schenectady, NY, USA) using a Thermo Delta Advantage mass spectrometer in continuous flow mode connected (via a ConFlo IV) to a Costech Elemental Analyzer for δ 13 C analysis or a Thermo TC/EA for δ 18 O analyses. Isotope values are presented as 1000 × (R sample /R standard − 1), where R sample and R standard are the ratios of heavy to light isotopes (e.g. 13 C / 12 C) and are referenced to VPDB and VSMOW for C and O, respectively. Carbon isotope data were corrected using sucrose (IAEA-CH-6, −10.449 ‰), acetanilide (in house, −37.07 ‰) and caffeine (IAEA-600, −27.771 ‰). Oxygen isotope data were corrected using sucrose (IAEA-CH-6, 36.4 ‰), cellulose (IAEA-C3, 31.9 ‰) and caffeine (IAEA-600, −3.5 ‰) with values from Hunsinger et al. (2010). Oxygen isotope standardization was further checked with the whole wood standards USGS54 and USGS56. The combined instrument uncertainty for δ 13 C (VPDB) is < 0.1 ‰ based on the in-house acetanilide standard and < 0.5 ‰ for δ 18 O (VSMOW) based on the cellulose standard (IAEA-C3). We performed isotope analyses on whole-plant tissue rather than on cellulose extracts. In living Sphagnum samples, there is a strong linear relationship between the isotopic composition of these two components for both δ 13 C (R 2 values 0. 89-0.96;Kaislahti Tillman et al., 2010;Ménot and Burns, 2001;Skrzypek et al., 2007b) and for δ 18 O (R 2 values 0. 53-0.69;Kaislahti Tillman et al., 2010;Jones et al., 2014). Focussing on whole-plant tissue allowed us to analyse a higher number of samples for this study, allowing larger numbers of sites and more replication.

Environmental variables
The modelled δ 18 O signal in meteoric water (precipitation) (Bowen and Wilkinson, 2002) was obtained from http:// www.waterisotopes.org (last access: 2 October 2017) as annual and monthly isotope ratio estimates at 10 arcmin resolution. These global estimates have shown to be highly accurate (R 2 = 0.76 for mean annual δ 18 O in precipitation) and are based on absolute latitude and elevation and account for regional effects on atmospheric circulation patterns (for details see Bowen, 2010Bowen, , 2017IAEA/WMO, 2015). To test which temporal period of δ 18 O values in precipitation showed the highest correlation with tissue δ 18 O values, we calculated annual (January-December), growing season (May-October) and winter-spring (January-April) mean isotope ratio. We calculated both unweighted and weighted means against precipitation for each month. Monthly precipitation (PRECTOTCORR), land evapotranspiration (EV-LAND) and surface air temperature (TLML) for each site and year of sampling (2013 or 2014) were retrieved from the NASA GESDISC data archive, land surface and flux diagnostics products (M2T1NXLND, M2TMNXFLX; resolution longitude 0.667 • , latitude 0.5 • ; Global Modeling and Assimilation Office, 2015a, b). Total precipitation and evapotranspiration (ET), and mean temperature, from April to October were used as predictors in the statistical models. As ET can be compensated for by precipitation, we used the ET/P quotient as a predictor for the effect of water loss. A high value (> 1) indicates a net loss of water to the atmosphere. Site altitude was retrieved from a global database using the R package elevatr (ver 0.1-2, Hollister and Shah, 2017).
The distance from the moss surface to the water table (HWT) was measured using water wells (commonly a PVC pipe, 2-5 cm in diameter and slotted or perforated along the sides) with a "plumper" (a cylinder on a string that makes a "plump" sound when it hits the water surface) or a "bubbler" (a narrow tube that makes bubbles when it hits the water surface while the user blows into it). HWT was measured in the spring and in the autumn and there was a strong correlation between the two time points (r = 0.74). As growth mainly Table 1. Sample sizes, standard variation and overall partitioning of measured variation for each species and response (δ 13 C and δ 18 O). N site is the number of sites and N obs the total sample size. Standard deviation of the responses is given for within and between sites, together with the proportion of total variance measured between sites and within sites. occurs in late summer to autumn in temperate and boreal regions, we used HWT at the end of the season as the proxy of relative HWT between sites.

Moss growth
Moss growth (or productivity, NPP) was measured with a modified version of the cranked wire method (see Clymo, 1970;Rydin and Jeglum, 2013 for details), with bristles from a paint brush spirally attached to a wire. These "brush wires" were inserted in the moss layer with the end of the wire protruding above the surface. Height increment (i.e. vertical growth) was measured over the growing season as the change in distance (to nearest millimetres) between the moss layer and the top of the wire. A minimum of three wires were inserted within a 1×1 m uniform area (same microhabitat, vegetation and general structure). To determine moss bulk density (kg m −3 ), we dried (24 h at 60-65 • C) the top 30 mm of the stems (area 78 cm 2 ) in our collected core (see Sect. 2.1). Biomass growth on an area basis (g m −2 yr −1 ) was calculated as height increment × bulk density.

Statistical analyses
To test and quantify the influence of environmental variables and species identity on isotope composition, we used linear mixed models in R (R core team, 2016), employing the R package lme4 ver 1.1-12 (Bates et al., 2015). Site dependence (i.e. multiple samples from the same site) was accounted for by adding site as a random factor. For tissue δ 13 C, we first fitted two separate models to test the independent effects of HWT, NPP and species identity (S. fuscum and S. magellanicum), and whether the HWT or NPP effect varied between species by fitting a species interaction term. To test the explanatory power of environmental variables (ET/P , precipitation, temperature, elevation) we first constructed a base model with HWT and NPP, as they were identified as the main predictors in literature. For simplicity we removed negligible interactions from this model. Each environmental variable and its interaction with species was then tested against the base model. For tissue δ 18 O, we first explored which temporal period of modelled δ 18 O precip (annual, growing season, winter-spring) had the highest explanatory power and whether the relationship varied between species. The identified best model was then used as the base model to separately test each environmental variable (HWT, ET/P , precipitation, temperature) and its interaction with the species. The proportion of variance explained by the predictors was calculated at the site level (Gelman and Hill, 2007) or as marginal R 2 (Nakagawa and Schielzeth, 2013; R package piecewiseSEM ver 1.1-4; Lefcheck, 2016). Although our study focused on explained variance by predictors, we also performed statistical tests of predictors and their interactions using type-2 (main effects tested after all the others in the model but without the interaction term) F tests, applying Kenward-Roger adjustments to the degrees of freedom, as implemented in the car package (ver. 2.1-3, Fox and Weisberg, 2011). Standard model checking was performed (e.g. residual analyses and distribution of random effects) to ensure compliance with model assumptions. Covariances between predictors were small (r < 0.15) or moderate (r = 0.40-0.50 between ET/P , precipitation and temperature) and this multicollinearity had minor impact on model estimates. (Table S2). Due to uncertainty in height increment measurement we recorded a few negative values resulting in negative NPP. These values were kept in the analyses. The means and standard deviations (in parenthesis) of height increment (HI, millimetres) were 14.3 (10.1) for S. fuscum and 19.5 (14.1) for S. magellanicum, and the means and standard deviations of bulk density (BD, kg m −3 ) were 17.8 (9.9) for S. fuscum and 10.2 (7.6) for S. magellanicum. Table 2. Results from linear mixed models for δ 13 C values. Statistical tests are based on type-2 F test using Kenward-Roger-adjusted degrees of freedom. The second model only included S. magellanicum. Elevation (m a.s.l.) and the three climatic variables (growing season sums and means: ET/P , P in millimetres, temperature in degrees Celcius) were tested one by one in the model including HWT (height above the water table in centimetres), species and NPP (g m −2 yr −1 ). For simplicity, the negligible HWT × NPP term was dropped from this model (P = 0.36). Estimated effects (± SEs) are only given for the main effects if interactions were considered negligible. These effects are slopes for continuous variables (all variables except species) and for species (categorical) the difference between S. magellanicum and S. fuscum (i.e. S. fuscum being the reference level). In the presence of an interaction between HWT and species, the species effect was estimated at mean HWT. R 2 site is explained between-site variance; R 2 marginal is explained total variance.

δ 13 C signal
Variation in Sphagnum tissue δ 13 C values was marginally greater within sites than between sites (Table 1). HWT predicted the δ 13 C values, but the relationship differed between the two species (Table 2, Fig. 2). Although δ 13 C values decreased with increasing HWT for both species, the slope was less steep for S. fuscum and this species had slightly higher δ 13 C values overall. In separate models for the two species, HWT for S. fuscum had near-zero explanatory power, while for S. magellanicum HWT explained 33 % of the betweensite variation and 17 % of the total variance (i.e. marginal R 2 ).
Measured δ 13 C values were related to moss NPP and δ 13 C values increased by 0.0023 ‰ (SE: 0.00048) for each milligram of biomass produced per square metres. NPP ex-plained 11 % of the between-site variation in δ 13 C and 7 % of the total variation. HWT and NPP explained 48 % of the between-site variation of δ 13 C in S. magellanicum and 24 % of the total variation. Corresponding values for S. fuscum were 6 % and 7 %. Of the additional environmental variables tested, we found weak evidence that ET/P and temperature were positively correlated with δ 13 C but only for S. magellanicum (Table 2).

δ 18 O signal
Sphagnum tissue δ 18 O values varied more between sites than within sites, and at similar magnitudes and proportions for both species (Table 1). Tissue δ 18 O values were predicted by the spatially explicit estimates of δ 18 O value isotope signature in precipitation (Fig. 3, Table 3). Annual mean δ 18 O precip Table 3. Results from linear mixed models for δ 18 O values. Statistical tests are based on type-2 F test using Kenward-Roger-adjusted degrees of freedom. Three time periods for modelled δ 18 O values (‰) in precipitation were tested individually: annual mean, growing season (April-September) and spring (January-April). The three climatic variables (growing season sums and mean: ET/P , P (mm), temp [ • C]) were tested one by one in a model including HWT (cm) and mean annual δ 18 O values). Estimated effects (± SEs) are only given for the main effects if interactions were considered negligible. These effects are slopes for continuous variables (all variables except species) and for species (categorical) the difference between S. magellanicum and S. fuscum (i.e. S. fuscum being the reference level). R 2 site is explained between-site variance; R 2 marginal is explained total variance.  (Fig. 3, Table 3).
HWT at the end of the growing season was, on average, 11 cm lower in S. magellanicum patches (wetter habitat) compared to S. fuscum (HWT = 33 cm) patches (F 1,224 = 131.9, P < 0.0001). However, we found only very weak support for the hypothesis that HWT predicts tissue δ 18 O values, as HWT explained < 1 % of the δ 18 O variation (Table 2). There was negligible influence of the additional environmental variables on δ 18 O values (Table 2). ET/P was associated with higher δ 18 O values in S. magellanicum and lower in S. fuscum (but not different from the zero effect), while increasing temperature was weakly associated with overall lower δ 18 O values.

Stable carbon isotope discrimination in Sphagnum
Our data were consistent with the hypothesis that moss growing closer to the water table (low HWT) has reduced carbon isotope fractionation, leading to greater fixation of 13 CO 2 and more 13 C-enriched tissue (Rice and Giles, 1996;Williams and Flanagan, 1996). Given that the water table position was measured in different places at different times and all are one-time measurements, this result is remarkably robust. For example refixation of 12 C-enriched substratederived CO 2 in living Sphagna (Turetsky and Wieder, 1999;Raghoebarsing et al., 2005) can potentially contribute to within-site variation in δ 13 C as it potentially affects both the ambient concentration of CO 2 as well as its isotopic composition. Interestingly, the strength of the δ 13 C -HWT relationship differed in the two species, with S. magellanicum exhibiting a greater reduction in δ 13 C in response to drier conditions (high HWT) than S. fuscum. The weaker effect of HWT on δ 13 C values in S. fuscum is likely a consequence of limited fluctuation in tissue water content, as this species is well known to store abundant water within capillary spaces  Table 1. and resist drying (Rydin, 1985), thus maintaining the waterfilm that results in reduced fractionation. Loader et al. (2016) reported a similar slope estimate for S. magellanicum in a single peatland, and several studies have confirmed the effects of contrasting microtopography (i.e. hummock -hollow differences) using multi-species comparisons (Price et al., 1997;Loisel et al., 2009;Markel et al., 2010). As such, our results suggest that species-specific differences in carbon isotope discrimination in Sphagnum are related to water retention capacity and, consequently, become more apparent under drier conditions. This supports the results of previous, smaller-scale studies (Rice, 2000). The influence of species identity on the relationship between δ 13 C values and water table position has important implications for palaeoenvironmental reconstructions based on δ 13 C values. The relationship between δ 13 C and HWT has been used in palaeoecological reconstructions of surface wetness (e.g. Loisel et al., 2009). In our data set the strength of the relationship was weaker than previously reported. For instance, Loader et al. (2016) reported R 2 = 54 % for S. magellanicum in a single site. Given the characteristics of our data (large-scale, circumpolar), the explanatory power (R 2 marginal = 17 %) can be considered acceptable and comparable to other proxies such as testate amoebae (16 % in Loader et al., 2016;Sullivan and Booth, 2011). Our results imply that isotopic signals of peatland wetness in hummockdwelling species (such as S. fuscum) may be weaker or ab-   sent compared to lawn species. It is therefore important that the same species or species type (e.g. lawn species as they likely have a broad HWT niche) is used if δ 13 C values are employed as a proxy to infer changes in HWT.
We also identified evidence that evapotranspiration (ET) and NPP modify δ 13 C values, although the effect of ET was weak and restricted to S. magellanicum. We expected a stronger relationship, as ET and temperature control δ 13 C by increasing water loss at the moss surface and reducing the diffusive resistance (i.e. reducing CO 2 limitation), enabling increased discrimination against 13 C (Williams and Flanagan, 1996). NPP only explained a small proportion of the variation in δ 13 C values, but the relationship was apparent across species. Several studies have proposed the use of δ 13 C values to infer Sphagnum productivity (e.g. Rice and Giles, 1996;Rice, 2000;Munir et al., 2017) and our study is the first to test this at the pan-Holarctic scale. Deane-Coe et al. (2015) investigated δ 13 C values across moss species (including Sphagnum) and years at one site and found a weak relationship between productivity and δ 13 C values (R 2 = 0.10 and 0.31). Similarly, Rice (2000) reported that the relative growth rate explained about 25 % of the variation in 13 C discrimination. We did not find as strong a relationship (R 2 < 0.12), but our study was geographically broader and less controlled; consequently, our results were likely influ-enced by more complex interactions among environmental factors that affect Sphagnum growth across our sites. Nevertheless, our results indicate independent effects of evaporation and productivity on δ 13 C values. The lack of a strong NPP pattern somewhat limits the ability to infer productivity of Sphagnum in palaeoecological studies.

Global patterns of δ 18 O values in Sphagnum
Modelled δ 18 O values in precipitation (Bowen, 2010) explained much of the variation in δ 18 O tissue values between sites (R 2 = 68 % for annual mean δ 18 O precip ). The percent variance explained was even higher if the spring period for modelled 18 O precip was used but was lower for the growing season average. This result does not necessarily mean that spring season water was utilized by the plants during the growing season. Between-site variation in 18 O precip values are much larger in the winter (Fig. S1), more effectively discriminating maritime and continental regions (Bowen, 2010). The better fit may simply be an effect of a more distinct separation of 18 O precip in the winter data. Although the δ 18 O tissue − 18 O precip relationship presented here is robust, a few δ 18 O values are less well predicted by the regression model and they originate from Northwest Territories (Canada) and West Siberia (Russia). Likely, this suggests that the 18 O precip model is less accurate in these interior regions with fewer precipitation collection stations.
In contrast, the data did not support a negative correlation between precipitation amount and δ 18 O tissue values and δ 18 O tissue values were only weakly affected by predictors associated with water loss (ET/P and/or temperature) and species identity. The indication of 18 O enrichment in S. magellanicum due to ET/P was expected as the lighter isotope 16 O needs less energy to vaporize. However, the opposite trend was suggested for S. fuscum, and surprisingly, higher surface temperatures decreased 18 O enrichment. Hence we conclude that climatic variables associated with water loss were weak predictors after controlling for δ 18 O precip values. This result may not be too unexpected as laboratory experiments have so far failed to relate 18 O enrichment in Sphagnum to differences in evaporation rates (Brader et al., 2010).
There have been few regional studies on moss δ 18 O tissue values that span gradients of δ 18 O precip values (Royles et al., 2016;Skrzypek et al., 2010) and most interpretations of moss δ 18 O tissue -climate relationships come from peat core studies (e.g. van der Knaap et al., 2011). In Antarctic non-Sphagnum peat banks variation in δ 18 O cellulose values tracked δ 18 O values in moss water across a latitudinal gradient (61-65 • S) despite a lack of difference in δ 18 O precip . This result led Royles et al. (2016) to suggest that moss water and tissue δ 18 O values are better temporal integrators of source water than point rainfall measurements. The authors interpreted site-to-site differences as relating to differential evaporative enrichment and other physio-chemical factors that affect 18 O exchange, fixation and biochemical synthesis. Skrzypek et al. (2010) explored variation in Sphagnum δ 18 O tissue values across a regional altitudinal gradient and found no consistent trend or significant relationship linking δ 18 O tissue values to altitude, where δ 18 O in source water is expected to differ. Although fractionation in source water caused by adiabatic cooling with altitude should lead to altitudinal effects, differences in precipitation amount can confound this pattern (Gat et al., 2000). Unfortunately, there are limited regional studies that have tested the effects of variation in source water on δ 18 O tissue values. The present study provides a much greater range of geographical and environmental variation and shows strong support for Sphagnum strongly tracking source water.
Interestingly, the relationship between δ 18 O tissue and δ 18 O precip values detected here is very similar to that proposed some time ago by Epstein et al. (1977); δ 18 O cellulose = 27.33 + 0.33 × δ 18 O precip (note that Jones et al., 2014, show high correspondence between δ 18 O cellulose and δ 18 O tissue values). However, our data suggest a slightly steeper slope and lower intercept, particularly for S. magellanicum. The species effect on δ 18 O suggests a difference in the degree of evaporation from the plant surface prior to the uptake of water. The lower δ 18 O values for S. magellanicum compared to S. fuscum (−0.83 ‰) is comparable to the results from bogs in Canada for the same species (−2.2 ‰, Aravena and Warner, 1992) and between a hollow and a hummock species in the Netherlands (−2 ‰, Brenninkmeijer et al., 1982). This suggests that the absorbed water in S. magellanicum was subject to less evaporation. In Sphagnum plants, surface water is largely affected by capillarity, water storage and reducing conductance with compact morphology. Plant traits that enhance these functions are more pronounced in species and individuals found at high HWT as these characteristics maintain high tissue water content (Hayward and Clymo, 1982;Laing et al., 2014;Waddington et al., 2015). Consequently, during droughts, Sphagnum species growing close to the water table will dry out quickly as the evaporative demand cannot be balanced, and simultaneously photosynthesis is shut down. Sphagnum species higher above the water table wick water from below and store water effectively, thereby remaining photosynthetically active while water is lost due to evaporation. This mechanism would result in 18 O enrichment being higher above the water table (Brenninkmeijer et al., 1982;Aravena and Warner, 1992) and explains the positive relationship between HWT and δ 18 O in S. magellanicum reported by Loader et al. (2016) along a 10 m transect. We found a weak positive relationship of δ 18 O with HWT, which suggests that HWT cannot entirely explain speciesspecific differences in 18 O enrichment. Instead, this can be attributed to lower water retention (i.e. higher evaporation at the same water deficit) in S. magellanicum compared to S. fuscum (Clymo, 1973;McCarter and Price, 2014). Although species differences in 18 O have been reported (Aravena and Warner, 1992;Zanazzi and Mora, 2005;Bilali et al., 2013), our study suggests that the species-specific δ 18 O signals may not simply be a consequence of growing at different HWT but can rather reflect distinct water retention capacity in these species.
The strong influence of δ 18 O precip values and, to a much lesser extent, environmental variables related to water loss, combined with a relatively small within-site variation in δ 18 O tissue values, suggest that macroclimatic drivers, such as precipitation inputs, largely determine the δ 18 O value of peatland moss tissue. These results are promising for the use of oxygen isotopes in large-scale palaeoecological reconstructions from peat cores (Ellis and Rochefort, 2006;Chambers et al., 2012;Daley et al., 2010), although a better understanding of O isotope fractionation within tissue components and their decay relationships would improve their utility. Moreover, the simple relationships presented here can potentially be utilized to trace changes in δ 18 O precip values that mirror climate variability.

Conclusions
Our study provides new insights into large-scale variation in Sphagnum tissue isotopic signature and suggests that isotopic composition can be used for climatic reconstructions. We show a close link between precipitation and tissue δ 18 O values and conclude that variation in δ 18 O values are mainly driven by the macroclimate, but species differences exist. In contrast, δ 13 C values were strongly related to local microtopography, while the influence of macroclimate was negligible. As suggested in earlier studies, δ 13 C values were also weakly associated with NPP. These conclusions were most strongly supported for the cosmopolitan S. magellanicum complex and species identity should be accounted for in future carbon isotope studies to avoid spurious conclusions.
Author contributions. SKR, GG and HR initiated the study and formulated the research objectives. All authors were involved in data collection and SKR, NB, KP, AV and DPG performed the isotope analyses. GG performed the statistical analyses and wrote the first draft with input from SKR and HR. All authors read and commented on the manuscript and approved the final version.
Competing interests. The authors declare that they have no conflict of interest.