Differences in leaf gas exchange strategies explain Quercus rubra and Liriodendron tulipifera intrinsic water use efficiency responses to air pollution and climate change

Trees continuously regulate leaf physiology to acquire CO2 while simultaneously avoiding excessive water loss. The balance between these two processes, or water use efficiency (WUE), is fundamentally important to understanding changes in carbon uptake and transpiration from the leaf to the globe under environmental change. While increasing atmospheric CO2 (iCO2) is known to increase tree intrinsic water use efficiency (iWUE), less clear are the additional impacts of climate and acidic air pollution and how they vary by tree species. Here, we couple annually resolved long‐term records of tree‐ring carbon isotope signatures with leaf physiological measurements of Quercus rubra (Quru) and Liriodendron tulipifera (Litu) at four study locations spanning nearly 100 km in the eastern United States to reconstruct historical iWUE, net photosynthesis (Anet), and stomatal conductance to water (gs) since 1940. We first show 16%–25% increases in tree iWUE since the mid‐20th century, primarily driven by iCO2, but also document the individual and interactive effects of nitrogen (NOx) and sulfur (SO2) air pollution overwhelming climate. We find evidence for Quru leaf gas exchange being less tightly regulated than Litu through an analysis of isotope‐derived leaf internal CO2 (Ci), particularly in wetter, recent years. Modeled estimates of seasonally integrated Anet and gs revealed a 43%–50% stimulation of Anet was responsible for increasing iWUE in both tree species throughout 79%–86% of the chronologies with reductions in gs attributable to the remaining 14%–21%, building upon a growing body of literature documenting stimulated Anet overwhelming reductions in gs as a primary mechanism of increasing iWUE of trees. Finally, our results underscore the importance of considering air pollution, which remains a major environmental issue in many areas of the world, alongside climate in the interpretation of leaf physiology derived from tree rings.


