Continental scale variability of foliar nitrogen and carbon isotopes in Populus balsamifera and their relationships with climate.

Variation across climate gradients in the isotopic composition of nitrogen (N) and carbon (C) in foliar tissues has the potential to reveal ecological processes related to N and water availability. However, it has been a challenge to separate spatial patterns related to direct effects of climate from effects that manifest indirectly through species turnover across climate gradients. Here we compare variation along environmental gradients in foliar N isotope (δ15N) and C isotopic discrimination (Δ13C) measured in 755 specimens of a single widely distributed tree species, Populus balsamifera, with variation represented in global databases of foliar isotopes. After accounting for mycorrhizal association, sample size, and climatic range, foliar δ15N in P. balsamifera was more weakly related to mean annual precipitation and foliar N concentration than when measured across species, yet exhibited a stronger negative effect of mean annual temperature. Similarly, the effect of precipitation and elevation on Δ13C were stronger in a global data base of foliar Δ13C samples than observed in P. balsamifera. These results suggest that processes influencing foliar δ15N and Δ13C in P. balsamifera are partially normalized across its climatic range by the habitat it occupies or by the physiology of the species itself.

Introduction

Spatial variation in foliar chemical traits represent phenotypic responses to environmental and genotypic factors. Because foliar chemicals regulate key biological functions, such as photosynthesis, spatial variation in foliar chemicals can also be used to understand the functioning of landscapes, including productivity and water relations of plants[1, 2]. Global gradients in foliar chemical traits, therefore, have the potential to reveal ecological controls on plant functioning. However, both intra and inter-specific variation contribute to relationships between foliar chemical traits and clin class="Disease">mate gradients[3], and these relationships are not necessarily congruent[4]. Analysis of how individual species respond to climate gradients could provide greater insight into ecological controls on plant functioning, possibly aiding the understanding of plant response to global change.

Of particular interest to understanding plant-environment relationships are the isotopic composition of nitrogen and n class="Chemical">carbon in leaves. Foliar δ15N, which serves as an indicator of terrestrial N cycling[5, 6], varies systematically between species[7-9], mycorrhizal associations[10, 11] and along gradients in climate[12], and, when measured over time, can be used to infer ecosystem response to disturbance and climate change[13-16]. Local to regional-scale anthropogenic impacts on N cycling can also be important influences on foliar δ15N[17-20]. Mechanistically, processes that fractionate the stable isotopes of nitrogen occur to a greater extent when N availability (i.e., supply relative to demand) is high[18, 21]. Thus, over time, processes associated primarily with gaseous N loss typically increase the δ15N value of the remaining inorganic N pool, causing 15N-enrichment in foliar tissues that draw from that pool. Conversely, systems experiencing low N availability have a more conservative N cycle and do not lose as much N with low δ15N values. Furthermore, when N availability is low, plants depend more strongly on N fixation, atmospheric deposition, and N from mycorrhizal fungi[6]. While the δ15N of N passed to plants through these processes differs, it is generally depleted in 15N relative to the enriched pool remaining following leaching and denitrification losses[22]. Therefore, continental- to global-scale patterns in foliar δ15N related to climate reflect climatic control on processes that influence N availability[6, 23], but these patterns must account for mycorrhizal association, and nevertheless, leave considerable small-grain variability related to local processes[24].

Carbon isotope discrimination (Δn class="Chemical">13C) represents an integration of the balance of CO2 and water fluxes in plants, and variation in Δ13C is controlled by both environmental and genetic factors[22, 25]. For example, there is a well-documented positive relationship between plant-moisture status and Δ13C across a variety of C3 plants e.g. refs 26–28. This pattern emerges because long-term intrinsic water-use efficiency (ratio of photosynthesis to stomatal conductance) and Δ13C are both influenced by Ci/Ca, the ratio of intercellular to ambient [CO2][29]. Dry conditions lead to reduced stomatal conductance (and Ci, the reservoir of CO2 available for photosynthesis) and a proportionately greater decrease in transpiration than photosynthesis[30, 31]. This relative decline in transpiration causes increased water-use efficiency and an increase in the fixation of 13C relative to 12C. Δ13C decreases with increasing elevation[28], likely related to a combination of atmospheric and edaphic factors[32], but the mechanisms are unclear and investigation continues[22, 33].

