Great Basin Bristlecone Pine Volatiles as a Climate Change Signal Across Environmental Gradients

Alpine treeline species, like Great Basin bristlecone pine (GBBP) (Pinus longaeva Bailey), have received attention for their potential as indicators of climate change. Most studies have focused on climate-induced changes to treeline position, but climate effects on the physiology and stress of treeline plants remain poorly understood. Volatile organic compounds (VOCs) could provide insights into plant responses to climate change since the quantity and blends of VOCs released by plants exhibit variation in response to the environment, and can convey detailed information about the status of the emitting plant. We collected and analyzed GBBP VOCs and within-needle chemistry along elevational gradients (lower treeline, upper treeline, and midway in between) near the northern and southern geographic limits of GBBP. Random Forest analysis distinguished elevation classes using VOCs with 83% accuracy and identified the compounds most important for classification. Ordination revealed that elevation, temperature, heat load index, and relative humidity were significantly correlated with VOCs. Within-needle chemistry provided less predictive value in classifying elevation class (73% accuracy) and was correlated only with elevation and heat load index. These findings suggest that GBBP VOCs are highly sensitive to the environment and could be used to assess and predict tree status and responses to environmental change. The potential effects of climate- and elevation-induced changes in GBBP chemistry on abiotic and biotic interactions are discussed.