| INTRODUC TI ON
Forest ecosystems represent 97% of the land carbon (C) sink in the conterminous United States and sequester 16% of U.S. national carbon dioxide (CO 2 ) emissions from fossil fuels into trees and soil (Lu et al., 2015;Woodbury et al., 2007). At the same time, forests also play a central role in the hydrologic cycle by regulating large fluxes of water through leaf stomata during transpiration (Schlesinger & Jasechko, 2014). Thus, environmental influences on trees can alter the balance between forest carbon uptake and water loss, with cascading effects over ecosystem function. Yet, these impacts can often be difficult to predict due to the complex, and sometimes counteracting, nature of environmental drivers on physiology (Fernández-Martínez et al., 2017;Guerrieri et al., 2019;Jennings et al., 2016;Mathias & Thomas, 2018). As forests currently contribute one of the largest sources of uncertainty in future climate projections (Friedlingstein et al., 2006), it is imperative to better understand the complexity behind the many factors simultaneously driving tree physiology to improve our predictions of future forest dynamics.
Indeed, reconstructions using tree rings have shown up to a ca. 40% increase in analogous tree-level intrinsic water use efficiency (iWUE; the ratio of net photosynthesis, A net , to stomatal conductance to water, g s ) since the early 1900s driven overwhelmingly by iCO 2 (Belmecheri et al., 2021;Frank et al., 2015;Guerrieri et al., 2019;Peñuelas et al., 2011). However, local environmental factors such as temperature (Way & Oren, 2010;Wu et al., 2011), precipitation (Beer et al., 2010;Belmecheri et al., 2021;Wu et al., 2011), and nitrogen availability (Horn et al., 2018;Thomas et al., 2010) may also independently influence A net and g s, thereby modifying the response of iWUE to iCO 2 .
Often overlooked, nitrogen oxides (NO x ) and sulfur dioxide (SO 2 ) can also have complex, nuanced impacts on forest ecosystems that differ among tree species (Horn et al., 2018;Maxwell et al., 2019;Thomas et al., 2010) and site environmental conditions (Smith et al., 2016), and are consequential to large-scale carbon (Fernández-Martínez et al., 2017) and water cycling (Lanning et al., 2019).
Moreover, these acidic air pollutants remain a major issue in many forested areas of the world (Aas et al., 2019;Akimoto et al., 2022) and are important to interpreting tree-ring records in North America  Savard et al., 2020). Sulfur pollution negatively impacts leaf physiology by reducing rates of photosynthesis (Choi et al., 2014) and stomatal conductance (Borer et al., 2005) which is reflected in tree iWUE (Rinne et al., 2010;Thomas et al., 2013). In contrast, nitrogen pollution may benefit iWUE through stimulation of A net resulting from increased nitrogen availability in both angiosperm and gymnosperm tree species (Gharun et al., 2021;Jennings et al., 2016;Leonardi et al., 2012;Mathias & Thomas, 2018). Indeed, reductions in N deposition in some regions of the world have reduced iWUE below the level otherwise expected under iCO 2 (Gharun et al., 2021;Mathias & Thomas, 2018), highlighting the need for a careful consideration of iWUE responses to acidic air pollution.
Dendroecology has gained substantial traction as an approach to examine the multi-decadal response of forest stands to the myriad of changing environmental factors, including in locations experiencing high acid pollutant loads (Boettger et al., 2014;Engel et al., 2016;Guerrieri et al., 2020;Kosiba et al., 2018;Lanning et al., 2019;Malcomb et al., 2020;Mathias & Thomas, 2018;Savard et al., 2014;Thomas et al., 2013). Indeed, the carbon isotope signature (δ 13 C tr ) imprinted within annual growth rings of trees records information related to leaf physiology during carbon fixation  and can be used to reconstruct leaf photosynthesis, stomatal conductance, and iWUE (Farquhar et al., 1989;Lavergne et al., 2019;Mathias & Thomas, 2018;Thomas et al., 2013). Moreover, δ 13 C p is useful in contextualizing the nature of iCO 2 pollution impacts to leaf gas exchange (Belmecheri et al., 2021;Saurer et al., 2004) which occur along a spectrum of active to passive responses and may result in a constant (1) leaf internal CO 2 concentration (C i ), (2) ratio of C i to atmospheric CO 2 (C i /C a ), or (3) difference between C a and C i (C a − C i ), respectively (Marshall & Monserud, 1996), with important implications to the carbon cycle.
The eastern United States is one of the most productive, heavily forested areas of the United States (Pan et al., 2013), and has experienced significant climate change over the past century (PRISM Climate Group, 2004). Furthermore, the eastern United States has been exposed to some of the highest levels of acidic air pollution in the country since the rise of industrialization (Driscoll et al., 2001).
Recent work has documented increases in conifer tree iWUE due to iCO 2 and reductions in pollutant emissions after adoption of the Clean Air Act (Mathias & Thomas, 2018;Thomas et al., 2013). The extent to which the leaf physiology of broadleaf deciduous tree species is impacted by acidic pollutants remains fully unresolved, but is nevertheless critical to making informed predictions of forest responses to environmental change across larger geographical areas (Bishop et al., 2015;Levesque et al., 2017;Sullivan et al., 2013).
In this study, we tease apart the combined effects of air pollution (i.e., iCO 2 , SO 2 , NO x ) and climate (i.e., temperature, precipitation) on the long-term iWUE, of two ecologically important K E Y W O R D S acidic air pollution, climate change, dendroecology, intrinsic water use efficiency, stable carbon isotopes, temperate deciduous forest broadleaf tree species in the eastern United States, northern red oak (Quercus rubra L.; hereafter Quru), and tulip poplar (Liriodendron tulipifera L.; hereafter Litu). We focus our analysis on four mature forest stands located in the mid-Atlantic region downwind from a large number of coal-fired power plants along the Ohio River basin. Consistent with experimental observations, one may hypothesize iCO 2 will stimulate iWUE Walker et al., 2019). Alternatively, climate may be a more important driver of physiology, as has been shown in the eastern United States (Levesque et al., 2017;Maxwell et al., 2019). Lastly, given the proximity to coal-fired power plants in the Ohio River basin, high levels of acidic pollution may moderate and potentially overwhelm any positive effects from iCO 2 or changing climate as has been seen in conifers in this region (Mathias & Thomas, 2018;Thomas et al., 2013), as well as other tree species globally (Fernández-Martínez et al., 2017;Gharun et al., 2021). Here, we ask the fol-  (USGS, 1946). PNR is the northernmost study location and is part of the Carnegie Museum of Natural History, while SCBI is the southernmost location and is part of the Forest Global Earth Observatory (ForestGEO) network (Anderson-Teixeira et al., 2015;Bourg et al., 2013), and WS10 and WS13 are reference watersheds located in the Fernow Experimental Forest (Trimble, 1977). All four study locations have remained free from disturbance since the early 1900s and contain Quru and Litu as canopy tree species. Elevation across the four study locations ranges from ~400 to ~750 meters above sea level, while historical (1940-2015) mean peak growing season temperature (May-September) ranges from ca. 17.5 to 20.5°C, and total peak growing season precipitation ranges from ca. 480 to 625 mm (Table S1). Soils at each study location consist of a relatively thin organic layer with litter pri- During the 2016 growing season, ≥20 trees of each species larger than 10 cm diameter at breast height (1.4 m aboveground) and older than ~80 years were cored at each study location using a 5.15 mm increment borer (Haglöf Inc.), except for those from SCBI which were sampled in 2010-2011 (Bourg et al., 2013). Sampled trees were selected across size classes to avoid any potential biases during analyses (Alexander et al., 2018;Nehrbass-Ahles et al., 2014), and cores were taken perpendicular to the slope to avoid compression and expansion wood. Following collection, increment cores were returned to the laboratory where they were allowed to air dry and were processed according to standard dendrochronological techniques (Stokes & Smiley, 1996). Increment cores were scanned at 1000 dpi to a digital image, after which each growth ring was measured to the nearest 0.001 mm and was assigned at the boundary between earlywood and latewood using WinDENDRO (Regent Instruments, Inc.).
Following identification of individual growth rings, chronologies for each respective species within a given study location were then separately crossdated, detrended using a cubic smoothing spline, and statistically confirmed all using dplR (Bunn, 2008;Bunn et al., 2015), after which a calendar year was assigned to each growth ring. All Quru and Litu trees across the four study locations were similarly F I G U R E 1 Study locations in the eastern United States containing the two study species, Quercus rubra and Liriodendron tulipifera. Study locations span across northeastern West Virginia, northwestern Virginia, and southwestern Pennsylvania and are second-growth, mature unmanaged forests and experience a temperate climate (Table S1). The color gradient denotes elevation (1 m) ranging from lowest (yellow) to highest (blue) and were obtained from the USGS. aged, with the mean across all sites being 94 years for both Quru and Litu (Table S2). We excluded the ~30 years of juvenile growth preceding 1940 from our analyses as to focus only on the physiology during the mature growth phase when the trees were canopy dominant and at a relatively stable height (Anderson-Teixeira et al., 2022;Brienen et al., 2017;Marshall & Monserud, 1996;McDowell et al., 2011). The mean expressed population signal (EPS), a measure of how well the sampled trees reflect a chronology developed from an infinite population (Wigley et al., 1984), was 0.91 for Quru and 0.90 for Litu across the four study locations for the period 1940-2015, respectively, with site-specific values listed in Table S2.