The existence of broad relationships between foliar isotopic composition and environmental factors suggest that these foliar chemical traits index the availability of nitrogen and n class="Chemical">water to plants[22]. However, different species grown in common environments exhibit variation in foliar δ15N and Δ13C due to distinct nutrient acquisition strategies and physiological response to the environment among species[19, 22]. For example, variation in nitrogen source (e.g., NO3 −, NH4 +, or dissolved organic nitrogen (DON))[22] and mycorrhizal association between species influences foliar δ15N[10]. Similarly, C3 plants vary in foliar Δ13C due to variation in carboxylation efficiency and stomatal response to environmental conditions[22]. However, it is an outstanding question if these sources of variability among species strengthen or weaken the global to continental scale variability in foliar δ15N and Δ13C imposed by climate. The assembly of species adapted to environmental conditions at any given site, might be expected to reduce the effect of climate on foliar chemical traits. Alternatively, gradients in plant lifeform along climate or elevational gradients, such as from deciduous angiosperms to evergreen gymnosperms that exhibit divergent mean foliar Δ13C, could strengthen effects of precipitation and elevation[28]. To advance our understanding of the relative effects of climate and plant traits (which vary by species) on foliar isotopic composition, widely distributed observations of foliar δ15N and Δ13C in a single species are required. Yet, there remains a paucity of examples where foliar δ15N and Δ13C were measured in more than a small number of specimens (e.g., 100) of the same plant species in its natural environment e.g.,[10, 19, 28]. Furthermore, because past work has analyzed foliar δ15N and Δ13C separately, the effect of N availability on Δ13C (e.g., through the influence of foliar nitrogen concentration on carbon assimilation) remains understudied[34].

Here we measure 755 specimens of the widely-distributed tree species Populus balsamifera for foliar [n class="Chemical">N], δ15N and Δ13C. Using these data, our objective was to explore the extent to which large gradients in climate-associated variation in foliar δ15N and Δ13C observed globally e.g., refs 10, 28, are represented in measurements of a single species. Furthermore, it has been proposed that intraspecific variation in foliar traits can either correlate with broad-scale trends along climate gradients, show no variation with climate, or exhibit an intermediate response[4]. While global interspecific variation in foliar traits is widely recognized as being valuable for understanding plant-nutrient relationships[35], few studies have sampled a sufficient number of individuals to conduct the same analysis using intraspecific variation. Because interspecific variation in plant traits is generally stronger than intraspecific variation, we hypothesize that the global variation in foliar [N], δ15N and Δ13C of P. balsamifera will exhibit a weaker response to climate than observed in global samples collected from many species. To evaluate the extent to which patterns in P. balsamifera chemical traits correspond to patterns in plant functioning, we subsequently used structural equation modeling to compare the effects of climate variables and foliar [N] on Δ13C.

Results

Foliar δ15N

Across the range of sites sampled (Fig. 1), δ15N in n class="Species">P. balsamifera ranged from −20.6‰ to 8.8‰. Climates for sampled P. balsamifera trees spanned ~12 °C of mean annual temperature (MAT) (−1.9 °C to 10.2 °C) and ~2100 mm of mean annual precipitation (MAP) (390–2566 mm). In a model of foliar δ15N with MAT, log(MAP), and log[N] as model effects, the addition of P. balsamifera observations (n = 755; enlarging the global data set by 7%) exhibited little influence on model estimates compared with using the original global data set (n = 9828) (Table 1). Foliar δ15N in P. balsamifera generally overlapped with foliar δ15N measurements in the global database, though within the restricted range of MAT and MAP exhibited by the P. balsamifera samples (Fig. 2). In general, the addition of P. balsamifera samples did not alter global patterns found previously: (1) foliar δ15N increased with increasing MAT (particularly for MAT >0 °C), (2) foliar δ15N decreased with increasing log-transformed MAP, and (3) foliar δ15N increased with increasing log-transformed [N] (R2 = 0.46 for both models; Table 1). Foliar δ15N of the P. balsamifera samples alone were generally more weakly related to the effects modeled than were plants from the global data set (Table 1). However, in contrast to the pattern observed for the global dataset, for P. balsamifera, foliar δ15N decreased with increasing MAT (P = 0.0005). P. balsamifera foliar δ15N showed no relationship with log-transformed MAP (P = 0.6). The explanatory ability of mean climate factors and foliar [N] was low (R2 = 0.05), leaving considerable small-grain variability remaining in the P. balsamifera δ15N observations.