INTRODUCTION
Treeline species that live in harsh, high-elevation environments are particularly susceptible to climate change and can serve as early signals of change (Körner, 1998(Körner, , 2012. This is thought to be because trees growing near the treeline are at, or near, their survival limits and so should be highly sensitive to environmental change and respond earlier than the rest of the forest (Smith et al., 2009;Körner, 2012). Components of climate change such as temperature, elevated concentrations of CO 2 , and altered precipitation can each affect alpine treeline species, with warming being the best studied Smith et al., 2009;Körner, 2012). For example, warming temperatures can alter alpine treelines by increasing tree mortality, promoting invasive plant establishment, changing forest fuels, plant community structure, and altering snowfall and snow melt patterns (Flannigan et al., 2000;Gibson et al., 2008;Balch et al., 2013). Warming temperatures can also increase risk of mortality to treeline species by stimulating abiotic and biotic disturbances such as mountain pine beetles (MPB, Dendroctonus ponderosae Hopkins), white pine blister rust (Tomback and Resler, 2007), and altering natural fire regimes (Gibson et al., 2008;Gray and Jenkins, 2017). To date, most studies that have examined climate effects on alpine treelines have focused on the advance or retreat of treeline position (Gehrig-Fasel et al., 2007;Paulsen and Körner, 2014;Schibalski et al., 2014;Millar et al., 2015;Bruening et al., 2017). Climate effects on the physiology of treeline plants have received much less attention and remains largely unknown, despite the fact that understanding the impacts of environmental variability on plant physiology could help predict how treelines will change.
Great Basin bristlecone pine (Pinus longaeva Bailey) (GBBP) is an alpine treeline, five-needle pine that has served as an important indicator of climate change (Rochefort et al., 1994;Körner, 1998;Salzer et al., 2009Salzer et al., , 2014. GBBPs are the oldest trees on earth, with some individuals exceeding 5000 years old. This longevity has allowed GBBP to play an important role in climate change research due to their extremely long tree-ring chronologies (Cook et al., 1995;Cook and Peters, 1997) with recent studies finding increased tree-ring growth at the highest elevations (Salzer et al., 2009(Salzer et al., , 2014. GBBP distribution is confined to the highest mountains (2,700-3,700 m) of the Great Basin in the western United States, where global climate models predict a mean 2-4 • C increase in annual temperatures within the next several decades [Scalzitti et al., 2016; IPCC Fifth Assessment Report (AR5)] accompanied by a likely decrease in precipitation (Cook et al., 2010). The physiological responses of GBBP to warming temperatures and decreases in precipitation are largely unknown. Moreover, such warming could add stress by increasing wildfire activity (Westerling et al., 2006) in these habitats with historically sparse fuel conditions, in part by increasing fine woody surface fuels (Gray and Jenkins, 2017).
Volatile organic compounds (VOCs) could be useful for detecting and measuring plant physiological responses to environmental change because the quantity and composition of VOCs emitted by a plant can be affected by the environment (Dudareva et al., 2006;Peñuelas and Staudt, 2010;Jaeger et al., 2016). Plant VOCs have many known ecological roles, such as attracting pollinators (Burkle and Runyon, 2016) and plant defense against herbivores, pathogens, and parasitic plants (Runyon et al., 2006;Huang et al., 2012;Gray et al., 2015;Turlings and Erb, 2018). VOCs emitted by GBBP are known to play important roles in defense, for example by repelling host-searching MPB (Gray et al., 2015). Plant VOCs can also protect against certain abiotic stresses, including high temperatures (Duhl, 2008), oxidative damage (Holopainen and Gershenzon, 2010), and cold stress (Cofer et al., 2018). In general, both biotic and abiotic stress increases VOC emissions from plants (Holopainen and Gershenzon, 2010), and trees are known to alter emission of VOCs in response to herbivory and changes in heat, light, precipitation, and season (Helmig et al., 2007;Trowbridge et al., 2014). Moreover, elevated temperatures typically increase VOC emissions (Tingey et al., 1980;Peñuelas and Llusià, 2003;Jamieson et al., 2017), and elevated CO 2 can increase emission of volatile terpenoids Yuan et al., 2009;O'Neill et al., 2010).
This study quantifies GBBP tree chemistry across elevational gradients near the southern and northern extent of GBBP's distribution as an approximation for future climate conditions. Elevation gradients are valuable surrogates for inferring broader climate change effects by providing variation in abiotic factors (Hodkinson, 2005;Körner, 2007;Beier et al., 2012). Moreover, because the elevation gradients at each latitude used in this study occurred over short distances (<2.5 km), this minimized the confounding effects of biogeographical differences such as the community of plants, herbivores, and pathogens that are present (Hodkinson, 2005). We address the following questions in this paper: • How does GBBP chemistry change with elevation, and is this correlated with air temperature, relative humidity (RH), and heat load index? • Do VOCs emitted from GBBP increase at lower elevations as a proxy for climate warming? • Can we apply understanding of GBBP response to elevational gradients to interpret potential susceptibility to threats like fire or MPB?

Study Sites
The Great Basin of the western United States has an arid climate and basin and range topography, with numerous high mountain peaks (sky islands) separated by low intermountain valleys. In the Great Basin, upper treeline occurs at ∼3,300 m and lower treeline at ∼2,100 m ( Thompson and Mead, 1982), although this varies with latitude. Mountains and basins create steep environmental gradients, which greatly influence the composition and structure of vegetative communities (Peet, 2000). For this study, we exploited these environmental gradients by sampling GBBP trees in two sky islands in Nevada near the northern and southern extremes of GBBP's geographic range. Sample transects were installed at Cave Mountain, NV (39.167, −114.616) and the Spring Mountains National Recreation Area, NV (36.293, −115.686) ( Figure 1A). Three discrete elevation classes, "low" "mid" and "high, " were sampled at each site. At each site and elevation, four GBBP trees of similar size were haphazardly selected that showed no obvious signs of stress (e.g., herbivory, pathogen attack). Mean height of sampled trees was 13.5 ± 0.8 m, and mean diameter at breast height (dbh) was 74.5 ± 7.8 cm (n = 24). Trees were sampled once each month from July to September in 2013, at each site, it took field crews 2 days to sample trees from the three elevations. VOCs were sampled three times daily (1,000, 12,00, and 1,400) in the first week of July, August, and September, needles for withinneedle terpene concentrations were collected from each tree at the end of the day (1,400). The July samples from the Spring Mtns were omitted from analysis because sampling was affected by the nearby Carpenter 1 fire (July, 2013, 36.25, −115.69, ∼4 km away) (total VOC collections: Cave Mtn = 108, Spring Mtns = 72). Other environmental variables measured at time of VOC sampling were temperature ( • C), relative humidity (RH in %), and topographic heat load index (HeatI). Heat load index, or the potential direct radiation at a site, is not symmetrical, as a slope with afternoon sun will be warmer than an equivalent slope with morning sun. Heat load index was calculated using the equation: HeatI = 1-cos(θ -45)/2, where θ = aspect in degrees east of north, as an approximation of heat, rescaling aspect to a scale of zero to one, zero being the coolest slope, and one having the most direct radiation (McCune and Keon, 2002). The only environmental variables measured for within-needle compounds were heat load index and elevation, as needles were collected once at the end of a sampling day.

Collection and Analysis of VOCs and Within-Needle Chemistry
We used GC-MS to examine chemicals emitted into the airspace surrounding plants (VOCs) and pulverized needles to examine the chemicals inside. In this paper, we use the term "withinneedle chemistry" to refer to analysis of pulverized needles and to distinguish from emitted VOCs. While GC-MS allows measurement of nearly all VOCs, this methodology did not permit us to measure all chemistry within-needles and amounts of larger within-needle compounds such as diterpenes, some phenolics, and sugars could not be determined in this study. Volatile emissions were collected by enclosing 50 cm of the apical end of one randomly selected branch with clear Teflon bags (50 cm wide x 75 cm deep, American Durafilm Co., Holliston, MA, USA). The apical branches were ∼1.5 m above the forest floor, following the methods of Page et al. (2012) and Gray et al. (2015). The same branch on each tree was marked with flagging and was used for all VOC sampling periods. Air was pulled from a side port (0.5 l min −1 ) of the Teflon bags through volatile traps containing 30 mg of the porous polymer adsorbent HayeSep-Q (Restek, Bellefonte, PA, USA) using portable volatile collection systems comprising automated vacuum pumps (Volatile Assay Systems, Rensselaer, NY, USA). VOCs were collected for 30 min and Teflon bags were removed between sample times. The volatile traps were stored on ice in a cooler, and then shipped immediately to the USFS RMRS in Bozeman, MT where they were processed within a week of sampling. VOCs were eluted from traps with 200 µl of dichloromethane, and 1,000 ng of n-nonyl acetate was added as the internal standard. After the last sampling period in September, the branches sampled for VOCs were clipped and needles removed and weighed.
To measure within-needle chemistry, ∼20 g of needles from randomly selected branches (1-2 m above the forest floor) on the same trees were removed at the end of each sampling date and stored on ice in a cooler (for 2-3 days) and then in a freezer at −80 • C until processed. Terpenoids were extracted from GBBP foliage following methods of Ormeno et al. (2009) and Page et al. (2014). For each sample, needles were thoroughly mixed and 2 g randomly selected and ground into a fine powder in liquid nitrogen using a mortar and pestle. Approximately 0.1 g of powdered needles were transferred into 2-ml FastPrep tubes (MP Biomedicals, Solon, OH), and 1.5 ml of cyclohexane was added and sonicated at room temperature for 20 min. Vials were then centrifuged at 13,000 g for 1 min and 200 ml of cyclohexane (top layer) was transferred to a gas chromatograph vial and 1,000 ng of n-nonyl acetate was added as the internal standard.
VOC and needle samples were analyzed using an Agilent 7890A gas chromatograph (GC) coupled with a 5975C mass spectrometer (MS) and separated on an HP-1 ms (30 m × 0.25 mm inside diameter, 0.25 µm film thickness) column, helium was used as the carrier gas. The GC oven was maintained at 35 • C for 3 min and then increased by 5 • C min −1 to 200 • C, then 25 • C min −1 to 250 • C. Quantifications were made relative to the internal standard using ChemStation software (Agilent Technologies, Wilmington, DE, USA). Identification of compounds were made using the NIST 08 Mass Spectral Search Program (National Institute of Standards and Technology, Gaithersburg, MD, USA) and confirmed by comparing retention times and mass spectra with commercial standards, when possible. Remaining unidentified compounds were labeled as unidentified monoterpenoids (MT), unidentified sesquiterpenoids (ST), unidentified benzenoids (B), or unidentified green leaf volatiles (GLVs). VOC emission rates were corrected for needle weight (ng hour −1 g −1 ) as were within-needle chemical concentrations which are reported on a fresh needle weight basis (µg g −1 ).

Analysis of Variance (ANOVA)
Statistical analyses for VOCs and within-needle compounds were performed using the non-parametric Kruskal-Wallis oneway analysis of variance (ANOVA) (Kruskal and Wallis, 1952) to identify compounds with significant (P ≤ 0.05) differences among elevations, and to test whether samples originate from the same distribution. To identify which elevation classes had significant differences in compounds, pairwise differences between all sampling periods were square root or log-transformed to normalize data and stabilize variance, and assessed with a Tukey's HSD means test (Sokal and Rohlf, 1995) and calculated using the Stats package in the R v.3.3.1 statistical software (R Development Core Team, 2016).

Random Forest
We used the Random Forest classification algorithm (Breiman, 2001) to investigate whether VOCs hold value for classifying elevation (low, mid, and high). Random Forest is a machinelearning algorithm that assigns samples to predefined groups and combines a consensus of multiple classification trees in numerous iterations and estimates the importance of each compound (Breiman, 2001). Constructing multiple classification trees using hundreds of bootstrapped training sets and averaging the resulting predictions significantly reduces the variance and increases prediction accuracy compared with a single tree. Additionally, when the number of variables is very high (in this case 42 compounds), applying only a randomly selected subset of variables for each model is computationally efficient. This statistical method has been applied in other ecological studies to classify VOC samples to reduce noise and correctly identify volatile signatures of plants (Ranganathan and Borges, 2010;Jaeger et al., 2016). Two Random Forest trials were carried out to classify elevation based on VOCs and within-needle chemistry. These analyses used random Forest package v4.6-12 (Liaw and Wiener, 2002) in R v.3.3.1 statistical software (R Development Core Team, 2016) utilizing the optional measure of importance of predictor variables, and proximity, a measure of the internal structure of the data used to detect outliers, with number of trees set at 1,000, and all other parameters were set as the defaults. Random Forest returns a confusion matrix that summarizes the accuracy of the classification as well as the variable importance. The importance of each VOC for classification was ranked using mean decrease in accuracy (MDA) which measures the accuracy in which the compound can be used to partition the data variables. Random Forest has two measures of variable importance. The first is based on mean squared error (MSE) and relates to the prediction accuracy of the out-of-bag portion of the data. The difference between the two MSEs is then averaged over all trees. The second measure describes variable importance based on the Gini impurity index, which is based on overfitted models (Breiman, 2001). The variables predicted to be important in the model help us to understand what variables are driving the differences in chemical signatures at different elevations.

NMDS Ordination
To characterize VOC dissimilarity along environmental gradients, a non-metric multidimensional scaling (NMDS) ordination based on a matrix of Euclidean dissimilarities (Kenkel and Orlóci, 1986;Dixon, 2003) was calculated on the rank order proportion of VOCs using most important compounds determined based on MDA from Random Forest. Ordination is the ordering of sample units along an axis, multiple axes, or a gradient , and is most often used in ecology for gradient analysis and is used to select the most important factors, separate patterns, and reveal unforeseen patterns and processes . The goal of NMDS is to collapse information from multiple dimensions into just a few, so that they can be visualized and interpreted. The iterative NMDS procedure defines position of variables in multidimensional space, constructs an initial configuration in 2-dimensions, regresses distances in the initial configuration against the observed distances, determines the stress, or the disagreement between 2-D configuration and predicted values from the regression, and if stress is high, repositions the points in 2 dimensions in the direction of decreasing stress. The scatter of points representing samples in NMDS ordination space is iteratively allowed to evolve until it resembles the observed dissimilarity matrix as closely as possible. This method is nonmetric because stress, the measure of closeness of fit, is based on the ranking of dissimilarity values rather than actual values of dissimilarity. The NMDS was based on a similarity matrix using the Bray-Curtis index (Clarke, 1993) running with 600 iterations to permit the NMDS algorithm to develop an ordination with minimal stress. NMDS ordinations were completed using the VEGAN package (Dixon, 2003) in R v.3.3.1 statistical software (R Development Core Team, 2016). The ordination space was then visualized by overlaying environmental variables with the ordination to specifically address how the compounds respond to spatial variation in the natural environment.
Finally, regression coefficients of several ratios of compounds of interest were calculated using generalized linear models (GLM) relating the compound ratios to temperature gradients. Regression coefficients of flammable within-needle compounds were calculated for elevational sites and heat load index. All statistics were completed using R v.3.3.1 statistical software (R Development Core Team, 2016).

GBBP VOCs
Gas chromatography-mass spectrometry (GC-MS) analysis identified 42 volatile compounds emitted by GBBP trees (21 select compounds based on MDA are shown in Table 1; all compounds are shown in Supplementary Table 1). The majority of compounds were monoterpenoids (29 compounds, ca. 99% of total VOC emissions) followed by sesquiterpenoids (10 compounds) and benzenoids (3 compounds) (Supplementary Table 1). These 42 compounds were emitted by all trees across all sample dates and elevations. However, volatile emissions varied quantitatively between study sites and across elevations for total VOCs and for some individual compounds ( Table 1). VOCs increased with decreasing elevation for several compounds (highlighted in bold in Table 1 and Supplementary Table 1) and this trend was common among the sesquiterpenoids. Amounts emitted for every sesquiterpenoid compound decreased with increasing elevation for at least one of the sites. Bornyl acetate was the only monoterpene to exhibit this inverse elevational trend at both study sites, while emission of three unidentified sesquiterpenoids (ST1, ST2, and ST4) decreased with increasing elevation at both sites. Additionally, GLM regression analysis revealed 12 volatile compounds that significantly increased (P < 0.001) with temperature (Supplementary Table 3).

Within-Needle Terpenoid Concentrations
In total, 78 within-needle compounds were identified using GC-MS, the majority of which were monoterpenoids (46 compounds, ca. 91% of total terpene concentrations) and sesquiterpenoids (27 compounds, ca. 7.5% of total terpene concentrations) (21 select compounds based on MDA are shown in Table 2; all compounds are shown in Supplementary Table 2). Overall, within-needle concentrations varied more by site than by elevation. Mean total within-needle terpene concentrations averaged across all sample periods were not significantly different among elevation classes, but the concentrations were significantly greater (P < 0.05) at the Spring Mtns site in the south than the northern Cave Mtn site. Total within-needle terpene concentrations had high variability among elevation classes and study sites. As with emitted VOCs, some within-needle terpene compounds showed trends along elevational gradients with concentrations increasing with decreasing elevation ( Table 2 and Supplementary Table 2). Bornyl acetate again exhibited decreases in concentrations with increase in elevation at both sites, along with two unidentified monoterpenes (MT8 and MT). Two unidentified sesquiterpenes (ST2 and ST3) were correlated with elevational gradients at both sites, however, these compounds decreased with increasing elevation at the Spring Mtns site, and increased with elevation at Cave Mtn site.

Random Forest Classification of Elevation Using VOCs
Random Forest classification correctly assigned high, mid and low elevation sites using VOCs in 82.7% of samples ( Table 3). When holding back 10% of the sample to validate model classification, the accuracy averaged 83.3% (10 trials), slightly better than the initial out-of-box (OOB) accuracy. Multidimensional scaling (MDS) of the Random Forest model using the 42 VOCs showed strong clustering by elevations class (Supplementary Figure 1). The most important compounds for differentiating elevation classes were (E)-β-ocimene, αfarnesene, B2, MT 8, geranyl acetate, tricyclene, ST 7, linalool, bornyl acetate, γ -terpinene, and caryophyllene ( Table 4) which was partly supported by the Kruskal-Wallis-tests (Table 1).
Random Forest classification OOB accuracy for within-needle terpene concentrations was not nearly as good (73% correctly classified), with 25% of the high elevation sites misclassified, 30% of the middle, and 25% of the low sites misclassified. Linalyl acetate, ST 21, ST 4, longifolene, terpinylacetate, MT 21, and α-terpineol (Supplementary Table 4) were the most important variables in classifying elevation using within-needle chemistry.

Ordination of VOCs and Within-Needle Compounds
We used NMDS ordination to visualize elevational effects on VOCs composition (Figure 2). NMDS ordination of the VOC data had a high linear fit (R 2 = 0.98) and a low stress value (0.075). A stress of ca. 0.05 provides an excellent representation in reduced dimensions, while a stress of ca. 0.1 provides a good representation (Clarke, 1993). The distances between points in Figure 2 are proportional to the dissimilarity between the blends of VOCs for each tree. While there is high variability among elevations, the first axis of the ordination plot (NMDS1) was highly correlated with heat load index [positive], whereas the vertical axis (NMDS2) was highly correlated with elevation [positive], relative humidity [positive], and temperature [negative]. Mean temperatures for Spring Mountains at low, mid, and high were 16.6, 13.1, and 11.9 • C, respectively. We overlaid NMDS with contours of environmental variables (Supplementary Figure 2). The NMDS axis dissimilarities are more pronounced along the NMDS1 axis and the environmental variables may have driven VOC amounts. The dissimilarity of compounds along the NMDS1 is perpendicular to the contours of the environmental variables humidity (Supplementary Figure 2A) and temperature (Supplementary Figure 2B). While elevation did not drive the ordination like humidity and temperature, there is a clear similarity of points along the elevation contours (Supplementary Figure 2C).   The NMDS ordination of the within-needle compounds (Figure 3) also had a high linear fit (R 2 = 0.93) and a medium stress value (0.142), indicating a good representation of the data (Clarke, 1993). Sesquiterpenoids clustered in the lower values of NMDS1 and along the heat load index vector (Figure 3). The NMDS ordination plot shows a pattern with within-needle chemistry of trees from the highest elevation class distributed along the NMDS1 axis (gray triangles and circles in Figure 3), with the mid elevations points clustered centrally along NMDS2, and the lowest elevation trees exhibiting more dispersion/variability (Figure 3). NMDS2 was negatively correlated with heat load index and elevation was not significant in the within-needle ordination.

DISCUSSION
Results showed that GBBP foliage chemistry, especially VOCs, changed with elevation in predictable ways. The Random Forest algorithm correctly identified the elevation (low, mid, high) at which GBBP trees were growing using VOCs with 83% accuracy. Within-needle chemistry was much less predictive of elevation (only 73% accurate) than VOCs, possibly because terpenes accumulate within conifer needles (Schönwitz et al., 1990;Martin et al., 2003) and these larger stored pools might be slower to change or less responsive to the environment in general. We don't know if changes in VOCs are due to emission of stored or newly synthesized compounds, or a mixture of both. In contrast, VOC emissions are highly responsive to the environment and plants have considerable control over the identity and amount of VOCs released into the surrounding atmosphere (Widhalm et al., 2015;Adebesin et al., 2017). Interestingly, VOCs across elevations differed quantitatively and not qualitatively, thus each tree emitted an identical suite of VOCs, but elevation strongly affected the quantities of the compounds that were emitted. It is surprising that distinct VOC differences occurred, and were detectable, over such small elevation changes, suggesting there can be a notable changes in micro climate even over these short distances. For example, the difference between low and high elevation sites at the Spring Mountains was just 520 m and only 225 m at Cave Mountain, suggesting that GBBP VOCs are highly sensitive to relatively small elevation changes and the concomitant environmental changes. Although we cannot rule out a role for genetic differences, which can influence VOC emissions (e.g., O' Reilly-Wapstra et al., 2011;Jaeger et al., 2016), our finding only a change in emitted amount and not a change in individual compounds (sharing the same "chemotype") is an indication that the plants are genetically rather similar. The two most important VOCs for classifying elevation were the monoterpene (E)-β-ocimene and sesquiterpene α-farnesene, and emissions of both are known to be highly responsive to the environment (e.g., Joó et al., 2011;Farré-Armengol et al., 2017).
What environmental factors could be affecting GBBP VOCs across elevations? We found that air temperature, heat load index, and relative humidity were strongly correlated with emitted VOCs (Figure 2). Not surprisingly, temperature consistently decreased with increasing elevation, with only three cases (out of 45) in which a higher elevation site was warmer than a lower elevation site at the same time of day. Higher temperatures did not consistently increase total VOC emissions: the greatest total VOCs were emitted by trees at low elevation just 27% of the time. This suggests that factors other than temperature are influencing VOCs and overriding the generally positive effect temperature has on plant volatiles (Jamieson et al., 2017) or that temperature is altering the blend of VOCs rather than total amounts. The effects of heat load index on plant VOCs is not studied, but at sample sites heat load index and elevation were not correlated (Pearson's r = 0.06). Relative humidity can also affect plant VOCs (Vallat et al., 2005) and RH tended to be lowest at the low elevation sites (67% of the time). Random Forest was particularly accurate at identifying trees at the low elevation sites using VOCs (only 12% were misclassified) ( Table 3), possibly indicating that environmental conditions at these low sites are more stressful resulting in altered VOCs. Within-needle terpenoids were correlated only with heat load index and not elevation, suggesting that microsite (e.g., aspect, openness, etc.) variables play an important role in these chemicals. Both VOCs and within-needle chemistry are likely affected by environmental variables we did not measure, including precipitation, soil moisture, C0 2 , UV, and ozone to name a few.
While there is uncertainty in how climate change will affect mountain regions of the Great Basin, extensive research has modeled probable effects (Ziaco et al., 2014;Barbero et al., 2015;Millar et al., 2015;Smithers et al., 2018). For example, dynamical downscaling of Global Climate Models (GCMs) suggest the Great Basin will experience a 2-4 • C increase in March-April mean temperatures, 60-100% decrease in mean April snow water equivalent (SWE), yet possibly a 10-20% increase in October-April precipitation by the end of the twenty-first century (Scalzitti et al., 2016). Assuming that the upper treeline of GBBP is temperature limited (Körner, 2012), it is reasonable to expect higher stands to more closely approximate the climate conditions that are currently at the mid and low elevations (predicted 2-4 • C warmer on average). The GBBP living at lower treeline might experience higher mortality rates due to heat stress, drought stress, competition pressure from other tree species, and increased likelihood of fire (Gray and Jenkins, 2017). GBBP have often been used to analyze climate change effects on conifers and reconstruct climate histories. For example, recent unprecedented tree-ring growth in bristlecone pine at the highest elevations and possible causes were presented by Salzer et al. (2009). GBBP chronologies from near treeline show an increase in precipitation sensitive tree rings (similar to lower elevation GBBP chronologies), yet tree rings rarely correlate well with temperature (Salzer et al., 2014).
Because of the high requirement of photosynthetic carbon for terpene synthesis (Loreto et al., 1996), terpenes can be relatively expensive for a tree to synthesize and store (Lewinsohn et al., 1991;Niinemets, 2004;Peñuelas and Staudt, 2010). Altered terpene production at lower elevations, which FIGURE 3 | Non-linear Multidimensional Scaling (NMDS) plot of the first and second dimensions for 21 within-needle terpenoid concentrations (ng per gram; n = 60) of Great Basin bristlecone pines growing at low, middle and high elevations. Green vectors are environmental gradients that are significantly correlated (α = 0.01) to the terpenoid ordination. Points are tree samples from each date; colors, and shape indicate elevation class.
approximates climate change induced heat stress, could reduce tree productivity. Examples of abiotic stresses include extreme temperatures, drought, and wildfires, all of which are factors to GBBP. Abiotic stress is an integral part of climate change which has a wide range of uncertain impacts on plants and trees. A stress event (e.g., extremely hot days) could deplete monoterpenes stored within-needles and exhaust reserves (Kravitz et al., 2016). Kravitz et al. (2016) suggest that temperature stress before VOC sampling may affect which terpenes are available, which in turn may provide a potential "early warning" of stress to the ecosystem. There was a period of high temperatures preceding the July sample period and lower temperatures preceding September (Figure 5). These preceding weather events could affect the VOC emissions and available within-needle terpenes, and might explain the monthly differences of VOC emissions observed in the sample (Supplementary Figure 3). The temperatures at our two sites (approximately the northern and southern extents of GBBP distribution) are highly correlated with each other (Figure 5). The southern Spring Mtn site was about 5 • C warmer than the Cave Mtn site. Considering terrestrial vegetation is one the largest sources of volatile organic compounds to the atmosphere (Misztal et al., 2015), there is a need to assess tree chemical responses to these environmental changes, and the connection between organism-level stress responses and broader atmosphere-scale studies. These responses are an instrumental gap between atmospheric chemistry and ecosystem studies (Kravitz et al., 2016), and can further the understanding of air pollution and ozone pollution that are linked to biogenic VOC emission (Kammer et al., 2019).
We developed an additional method to analyze future climate change effects on GBBP, rather than relying only on tree rings/dendrochronology to understand historical precipitation-based climate changes. NMDS ordination showed that the environmental gradients of elevation, temperature, RH, and heat load index all are statistically correlated with VOC emissions (Figure 2 and Supplementary Figure 2). Our Random Forest classification model was highly accurate at predicting elevation classes from VOCs, demonstrating a distinct chemical signature of trees growing at different elevations. We also improved methods to reliably assess and predict tree resiliency with climate change. These changes along environmental gradients may be a response to biotic or abiotic stress. This research improves understanding of VOC emissions as a physiological tree response to environmental gradients and helps inform land managers about forests threatened by native and non-native pests and pathogens, increased threat of fire, and changing distribution patterns from climate change.

AUTHOR CONTRIBUTIONS
CG, JR, and MJ performed VOC collections and field measurements. JR analyzed the VOC and foliage chemistry. CG performed data analysis and wrote the first draft of the manuscript. All authors contributed to manuscript writing and revision.

FUNDING
This project was supported by funding from USFS EM Project INT-EM-F-10-02, Utah Agricultural Experiment Station (UAES) UTA Project 1070, USDA National Needs Graduate Fellowship Competitive NIFA Grant, the Ecology Center at Utah State University, and PECASE (President's Early Career Award in Science and Engineering) to JR. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

ACKNOWLEDGMENTS
The authors thank Andrew Giunta, Chelsea Toone, Scott Frost, and Sarah Null who helped collect volatiles. The staff at Cave Lake State Park, NV provided site access and camping. Some portions of this paper have appeared previously in the author's dissertation (Gray, 2017). We would also like to thank the two reviewers whose comments greatly improved this manuscript.