| Carbon isotope analysis of tree rings
We prepared samples for C isotope analysis according to standard methods following Mathias and Thomas (2018). After development of absolutely dated chronologies (Helcoski et al., 2019), a subset of five increment cores were randomly selected from both species within each site for tree-ring C isotope signature analysis across the mature phase of the chronology . Increment cores were dissected under a stereo zoom microscope into individual growth rings at the boundary of earlywood and latewood, ground to a fine powder, and 1.0 ± 0.2 mg of ground whole wood from each individual tree ring (i.e., samples were not pooled) was packed into tin capsules (CE Elantech).
Carbon isotope composition of each tree ring was analyzed using a ThermoFisher Delta V+ isotope ratio mass spectrometer (ThermoFisher) interfaced with a Carlo Erba NC 2500 Elemental Analyzer (Carlo Erba) at the University of Maryland Center for Environmental Science (UMCES) Appalachian Lab. Tree-ring stable C isotope composition ( 13 C tr , ‰) was reported in standard delta notation using where R sample is the ratio of 13 C: 12 C in the wood sample and R standard is the ratio of 13 C: 12 C in the standard PeeDee belemnite (PDB) from the PeeDee River Formation in Hemingway, South Carolina. We then calculated leaf-corrected carbon isotope discrimination (∆ 13 C, ‰) from δ 13 C tr using where δ 13 C atm is the carbon isotopic signature of atmospheric CO 2 and d accounts for the carbon isotope fractionation from leaf to wood (Table S3; Quru: −0.21‰, Litu: −0.62‰) (Badeck et al., 2005). From ∆ 13 C, we calculated the leaf intercellular CO 2 concentration (C i , ppm), taking into account photorespiration and post-photosynthetic fractionation effects (Lavergne et al., 2019) using where f is the isotopic fractionation associated with photorespiration (12‰) (Ubierna & Farquhar, 2014), pC a is the partial pressure of atmospheric CO 2 (Pa), a is the fractionation factor (4.4‰) associated with the diffusion of CO 2 through the stomata , b (−28‰) is the fractionation associated with carboxylation by Rubisco (Cernusak & Ubierna, 2022), and Γ* is the CO 2 compensation point in the absence of dark respiration (Pa). We assumed no substantial changes in leaf structure or relative canopy position over the study period and calculations used records of C a and 13 C atm compiled by Belmecheri and Lavergne (2020) in the calculation of C i . Last, we calculated iWUE according to Farquhar et al. (1989) using where 1.6 is the constant ratio for the diffusivity of CO 2 and H 2 O vapor in the air. All carbon isotope calculations described were performed using the R package "isocalcR" (Mathias & Hudiburg, 2022).