Figure 1

Locations of foliar P. balsamifera samples acquired for this study mapped over mean annual precipitation and mean annual temperature for North America. Inset shows typical fine-grain variability in foliar Δ13C and δ15N. This map of North America was produced in ArcGIS 10.2.2 (http://www.esri.com/) using Bioclim mean annual temperature and annual precipitation layers[50]. The bivariate legend was produced using custom Python code (https://www.python.org/) developed by SMG.

Table 1

Models of foliar δ15N using each of the three data sets.

	Effect	SS1	Estimate ± SE2	P
Global data; R2 = 0.46	Intercept	3721	5.38 ± 0.27	<0.0001
	logMAP3	10561	−3.26 ± 0.10	<0.0001
	MAT4	65767	0.24 ± 0.003	<0.0001
	Log[N]5	24933	7.67 ± 0.15	<0.0001
Global data with P. balsamifera; R2 = 0.46	Intercept	3885	5.45 ± 0.27	<0.0001
	logMAP	11103	−3.30 ± 0.10	<0.0001
	MAT	70335	0.25 ± 0.003	<0.0001
	Log[N]	25031	7.61 ± 0.15	<0.0001
P. balsamifera data; R2 = 0.05	Intercept	16.1	−3.51 ± 2.69	0.1916
	logMAP	2.39	0.49 ± 0.97	0.6148
	MAT	115	−0.21 ± 0.06	0.0005
	Log[N]	185	5.17 ± 1.17	<0.0001

1Sum of squares. 2Standard error. 3Base 10 logarithm of Mean Annual Precipitation (mm). 4Mean annual temperature (°C). 5Base 10 logarithm of foliar nitrogen concentration.

Figure 2

Relationship between residual δ15N with mean annual temperature (MAT), mean annual precipitation (log transformed MAP), and nitrogen concentration ([N]). Residual δ15N is calculated relative to a model including the two remaining model effects. Globally distributed observations from 1273 species in black (n = 9828); observations of P. balsamifera (n = 755) in grey.

Relationship between residual δ15N with mean annual temperature (n class="Disease">MAT), mean annual precipitation (log transformed MAP), and nitrogen concentration ([N]). Residual δ15N is calculated relative to a model including the two remaining model effects. Globally distributed observations from 1273 species in black (n = 9828); observations of P. balsamifera (n = 755) in grey.

The unique patterns in relationships between foliar δ15N and clin class="Disease">mate for P. balsamifera were also present when subsetting the global data to the same climate envelope as P. balsamifera. Repeated random samples of foliar δ15N, N, MAT and MAP from the global database within the climate envelope defined by the P. balsamifera samples resulted in model estimates for each effect on δ15N that were different than the model estimates using P. balsamifera alone (Fig. 3). The P. balsamifera samples alone resulted in estimates of −0.21‰ °C−1 for MAT, 0.49‰ for logMAP, and 5.17‰ for log[N] (Table 1). After accounting for variation imposed by the differing mycorrhizal fungi associations of species, mean estimates resulting from 1000 repeated samples from the global dataset were −0.03 ± 0.001‰ °C−1 for MAT, −0.77 ± 0.02‰ for logMAP, and 6.71 ± 0.01‰ for log[N]. These effects of MAT, logMAP, and log[N] were significant (P < 0.05) 7.8%, 20.3%, and 100% of the time, respectively. Therefore, relative to the global data set, the effect of MAT on the P. balsamifera samples was more negative (P < 0.001), the effect of logMAP was removed (i.e., closer to zero; P = 0.6), and the effect of log[N] was reduced, although still significant (P < 0.001).

Figure 3

Distribution of model estimates on foliar δ15N for MAT, logMAP, and log[N] across 1000 samples of 755 observations in the global N database. Models were constructed as in Fig. 2. Samples were constrained to be within the climate range of P. balsamifera samples; the estimates of each effect using only P. balsamifera samples is represented by the vertical dotted line (also provided in Table 1).

Distribution of model estimates on foliar δn class="Chemical">15N for MAT, logMAP, and log[N] across 1000 samples of 755 observations in the global N database. Models were constructed as in Fig. 2. Samples were constrained to be within the climate range of P. balsamifera samples; the estimates of each effect using only P. balsamifera samples is represented by the vertical dotted line (also provided in Table 1).

Foliar Δ13C

Across the range of sites sampled, foliar Δ13C in n class="Species">P. balsamifera ranged from 18.2‰ to 25.7‰. Elevation for sampled P. balsamifera trees ranged from sea level to 2092 m, which was 60% of the elevation range spanned by the global data set (sea level to 3500 m). After accounting for elevational differences between samples, residual Δ13C was plotted against logMAP (Fig. 4) and it was observed that foliar Δ13C of P. balsamifera tended to be higher than global samples with the same MAP. Mean P. balsamifera foliar Δ13C was 22.2‰ and mean global foliar Δ13C for the same MAP range was 20.6‰ (21.2‰ for angiosperms and 19.1‰ for gymnosperms). Further, after accounting for variation with MAP, high elevation P. balsamifera sites exhibited elevated residual Δ13C relative to the global samples at the same elevation. As a result, models using P. balsamifera samples attribute a smaller fraction of model variance to logMAP and sqrt(elevation) (Sum of Squares values in Table 2) for their effects on Δ13C compared with these effects calculated using the global data set.

Figure 4

Relationship between residual Δ13C and mean annual precipitation (log transformed MAP) and the square root of elevation. For the plot with mean annual precipitation, residual Δ13C is calculated relative to a model that only includes the effect of elevation. Similarly, for the plot with elevation, residual Δ13C is calculated relative to a model that only includes the effect of precipitation. Global observations from previous work[28] in black, P. balsamifera observations (this study) in grey.

Table 2

Models of foliar Δ13C using each of the three data sets.

	Effect	SS1	Estimate ± SE2	P
Global data	Intercept	258	8.98 ± 0.89	<0.0001
	logMAP3	685	4.28 ± 0.26	<0.0001
	sqrt(Elevation)4	174	−0.05 ± 0.01	<0.0001
Global data with P. balsamifera	Intercept	1435	15.50 ± 0.67	<0.0001
	logMAP	379	2.44 ± 0.20	<0.0001
	sqrt(Elevation)	348	−0.06 ± 0.01	<0.0001
P. balsamifera data	Intercept	445	18.40 ± 1.20	<0.0001
	logMAP	23.1	1.36 ± 0.39	<0.0001
	sqrt(Elevation)	13.1	−0.02 ± 0.01	<0.0001

1Sum of squares. 2Standard error. 3Base 10 logarithm of Mean Annual Precipitation (mm). 4Square root of elevation (m).

Relationship between residual Δ13C and mean annual precipitation (log transformed MAP) and the square root of elevation. For the plot with mean annual precipitation, residual Δn class="Chemical">13C is calculated relative to a model that only includes the effect of elevation. Similarly, for the plot with elevation, residual Δ13C is calculated relative to a model that only includes the effect of precipitation. Global observations from previous work[28] in black, P. balsamifera observations (this study) in grey.

To examine how climate and n class="Chemical">N availability jointly influence Δ13C in P. balsamifera, we built a structural equation model (Fig. 5). Results exhibited a positive direct effect of log(MAP) (standardized estimate = 1.74; P < 0.001) and a negative direct effect of log(foliar [N]) on Δ13C (standardized estimate = −3.23; P < 0.001). A positive direct effect of log(MAP) (standardized estimate = 0.07; P = 0.025) and a negative direct effect of MAT (standardized estimate = −0.01; P = 0.001) on log(foliar [N]) represent indirect effects of climate on Δ13C through foliar [N]. The direct effect of MAT on Δ13C was not significant, therefore, the only effect of MAT on foliar Δ13C is indirectly through foliar [N].

Figure 5

Path diagram illustrating standardized effects (either positive or negative) of mean annual temperature (MAT), mean annual precipitation (log transformed MAP), and foliar N concentration (log transformed) on carbon isotope discrimination in leaves (Foliar Δ13C). Arrow widths are proportional to effect sizes; black arrows denote significant effects.

Path diagram illustrating standardized effects (either positive or negative) of mean annual temperature (MAT), mean annual precipitation (log transformed MAP), and foliar n class="Chemical">N concentration (log transformed) on carbon isotope discrimination in leaves (Foliar Δ13C). Arrow widths are proportional to effect sizes; black arrows denote significant effects.

Discussion

Foliar δ15N response to climate

Through a comparison of the effect of mean climate variables on δn class="Chemical">15N between a global data set, including hundreds of different plant species, and a single widely distributed species, P. balsamifera, we found that broad scale patterns in the isotopic composition of foliar tissues observed globally is largely weakened when observed in a single species. Across global gradients in climate, foliar δ15N measured in a wide variety of species has been shown to correlate with MAT, MAP, and [N], and mycorrhizal fungi association[10]. Demonstrating the same for a single species is often difficult due to the limited climatic variability encompassed by the range of most species. For the widely distributed P. balsamifera, we found weak, but significant effects of MAT and foliar [N], but not MAP (Table 1). The small effect of climate on foliar δ15N in P. balsamifera is consistent with the idea that the global patterns in foliar δ15N published previously were a result of geographic variation in species presence collinear with geographic variation in climate. However, the more limited climate range occupied by P. balsamifera complicates this interpretation.

Random samples from the global database within the climate envelope defined by the range of n class="Species">P. balsamifera (a sample pool of 1937 observations) exhibited significant effects of MAT, MAP and [N] on δ15N 7.8%, 20.3%, and 100% of the time (respectively, across 1000 resampling). Comparing global and P. balsamifera samples, the largest difference is the effect of MAT, which had a significant negative effect on foliar δ15N in P. balsamifera. These patterns afford some speculation as to the mechanisms behind global foliar δ15N patterns. For example, increasing foliar δ15N with increasing MAT between −0.5 and 30 °C is thought to represent a gradient in the dominant pathway responsible for the loss of N from ecosystems, from DON leaching (a low fractionation pathway) at low MAT to denitrification (a high fractionation pathway) at high MAT[10]. For P. balsamifera, occupying a MAT range from −1.9 to 10.2 °C, the enrichment of 15N decreased with increasing MAT, which is counter to the broad global pattern, but consistent with the 7.8% of models from the global database that were significant across this temperature range. For P. balsamifera, therefore, one possible interpretation is that N availability decreases with increasing temperature. However, it is also possible that P. balsamifera becomes more reliant on NO3 - as opposed to NH4 + with increasing temperature or is more likely to occupy sites with lower nutrient availability, leading to a small decline in foliar δ15N with increasing temperature[22].

After accounting for MAT, the effect of MAP on foliar δn class="Chemical">15N was weaker for P. balsamifera than for the global samples. The global pattern of decreasing foliar δ15N with increasing MAP has been explained as a response to higher gaseous N loss in more xeric environments[36] and low N availability (potentially due to high nitrate leaching) in wet environments[23, 37]. For the relatively mesic environments studied here, there was no effect of MAP, suggesting that precipitation amount is potentially not a driver of denitrification for the sites inhabited by P. balsamifera across this gradient. This might be a consequence of the habitat preferences of P. balsamifera, which is for wet locations such as riparian zones and the edges of wetlands, which might act to normalize the effect of precipitation on δ15N across its range.

The effect of foliar [N] on foliar δn class="Chemical">15N was smaller for P. balsamifera than it was for the global samples. After accounting for MAT and MAP, for every 1% increase in foliar [N] δ15N increased 5‰. Considering a foliar δ15N range of ~30‰, N availability appears to be highly variable in P. balsamifera, possibly associated with variation in productivity across its range. However, it is also clear that there is considerable site-to-site variability in N availability unexplained by the mean climate factors studied here. This variability is likely due to any number of factors, such as site soil moisture conditions, variable anthropogenic inputs, and site disturbance history[19].

Joint influence of precipitation and foliar nitrogen on Δ13C

Adding P. balsamifera samples to the global database reduced the effect of MAP and elevation on foliar Δn class="Chemical">13C. The effect of elevation on Δ13C has been discussed as related to temperature, the partial pressure of atmospheric CO2 or O2, irradiance, or edaphic factors[28], but is still a poorly understood process[33]. There is also a trend towards greater dominance of evergreen gymnosperms at higher elevations, which have lower Δ13C values than angiosperms. Diefendorf et al.[28] calculated that at the site level, evergreen gymnosperms on average exhibited 2.7‰ lower Δ13C than deciduous angiosperms. Therefore, not accounting for plant functional type in models will likely exaggerate the effect of elevation on Δ13C, but this is just one possible influence of elevation on Δ13C. By considering only P. balsamifera samples we explicitly account for interspecific variation, effectively removing most of the effect of elevation on Δ13C (Table 2). Beyond all coming from the same species, the samples of P. balsamifera also have the potential to normalize for habitat conditions typical of locations occupied by this species. Therefore, the reduced sensitivity of P. balsamifera Δ13C to elevation lends support to the idea that species turnover and variation in the site conditions encountered along elevation gradients explain the mechanisms behind the effect of elevation on Δ13C more than atmospheric differences.

The positive effect of MAP on Δ13C is understood as a result of reduced ston class="Disease">matal conductance at low MAP, which decreases water loss through transpiration more effectively than it decreases assimilation, thus increasing water use efficiency. The fact that we see this effect in P. balsamifera samples is consistent with the idea that trees from dryer locations exhibit greater water use efficiency. However, δ15N, and therefore N availability, trended lower with increasing MAP (Fig. 2), which has a similar effect on Δ13C. This is commonly seen in terrestrial plants: plants with high N availability exhibit increased assimilation, which reduces carbon discrimination (Δ13C) at the same stomatal conductance[34]. For a single tree species with presumably similar nitrogen and water acquisition strategies across its range, the effects of increased N availability and reduced stomatal conductance at low precipitation sites have the effect of reducing Δ13C. We see these direct and indirect effects on Δ13C in our structural equation model (Fig. 5), which shows significant direct effects of foliar [N] and logMAP on Δ13C. MAT only influences Δ13C indirectly through its effect on foliar [N]. The explanatory power of this model using P. balsamifera samples alone is small, however, suggesting that differences among species dominate the global pattern.

In summary, we generally found different responses to climate in plant foliar chemical traits measured in a single widely distributed species than the same measured in many species distributed globally. Foliar chemical traits in n class="Species">Populus spp. have been found to vary across distinct genotypes[38-40]. However, we find that this variation is not always congruent with variation due to species and plant functional type distributions. This is not to say that global variation in foliar isotopes, driven by species turnover across climate and elevational gradients, is in some way misleading. For many applications (e.g., interpreting paleo records[26]), global variation across many species is exactly what is needed. By studying a single species, we control for one source of variability on foliar isotopes, providing insight into mechanisms that are unique to that species or the habitat it occupies. Specific to δ15N, measurements of P. balsamifera exhibit a negative effect of MAT that was not measured across a similar climate range in globally distributed foliar samples. Among the possible interpretations is that N availability to these trees is reduced along the southern range edge of this species. Likewise, for foliar Δ13C, the effect of elevation was reduced in P. balsamifera, which suggests that the habitat’s P. balsamifera occupies somehow normalize for the effect of elevation on foliar Δ13C. It is difficult to use these observations of P. balsamifera to generalize to other species or mechanisms driving geographic variability. However, it seems clear that global variation in foliar δ15N and Δ13C is driven by variation between species, each with its own unique relationships to climate and elevation. The importance of this finding will only be realized as further research attempts rectify intraspecific and interspecific isotopic patterns in a variety of species, ideally leading to stronger generalities of how climate affects resource availability.

Methods

Populus balsamifera is a dominant deciduous tree across most of the boreal forest of n class="Chemical">North America (Fig. 1). In the southern edge of its extent (throughout the forests of New England and the upper Midwest, USA), it occurs in low density, preferring wetter sites and sites exhibiting recent disturbance. P. balsamifera generally exhibits ectomycorrhizal associations, however, there is some variation in mycorrhizal association across the Populus genus[41]. In the northern portion of its range, it is often the only deciduous tree or is codominant with Populus tremuloides. Studies of genetic variation in P. balsamifera show that the structure of population diversity reflects broad patterns of geographic expansion following the last glacial maximum, with three main demes evident located in the northwest, center, and eastern portions of its range[42, 43]. If intraspecific differences control variation in chemical foliar traits in P. balsamifera then we expect δ15N and Δ13C values to cluster into similar geographic regions. However, such patterns could also be superseded by global patterns in foliar traits driven by environmental gradients[10].

We sampled 755 specimens of P. balsamifera during the 2015 growing season. Samples were spaced >20 km apart along navigable roads to facilitate access, which also lead to a bias towards preferentially sampling within ~100 m of roads and across the southern range edge of this species (where there is greater access) (Fig. 1). While it is known that traffic emissions can influence foliar δn class="Chemical">15N in road-side plants[17, 44], we had no data on traffic patterns and ultimately chose not to investigate this effect in our samples. If an effect of traffic is present in our data, we suspect it would not be collinear with climate or elevation. For each sample, a hand-held global positioning system receiver was used to record the geographical location (including elevation) of the sample with a locational uncertainty of <30 m. Sites spanned nearly 23 degrees of latitude (39.87°N to 62.74°N) and 90 degrees of longitude (60°W to 150°W) (Fig. 1). Foliar samples were selected from sun leaves on each tree, acquired via a pruning pole or equivalent. Each sample was placed in a plastic bag with 10g of silica gel to act as a desiccant during transport to the laboratory where foliar samples were then dried at 60 °C for 48 hours. Leaves were homogenized, and approximately 2 mg of foliar tissue from each sample was analyzed for [C], [N], δ13C and δ15N using a Carlo Erba NC2500 elemental analyzer (CE Instruments, Milano, Italy) interfaced with a ThermoFinnigan Delta V + isotope ratio mass spectrometer (IRMS; Bremen, Germany). A MgClO4 trap was used to remove water vapor prior to the transfer of sample gases to the IRMS. The δ13C and δ15N data were normalized to the VPDB and AIR scales, respectively, using a two-point normalization curve with internal standards calibrated against USGS40 and USGS41. The long-term analytical precision (1σ) of an internal leaf standard analyzed alongside samples was 0.28‰ for δ13C and 0.24‰ for δ15N. Carbon isotope discrimination against 13C (Δ13C) was calculated according to Farquhar et al.[29] as:where the δ13Cair value used was the δ13C measured in the atmosphere at Mauna Loa, Hawaii in June 2015 (δ13Cair = −8.6‰). Elevated Δ13C signifies greater increases in intercellular CO2 concentrations (Ci) than atmospheric CO2 concentrations (Ca), which is caused by increased isotopic fractionation due to either greater stomatal conductance and/or lower photosynthetic assimilation rates in C3 plants[45]. Foliar [N] was log transformed to achieve a normal data distribution for comparison with δ15N and climate parameters.

We attained gridded climate data from DAYMET 1980–2015, which is an interpolated product providing minimum and maximum temperature and precipitation totals at daily intervals[46, 47]. Daily mean temperatures were calculated as Tmean = (Tmin + Tmax)/2. From the daily means, mean annual temperature (n class="Disease">MAT) was calculated by first calculating annual means 1980–2015, which were then averaged. Similarly, mean annual precipitation (MAP) was calculated as the mean annual sum of daily precipitation, 1980–2015. MAP was log10 transformed to achieve a normal data distribution. The square root of elevation was calculated to achieve a normal data distribution.

A global data set of foliar [N] and δn class="Chemical">15N (n = 9828) was acquired from Craine et al.[10]. Therefore, the number of measurements of δ15N of P. balsamifera was 7.7% of the global data set. The global data set included data on mean annual temperature (MAT), mean annual precipitation (MAP) both extracted from[48], geographic location, and species. The global data set represented foliar samples from across the full spectrum of global climates, whereas the P. balsamifera data were constrained to climate conditions in northern North America. A global data set of foliar Δ13C in trees (n = 570) was acquired from Diefendorf et al.[28]. These data represent species mean Δ13C at each site (between 1 and 227 individuals were sampled, n = 3310), rather than individual observations, thus removing within-species variability at each site. This data set included site MAP and elevation.

After completing a general survey of the data and evaluating the geographic and climatic range represented by the data, we followed procedures described by Craine et al.[10] and Diefendorf et al.[28] to build models relating δn class="Chemical">15N and Δ13C to climate factors. We combined the global and P. balsamifera datasets and built ordinary least squares models of foliar δ15N and Δ13C with a selection of MAT, logMAP, log[N] (for δ15N), and logMAP and sqrt(elevation) (for Δ13C) as model effects. Models using the (1) global data sets, (2) the global data sets combined with the new P. balsamifera observations, and (3) the P. balsamifera observations alone were compared to investigate the impact of the addition of the P. balsamifera samples on the global relationships. To facilitate visualization, for each of the models using the combined data sets, we iteratively built models that left one predictor out and then plotted the model residuals against the remaining predictor. For δ15N, the global data set contained a sufficient number of observations within the climate range exhibited by the P. balsamifera data (n = 1937) to estimate the effects of MAT and MAP on foliar δ15N across samples spanning a similar range in these predictors. Therefore, we repeatedly (x1000) sampled 755 random observations from the global database constrained to have the same minimum and maximum climate values as exhibited by the P. balsamifera samples. To account for variation in mycorrhizal association between the many species represented in the global foliar δ15N database, we used normalization factors published by Craine et al.[10] and applied an offset to each measurement so that mean foliar δ15N was the same for ericoid (+2.7‰), ectomycorrhizal (0‰), arbuscular (−1.2‰), and non-mycorrhizal (−3.2‰) plants. Using the 755 selected observations, we ran a model with MAT, logMAP, and log[N] as model effects and recorded the estimates for each predictor. We then repeated the model on each of the 1000 random resampled data sets and calculated the mean and standard deviation of the resulting model estimates for comparison with estimates of the same predictors on P. balsamifera foliar δ15N alone. There were insufficient samples in the global database to complete the same random resampling procedure for Δ13C.

The impact of global patterns of N availability and resulting foliar n class="Chemical">N on water relations in plants (and foliar Δ13C) is an outstanding question, particularly among widely distributed plant species that exist across a wide range of climate conditions. However, such relationships have rarely been evaluated for a single, widely distributed species. Therefore, we explored the importance of direct and indirect effects of MAT, MAP, and log[N] on Δ13C using a structural equation model (SEM). The model we used contained no latent variables, but did include all possible pathways between climate variables and the foliar chemical variables. This model was implemented in the R programing language using the Lavaan package version 0.5–22[49].