| Coupling photosynthesis measurements with ∆ 13 C to simulate A net and g s
To simulate historical A net and g, we measured A-C i curves on Quru (N = 47) and Litu (N = 32) trees during June, July, and August of 2016 across PNR, WS10, and WS13. Measurements were taken between 1000 and 1600 EST on both sun and shade leaves to capture a canopy-integrated A-C i response curve covering the range of CO 2 in our study ( Figure 2). First, branches containing fully intact leaves with no visible sign of damage were collected and placed into floral water picks to maintain the transpiration stream. After a 5-min equilibration period, leaves were placed into the cuvette of a LI-6400XT (Li-Cor, Inc), an open-flow gas exchange system with red and blue LED lights at saturating light (1500 μmol photons m −2 s −1 ) where they were allowed to equilibrate for 5 additional minutes. We then measured rates of photosynthesis at six ambient CO 2 concentrations between 50 and 450 ppm beginning at ambient CO 2 (400 ppm) at the time of the study.
We coupled contemporary measurements of leaf physiology (A-C i curves) with ∆ 13 C-derived C i in the reconstruction of A net and g s throughout the chronology following methods and theory described in Farquhar and Sharkey (1982) and used in Thomas et al. (2013) and Mathias and Thomas (2018). Briefly, for each year we used C i calculated using Equation (3) to compute A net across each respective chronology (N = 5 cores species −1 site −1 ) using each unique relationship between A net and C i measured for each of the two species (Figure 2; Figure S1). As such, each year for a given tree yielded 47 unique estimates of A net and g s for Quru and 32 unique estimates for Litu, which captured not only the range of mean physiological parameters measured for each species, but also incorporates the variability associated with how iCO 2 impacts leaf level processes. For the reconstruction of A net , we assumed the response of A net to C i was conserved throughout the chronology as CO 2 increases, an assumption supported by repeated measurements from FACE experiments in the eastern United States (Springer et al., 2005) and Europe (Gardner et al., 2021). For each estimate of A net , the supply function of photosynthesis was used to calculate stomatal conductance following methods from Farquhar and Sharkey (1982) using

| Statistical analyses
In the present study, we focus our analyses on changes in isotopederived indices of tree physiology, but see Helcoski et al. (2019) for a review of the drivers of interannual Quru and Litu tree growth in this region. We examined the temporal changes in ∆ 13 C, iWUE, A net , and g s across the four study locations for each species using linear mixed effects (LME) models. For each case, we used the lme() function from the "nlme" R package (Pinheiro et al., 2018) with Year as the sole fixed effect and Tree nested within Site as a random effect, including an AR(1,0) autocorrelation structure, and fit via restricted maximum likelihood (REML). We used the "segmented" R package (Muggeo, 2008) to identify statistical breakpoints in chronologies of iWUE.
We determined the underlying contribution of A net and g s to changes in iWUE for each year from 1940 to 2015 for each species within each site individually. To do this, we first fit generalized additive models (GAMs) for each chronology (i.e., A net , g s ) using the gamm() function from the "mgcv" R package (Wood, 2011). We fit the models with an AR(1,0) autocorrelation structure, fit via REML, with the basis dimension for smooths, k, set to 15. We then used finite differences to compute the first derivative, and associated 95% confidence interval (CI), of the relationship between model predicted values for a given year for each variable using the derivatives() function from the "gratia" R package (Simpson, 2021). Trends in any variable (i.e., A net or g s ) were significant when the 95% CI of the first derivative was different from 0.
We used model averaging of LME models and hierarchical partitioning to examine the influence and importance of, and potential interactions among, environmental factors on Quru and Litu iWUE for the period 1940-2015. LME models were constructed using the lme() function from the R package "nlme" (Pinheiro et al., 2018). We used the dredge() function from the R package "MuMIn" (Bartoń, 2017) to examine all possible model combinations (N = 40,043), up to order two, which we ranked by Akaike's corrected information criterion (AIC c ). We then averaged the model parameters across all models in which a given parameter was significant for models with a ∆AIC c < 2 from the model with the lowest AIC c . We mean-centered each continuous predictor variable to avoid issues related to multicollinearity among predictors (Schielzeth, 2010), included tree nested within site as a random factor, and initially fit the models via maximum likelihood to allow comparisons among different fixed effects. We constructed LME models separately for Quru and Litu to examine speciesspecific responses to environmental change. We then used the hier.part() and rand.hp() functions from the R package "hier.part" (Mac Nally & Walsh, 2004;Olea et al., 2010) to identify the sole proportional contribution of each environmental factor to the total variance in iWUE for each species.
Together, this resulted in a 16.1% (0.21% year −1 ) increase in iWUE for Quru trees across the 76-year chronologies and a 24.7% (0.33% F I G U R E 2 Relationship between net CO 2 assimilation rate (A net , μmol CO 2 m −2 s −1 ) and leaf internal CO 2 concentration (C i , ppm) for Quercus rubra (a, purple) and Liriodendron tulipifera (b, green). Measurements were taken on sun and shade leaves during June, July, and August of 2016 to capture the relationship between A net and C i for the historical reconstruction of A net and g s to explain underlying changes in tree intrinsic water use efficiency.
year −1 ) increase in iWUE for Litu trees, with site-level trends listed in  ( Figure 3c). Alternatively, in Litu iWUE reflected a constant leaf C i from 1940 to 1957, after which point a constant C i :C a was maintained ( Figure 3d).

| Reconstruction of historical A net and g s
Using contemporary A-Ci curves to reconstruct historical leaf physiology, we found seasonally integrated Quru A net increased by 0.57 μmol CO 2 m −2 each decade over the time period (p < .001), while Litu A net increased, on average, by 0.54 μmol CO 2 m −2 each decade since 1940 (p < .001) (Figure 4a,b). This translated to a ca.

| Attribution of changes in iWUE to A net and g s
Using finite differences to calculate the annual rate of change in modeled A net and g s based on the GAM fit over time, we found A net increased either alone, or in combination with changes in g s , in 79% of years examined for Quru and 86% of cases for Litu, respectively from 1940 to 2015 ( Figure 5; Figure S3; Table S4). Conversely, we found reduced g s , either alone or in combination with changes in A net , were attributed to 21% and 14% of years examined, in turn, for Quru and Litu trees across all sites during the study period ( Figure 5; Figure S3; Table S4). Increases in A net were more consistent in both species beyond ca. 1970, while changes in g s were more variable throughout the 76 years since 1940 ( Figure S3).

| Environmental drivers of iWUE
Our analysis using LME models revealed positive main effects of VPD and NO x and a negative effect of T mean on Quru iWUE. We

F I G U R E 4
Standardized net photosynthesis (a,b) and stomatal conductance to water (c,d) for the period 1940-2015 for Quercus rubra (a,c; purple) and Liriodendron tulipifera (b,d; green) across all four sites in the Central Appalachian Mountains. For each year A net and g s were reconstructed by coupling ∆ 13 C-derived C i with leaf physiological measurements (see Section 2). The solid trendline in each panel represents the average trend for each variable from 1940 to 2015. Each chronology was standardized by subtracting the value of a given year from the mean value across the entire chronology for each respective tree to remove differences in absolute values of A net and g s across chronologies while maintaining the nature of changes over time. As such, values shown represent either A net or g s relative to the mean value depicted in the lower right panel for each case (e.g., an A net value of 0 for Quru would correspond to 10.25 μmol CO 2 m −2 s −1 ).
found interactions among T mean and NO x with CO 2 on Quru iWUE, such that increasing T mean and NO x diminished the CO 2 effect, while NO x diminished the T mean effect. We also found increasing SO 2 enhanced the negative effect of T mean ( Figure S4; Table S5). Last, we found SO 2 enhanced both the positive effects of VPD and NO x on Quru iWUE. Similarly, we found positive main effects of CO 2 , VPD, and NO x , and negative main effects of T mean and SO 2 on Litu iWUE from 1940 to 2015 ( Figure S4; Table S5). We also found interactions among environmental factors in Litu iWUE, whereby increasing NO x reduced the CO 2 effect, increasing SO 2 reduced the NO x effect, and increasing NO x enhanced the positive effect of VPD ( Figure S4; Table S5).
Hierarchical partitioning reinforced LME model results and revealed together NO x and CO 2 accounted for more than 85% of the explained variance in Quru iWUE, with relatively minor contributions by T mean , SO 2 , and VPD ( Figure 6; Table S6). For Litu, increasing CO 2 contributed to nearly half of the explained variance in iWUE, followed by declining SO 2 and NO x emissions together explaining ca. 40%, with T mean , and VPD playing relatively minor roles ( Figure 6; Table S6).

| DISCUSS ION
In this study we combine tree ring-derived carbon isotope sig- Importantly, we identify a ca. 43%-50% stimulation in A net as the dominant mechanism driving increases in iWUE of both species across ca. 79%-86% of the 76-year chronology, particularly beyond ca. 1970 (Figure 4; Figure S3), adding to the evidence that g s may not always significantly decline as expected under iCO 2 (Guerrieri et al., 2019;. Moreover, we document differences in interspecific leaf gas exchange characteristics that likely resulted in a saturation in Quru iWUE, but a F I G U R E 5 Simultaneous changes in reconstructed seasonally integrated photosynthesis (a) and stomatal conductance to water (b) as underlying drivers of intrinsic water use efficiency for the period 1940-2015 for Quercus rubra (a, purple) and Liriodendron tulipifera (b, green) at four sites in the Central Appalachian Mountains. Each data point represents the first derivative of the generalized additive model fit for each variable (i.e., the instantaneous change) for a given year (see Section 2). Values falling on the abscissa and ordinate represent cases where the annual change in A or g s , in turn, were not different from 0.

F I G U R E 6
Proportion of the contribution of each environmental factor to the variance in observed changes in Quru and Litu intrinsic water use efficiency (iWUE) during the period 1940-2015 across the four study locations. Lighter shading with dashed lines indicates a negative effect of a given variable on iWUE, while darker shading and solid lines indicate a positive effect. The individual contribution of each factor to the total explained variance in iWUE was determined using hierarchical partitioning, which avoids issues related to multicollinearity among predictor variables (see Section 2). The total variance explained, in turn, by the Quru and Litu model is 20% and 31%, and full details can be found in Table S6.
sustained increase in Litu iWUE in recent years ( Figure 3). Finally, we highlight the relatively greater impact of air pollution (iCO 2 , NO x , and SO 2 ) than climate in explaining changes in tree iWUE at our study locations in the mid-Atlantic region ( Figure 6).
We found iCO 2 was one of the most important factors in explaining increases in tree iWUE of Quru and Litu over the last seven decades (Figure 6), supporting previous studies using tree-ring chronologies (Belmecheri et al., 2021;Frank et al., 2015;Keller et al., 2017;Peñuelas et al., 2011;Saurer et al., 2004;Thomas et al., 2013) and consistent with leaf physiology theory (Frank et al., 2015;Lavergne et al., 2019). However, we found several environmental factors interacted with iCO 2 to moderate the integrated response of iWUE in the two study species. For both Quru and Litu, we found increasing NO x stimulated iWUE, but to a greater extent at low [CO 2 ] ( Figure S4), likely through increased N deposition resulting in higher plant available N stimulating A net (Gharun et al., 2021;Jennings et al., 2016;Leonardi et al., 2012). We additionally found iCO 2 stimulated Quru iWUE less when it was warm ( Figure S4; Table S5), pointing to the interactive effects between atmospheric pollutants and climate on tree physiology , but also suggesting leaf respiration may increase proportionally more than photosynthesis. This may lead to Quru trees at relatively cooler, higher latitudes experiencing larger increases in future iWUE as CO 2 continues to rise. We highlight positive interactions between VPD and SO 2 and NO x , in turn, on Quru and Litu iWUE. While increasing VPD led to increasing iWUE in both species, likely through reductions in g s (Grossiord et al., 2020;Novick et al., 2016), this effect is likely compounded by SO 2 on Quru (Choi et al., 2014;Mathias & Thomas, 2018;Ooi et al., 2018), whereas NO x likely stimulated A net in Litu (Jennings et al., 2016;Leonardi et al., 2012;Mathias & Thomas, 2018 (Belmecheri et al., 2021), and at broader scales globally (Guerrieri et al., 2019;Levesque et al., 2017;. Thus, while our study shows increasing CO 2 is a strong driver of increased water use efficiency in trees, it is clear this is context dependent, and modulated by other environmental factors. Evidence for a spectrum of physiological behavior and responses at the leaf and tree level have been observed across a range of environmental drivers (Ainsworth & Rogers, 2007;Belmecheri et al., 2021;Brzostek et al., 2014;Horn et al., 2018;Lin et al., 2015;Saurer et al., 2004;Thomas et al., 2010).
Our results here indicate divergent leaf gas exchange strategies in provides evidence for a broader physiological phenomenon across unique plant functional types, but also lends support for a range physiological control under in tree species experiencing similar environmental conditions (Bryant et al., 2022;Brzostek et al., 2014;Roman et al., 2015), indicating the importance of species-level responses to environmental change.
Our analysis of changes in modeled A net and g s reveals that enhanced A net overwhelmed reductions in g s in driving Quru and Litu iWUE throughout the majority of the last seven decades (Figures 4 and 5; Figure S3; Table S4). In fact, cases where reduced g s contributed to increasing iWUE occurred almost exclusively before ca 1970 ( Figure S3) when mean precipitation was lower ( Figure S5).
Importantly, our findings are robust against changes in the nature the leaf A-C i response to iCO 2 ( Figure S6). While iCO 2 independently increases A net and decreases g s , thereby stimulating iWUE (Ainsworth & Long, 2004;Ainsworth & Rogers, 2007), local site conditions are known to modify this response (Belmecheri et al., 2021;Fernándezde-Uña et al., 2016). Indeed, observations of increasing water availability have been linked to similar patterns in g s in trees throughout North America despite simultaneous increases in iWUE (Guerrieri et al., 2019;Levesque et al., 2017), again underscoring the importance of relatively greater increases in stimulated A net . The lack of a strong reduction in g s over time for Quru and Litu trees in our study is in agreement with recent analyses using a dual carbon-oxygen stable isotope technique from trees around the globe (Guerrieri et al., 2019;, despite suggested caveats using the dual isotope method (Guerrieri et al., 2022;Lin et al., 2022), and is likely a product of generally mesic conditions at our study locations. Finally, our findings of increasing A net stimulating iWUE are consistent with observations of increasing carbon uptake globally (Campbell et al., 2017;Cernusak et al., 2019;Haverd et al., 2020), and highlight the nuanced, context-dependent response of leaf physiology to environmental change.
A careful consideration of the complexity behind annual tree growth and physiology is required to attribute the role of any given environmental factor on long-term trends in tree iWUE. While we highlight some of the environmental drivers and their interactions that affect iWUE in this study, there are factors such as competition (González de Andrés et al., 2018), potential legacy effects of SO 2 and NO x (Likens et al., 1996), and light environment ) that we are unable to consider and may contribute some uncertainty to our analysis. Moreover, while iWUE increases with tree height (McDowell et al., 2011), light environment can exert stronger controls over iWUE , complicating this response. Nevertheless, our study using mature, canopy-dominant trees contextualizes the response of Quru and Litu tree iWUE to changes in multiple environmental factors and is in direct agreement with assessments from regional (Mathias & Thomas, 2018) and global analyses (Belmecheri et al., 2021;Frank et al., 2015;Guerrieri et al., 2019;. The extent to which tree iWUE will continue to increase remains an open question. Our study clearly highlights the importance of local environmental factors beyond climate, many of which interact ( Figure S4; Table S5), that have large influences over iWUE ( Figure 6).
Indeed, there is growing evidence for stimulated A net persisting under iCO 2 (Ainsworth & Long, 2004;Ainsworth & Rogers, 2007;Gardner et al., 2021;Springer et al., 2005) which would be reflected in iWUE, so long as nutrient demand is met (Luo et al., 2004) and temperature increases are moderate . However, our data from Quru and Litu trees in the eastern United States also pinpoint the impact of acidic air pollution (i.e., SO 2 ) on leaf physiology. For example, it is likely the increase in g s observed at some sites in recent years ( Figure 4; Figure S3) is linked to reductions in SO 2 emissions (Borer et al., 2005;Ooi et al., 2018), which reinforces the importance of air quality on forest function (Engel et al., 2016;Mathias & Thomas, 2018). Future work will need to continue to directly assess the underlying mechanisms impacted by air pollution outlined here. Last, although we did not find precipitation to impact interannual variability of iWUE in our study, it is clear changes in water availability impacts tree iWUE-potentially saturating Quru iWUEand should not be overlooked (Belmecheri et al., 2021;Levesque et al., 2017).
Accurately resolving the contribution of concurrently changing environmental factors on tree physiology is not trivial. While in this study we document moderate increases in the iWUE of Quercus rubra and Liriodendron tulipifera trees located in the central Appalachian Mountains, how other forest ecosystems will respond to future climate change still requires attention. Indeed, we found iCO 2 to be an important factor in the proximate increases in iWUE for both species, but there is evidence for an attenuation of this response . Moreover, it has become clear considering the impact of air pollutants (i.e., NO x , SO 2 ) is important in the interpretation of tree iWUE, especially over decadal to century time scales (Gharun et al., 2021;Jennings et al., 2016;Leonardi et al., 2012;Mathias & Thomas, 2018), and must be taken into account moving forward, particularly in highly industrialized areas. The interactions among environmental factors on leaf physiology uncovered here reinforce previous work at larger scales across a range of tree species and site environmental conditions (Fernández-Martínez et al., 2017).
While it is clear the water cost of carbon uptake has declined over the last few decades, and in our study locations primarily due to increasing A net , the extent to which these findings are observed more broadly is critical to resolve as changes in climate and atmospheric composition continue to accelerate.

CO N FLI C T O F I NTE R E S T S TATE M E NT
The authors declare no conflict of interest.

DATA AVA I L A B I L I T Y S TAT E M E N T
All data and accompanying code are publicly available and found on