Temperature-sensitive biochemical 18 O-fractionation and humidity-dependent attenuation factor are needed to predict δ 18 O of cellulose from leaf water in a grassland ecosystem

of leaf cellulose ( δ 18 O cellulose ) in a drought-prone, temperate grassland ecosystem. A new allocation-and-growth model was designed and added to an 18 O-enabled soil – vegetation – atmosphere transfer model (MuSICA) to predict seasonal (April – October) and multi-annual (2007 – 2012) variation of δ 18 O cellulose and 18 O-enrichment of leaf cellulose ( Δ 18 O cellulose ) based on the Barbour – Farquhar model. (cid:1) Modelled δ 18 O cellulose agreed best with observations when integrated over c. 400 growing-degree-days, similar to the average leaf lifespan observed at the site. Over the integration time, air temperature ranged from 7 to 22 ° C and midday relative humidity from 47 to 73%. Model agreement with observations of δ 18 O cellulose ( R 2 = 0.57) and Δ 18 O cellulose ( R 2 = 0.74), and their negative relationship with canopy conductance, was improved signiﬁcantly when both the biochemical 18 O-fractionation between water and substrate for cellulose synthesis ( ϵ bio , range 26 – 30‰) was temperature-sensitive, as previously reported for aquatic plants and heterotrophically grown wheat seedlings, and the proportion of oxygen in cellulose reﬂecting leaf water 18 O-enrichment (1 – p ex p x , range 0.23 – 0.63) was dependent on air relative humidity, as observed in independent controlled experiments with grasses. (cid:1) Understanding physiological information in δ 18 O cellulose requires quantitative knowledge of climatic effects on p ex p x and ϵ bio .


Introduction
The oxygen isotope composition of plant cellulose (δ 18 O cellulose ) and its enrichment above source water (Δ 18 O cellulose ) are thought to record environmental and physiological information of great interest to a range of scientific disciplines, including functional plant ecology and climate change biology (e.g. Barbour, 2007;Battipaglia et al., 2013;Gessler et al., 2014). In particular, δ 18 O cellulose or even more so Δ 18 O cellulose has been discussed as an integrated proxy of past stomatal conductance (e.g. Farquhar et al., 1998;Scheidegger et al., 2000;Barbour et al., 2000b) that can provide information about environmental or climate change effects on the processes regulating the water-use efficiency of C 3 plants. However, empirical evidence for a direct link between Δ 18 O cellulose and stomatal conductance is currently incomplete, which prevents the use of δ 18 O cellulose time series to interpret changes in plant water use efficiency. Fundamentally, all of the oxygen in cellulose is derived from water (DeNiro & Epstein, 1979;Liu et al., 2016). Oxygen exchange with water can occur in multiple steps in the photosynthetic carbon cycle up to the formation of sucrose (or other transport sugars), and during metabolism of these sugars in sink tissues, when cellulose is being formed (Hill et al., 1995;Farquhar et al., 1998;Barbour et al., 2005). Therefore, δ 18 O cellulose depends on the oxygen isotope composition of water in sink tissues, often approximated by the δ 18 O of plant source water taken up by the roots (δ 18 O source ), and on the δ 18 O of leaf lamina water (δ 18 O leaf ), respectively. During photosynthesis, leaf transpiration enriches isotopically leaf lamina water above source water (Δ 18 O leaf ≈ δ 18 O leaf − δ 18 O source ). δ 18 O source can vary dynamically and is primarily controlled by the δ 18 O of meteoric inputs (δ 18 O rain ; Dansgaard, 1964;Bowen & Wilkinson, 2002), their mixing with soil water and the intensity of soil evaporation (Barnes & Allison, 1988). In addition, δ 18 O source is determined by the depth distribution of Ó 2020 The Authors New Phytologist Ó 2020 New Phytologist Trust New Phytologist (2020) 1 www.newphytologist.com This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited. roots and their specific uptake intensities as well as by the effect of transpiration and root water uptake on soil water emptying and refilling dynamics (Brinkmann et al., 2018;Hirl et al., 2019).
Arguably, the oxygen isotope enrichment of cellulose above source water (Δ 18 O cellulose ≈ δ 18 O cellulose − δ 18 O source ) has the most important physiological interest because of its link with Δ 18 O leaf and stomatal conductance. According to the Barbour-Farquhar model (Barbour & Farquhar, 2000), The product p ex p x reflects the proportion of oxygen in cellulose derived from source water, where p ex is the proportion of exchangeable oxygen in the intermediates formed during cellulose synthesis from sucrose exported from leaves, and p x is the proportion of source water at the site of cellulose formation. In other words, 1 − p ex p x represents the proportional effect of leaf lamina water 18 O-enrichment on the 18 O-enrichment of cellulose. ϵ bio is the average biochemical fractionation between the organic substrate for cellulose synthesis and water (Sternberg et al., 1986;Barbour, 2007).
Equation 1 relies on the assumption that Δ 18 O leaf represents the 18 O-enrichment of water in isotopic equilibrium with leaf sucrose, as supported by two laboratory studies on castor bean (Barbour et al., 2000a;Cernusak et al., 2003). In those studies, conducted in climate-controlled, steady-state conditions, the 18 O-enrichment of sucrose above source water (Δ 18 O sucrose ) was well approximated by Δ 18 O leaf + ϵ bio , with ϵ bio set to 27‰. However, a recent study by Lehmann et al. (2017) on two C 3 grass species indicated that newly formed sucrose might not always be in isotopic equilibrium with average leaf lamina water. Moreover, in submerged aquatic plants (where Δ 18 O leaf is c. 0) ϵ bio was found to be inversely related to growth temperature, with a particularly strong temperature dependence below about 20°C (Sternberg & Ellsworth, 2011). Virtually the same temperaturesensitivity of ϵ bio was found in wheat seedlings during heterotrophic growth (Sternberg & Ellsworth, 2011). Although there is a likelihood that this result applies to autotrophic terrestrial plants, it is difficult to prove because it cannot be verified, as p ex p x and ϵ bio cannot be measured simultaneously.
Probably, the least contentious parameter in Eqn 1 is p x , at least in nontranspiring tissue such as the leaf-growth-and-differentiation zones of grasses or the developing cells of tree trunks, where p x has been shown to stay close to 1 Liu et al., 2017a). Sucrose is the most common carbohydrate transported within plants (Lalonde et al., 2003) and the main substrate from which UDP-glucose, the immediate precursor of cellulose synthesis, is formed in heterotrophic tissue (Verbancic et al, 2018), such as the leaf growth zone of grasses (Baca Cabrera et al., 2020). Based on theoretical considerations of oxygen exchange between the metabolites of sucrose and water during cellulose formation and on observational data, the value of p ex is often assumed to be around 0.4 (Sternberg et al., 1986;Roden & Ehleringer, 1999;Barbour & Farquhar, 2000;Cernusak et al., 2005), but considerable variation in p ex has been suggested from other observational studies (Gessler et al., 2009;Song et al., 2014;Cheesman & Cernusak, 2017). As it cannot be measured in vivo, estimates of p ex must be obtained as a fitted parameter, whilst assuming values for p x and ϵ bio . Such estimations are sensitive to errors, including any effect of isotopic disequilibrium between Δ 18 O sucrose and Δ 18 O leaf (Lehmann et al., 2017). Higher isotopic disequilibria between leaf lamina water and sucrose have also been found when air relative humidity (RH) was lower (Lehmann et al., 2017). When interpreted with Eqn 1, this latter result suggests that p ex p x may increase with increasing RH. A positive relationship between RH and p ex p x was also suggested by the relationship between Δ 18 O cellulose and Δ 18 O leaf in the studies of Liu et al. (2016) and Helliker & Ehleringer (2002a) for a range of C 3 and C 4 grasses, perhaps pointing to isotopic disequilibria between Δ 18 O sucrose and Δ 18 O leaf also in these cases. Clearly, there remain important knowledge gaps on the effect of environmental conditions on ϵ bio and p ex p x and, hence, their implication for physiological interpretation of δ 18 O cellulose or Δ 18 O cellulose from field studies.
The temporal integration of all processes involved in the making of leaf cellulose (i.e. assimilation, mobilization of stored substrate, allocation of substrate to growth, and cellulose synthesis) represents an additional challenge when interpreting δ 18 O cellulose or Δ 18 O cellulose (Hemming et al., 2001;Damesin & Lelarge, 2003;Gessler et al., 2009Gessler et al., , 2014Royles et al., 2013;Liu et al., 2017b). In intensively managed grassland, shoot biomass is mostly vegetative and consists of short-lived leaves (Lemaire et al., 2000). Leaf production during the growing season is continuous, with new leaf production occurring simultaneously with the senescence of older leaves . Accordingly, the mean live leaf age (measured in growing-degree-days (GDD)) changes relatively little during the course of the vegetation period (Lemaire et al., 2000;Schleip et al., 2013). In a controlled environment study with Cleistogenes squarrosa (a perennial C 4 grass), Liu et al. (2017b) found a close linear relationship between the fraction of remaining leaf elongation and the fraction of oxygen in cellulose assimilated following a change of RH in the growth environment, as inferred from 18 O-abundance measurements. Similar results were previously obtained by Helliker & Ehleringer (2002b), highlighting the potential of grass leaves as recorders of environmental signals in δ 18 O cellulose . In light of this, allocationand-growth models, coupled with isotope-enabled process-based models of soil-vegetation-atmosphere CO 2 and H 2 O exchange, may be very useful tools for analysing the mechanisms and dynamics controlling Δ 18 O cellulose and δ 18 O cellulose in natural environments. Such models have already been applied to treerings (Roden et al., 2000;Barbour et al., 2002;Ogée et al., 2009;Keel et al., 2016;Lavergne et al., 2017), but are presently unavailable for grassland.
The present study explores the effects of environmental drivers on the parameters of the Barbour & Farquhar (2000) (Hirl et al., 2019). Here, we first validated the combined MuSICA and allocationand-growth model by testing its ability to predict the labelling kinetics of autotrophic respiration observed by Gamnitzer et al. (2009), observed root : shoot C allocation (Schleip, 2013), and the observed δ 18 O cellulose of leaves at the study site, using integration times consistent with leaf growth dynamics observations. We then applied the model to explore the three questions developed above.

Experimental site and sampling
The study was conducted inside pasture paddock no. 8 of Grünschwaige Grassland Research Station near Freising, Germany (for details on site, vegetation and grazing management see Schnyder et al., 2006). Average air temperature in the study years 2007-2012 was 9.3°C, and mean annual precipitation was 753 mm (recorded at the Munich airport meteorological station). The mineral topsoil has a low water-holding capacity (66 mm plant-available field capacity) causing frequent and prolonged drought periods (Hirl et al., 2019). The pasture was continuously grazed by Limousin suckler cows during the growing seasons (from mid-April to beginning of November). Animal stocking density was adjusted periodically to balance grass production and consumption by the cattle (Lemaire et al., 2009), so that mean sward height was maintained at about 7 cm.
Two replicate leaf samples were collected at around midday (between 11:00 and 16:00 h Central European Summer Time) at approximately fortnightly intervals during the vegetation periods 2007-2012 (Hirl et al., 2019). Each sample consisted of a mixed-species collection of the codominant species: four C 3 grasses (Lolium perenne L., Poa pratensis L., Phleum pratense L., Dactylis glomerata L.), one rosette dicot (Taraxacum officinale F.H. Wigg.) and one legume (Trifolium repens L.). The sample included only the green, nonsenescing, fully expanded leaf blades, as well as the exposed part of the growing leaf (cf. Fig. 1 of Liu et al., 2017a) from 16 vegetative tillers of L. perenne, P. pratensis and P. pratense and two vegetative tillers of D. glomerata, as well as one half of a leaf blade of T. officinale (with the midvein removed) and two trifoliate leaves of T. repens. Water in these leaf samples, along with source, soil, atmospheric humidity and rainwater samples had previously been analysed for δ 18 O, and presented and discussed in Hirl et al. (2019).
Maintenance respiration (F maint,resp , R pool ) Assimilation (F assim , R assim ) Allocation to growth (F growth , R pool ) Shoot growth and respiration (F shoot,growth + F shoot,growth,resp , R pool ) Root growth and respiration (F root,growth + F root,growth,resp , R pool ) Fig. 1 Scheme of the allocation-and-growth model for predicting carbon fluxes and cellulose isotope compositions in grassland. Photosynthetic assimilation (F assim ) supplies substrate with specific isotopic composition (R assim ) to a well-mixed metabolic pool (W pool , with isotope ratio R pool ). Maintenance respiration (F maint,resp ) and allocation of substrate to growth (F growth ) are both effected from the metabolic pool and thus carry the isotope signal of the pool (R pool ). Partitioning between shoot and root growth (F shoot,growth and F root,growth , which involve growth respiration F growth, resp = F shoot,growth,resp + F root,growth,resp ) is governed by xylem water potential (Ψ xylem ). The shoot and root structural biomass pools are represented as layered (or stacked, nonmixing) pools. The integration time (d) represents the maximum age of structural leaf biomass in a sample. The size (or 'thickness') and isotopic composition of a given daily layer (or stack) in the shoot or root is determined by shoot and root growth rates and by the isotope ratio of the pool on that day (for details see Materials and Methods section and Supporting Information Methods S2).

Cellulose extraction and isotopic analysis
Following the procedure of Brendel et al. (2000) as modified by Gaudinski et al. (2005), α-cellulose was extracted from a subsample (50 or 25 mg) of ground plant material (for details see Liu et al., 2017b). After redrying of the cellulose at 80°C for 24 h, 0.7 mg aliquots were weighed into silver cups (size: 3.3 × 5 mm, IVA Analysentechnik e.K., Meerbusch, Germany) and stored above Silica Gel orange (2-5 mm, ThoMar OHG, Lütau, Germany) in exsiccator vessels. Samples were pyrolysed at 1400°C in a pyrolysis oven (HTO, HEKAtech, Wegberg, Germany), equipped with a helium-flushed zero blank autosampler (Costech Analytical technologies, Valencia, CA, USA), interfaced (ConFlo III, Finnigan MAT, Bremen, Germany) to a continuous-flow isotope ratio mass spectrometer (Delta Plus, Finnigan MAT). A solid internal laboratory standard (cotton powder) was included after every third or fourth sample and used for V-SMOW scaling and instrument drift correction. Every sample was analysed in duplicate. All samples and the laboratory standard were measured against a laboratory working reference carbon monoxide gas, which had previously been calibrated against a secondary isotope standard (IAEA-601; accuracy of calibration AE 0.25‰ SD). The precision for the laboratory standard was < 0.3‰ (SD for repeated measurements). Oxygen isotope composition is expressed in per mil (‰) as δ 18 O = (R sample /R standard − 1), with R sample the 18 O : 16 O ratio of the sample and R standard that in the V-SMOW standard.

Model
The integral model was composed of the 18 O-enabled soil-vegetation-atmosphere transfer model MuSICA (Ogée et al., 2003(Ogée et al., , 2009Wingate et al., 2010;Gangi et al., 2015), parameterized as in Hirl et al. (2019), and a new allocation-and-growth model ( Fig. 1). MuSICA is a multilayer multileaf model that simulates CO 2 , water and energy redistribution and isotopic exchange processes in an ecosystem based on current mechanistic understanding (see Supporting Information Methods S1). MuSICA requires half-hourly climate data at a reference level above the vegetation (0.5 m here) as well as the isotopic composition of rainfall, CO 2 and water vapour at the same reference level. Meteorological variables were obtained from the meteorological station at Munich airport at about 3 km from the study site (wind speed, precipitation, air temperature, RH, air pressure), from two other meteorological stations at 10 and 12 km distance (radiation), and from an eddy flux station installed at the experimental site (CO 2 concentration) (Hirl et al., 2019). For the precipitation and water vapour isotopic input data (δ 18 O rain and δ 18 O vapour ), we used data collected at the site whenever available, and gap-filled the data with offset-corrected IsoGSM data (Yoshimura et al., 2011) as detailed in Hirl et al. (2019). δ 18 O and δ 13 C data of atmospheric CO 2 (δ 18 O CO2 and δ 13 C CO2 ) were obtained from NOAA/CMDL latitudinal products (J. Miller, pers. comm.). MuSICA parameterization included soil and vegetation properties that described the pasture system in terms of structural and hydrological characteristics, leaf gas exchange, and root distribution and hydraulics. MuSICA, in its 'standard parameterization' (Hirl et al., 2019), performed well in predicting δ 18 O leaf , δ 18 O source (termed δ 18 O stem in Hirl et al., 2019) and δ 18 O soil at different soil depths throughout seven growing seasons (2006--2012). The model also predicted the dynamics of transpiration, canopy conductance (g canopy ), root water uptake and plant-available soil water (PAW) (Hirl et al., 2019). MuSICA parameters and parameters of the allocation-and-growth model were based on singular measurements (e.g. parameters of the soil water retention curve) or on measurements performed at intervals during the study period (e.g. canopy height, leaf area index). In the latter case, average observed values were used for the standard simulation, in agreement with the constant sward-state objective that ruled pasture management (see above). Missing parameters were taken from the literature (e.g. photosynthetic parameters) (see table S1 in Hirl et al., 2019).
A detailed description of the new allocation-and-growth model ( Fig. 1, Table 1) is provided in Methods S2. Briefly, the model includes three compartments: one well-mixed metabolic pool (W pool ) and two structural (i.e. nonmetabolic) compartments representing the aboveground shoots (W shoot ) and belowground roots (W root ), similar to Ostler et al. (2016). Current assimilates replenish W pool with rate (F assim , gross primary production) and 13 C and 18 O signatures (δ 13 C assim and δ 18 O assim , see below) predicted by MuSICA (Fig. 1). Maintenance respiration (F maint,resp ) is directly supported by W pool , while growth respiration (F growth, resp ) is a constant fraction of growth and included in the flux of substrate allocated to the growth of shoots and roots (F growth ). F maint,resp is an exponential function of soil surface temperature and is assumed to be proportional to total plant carbon mass and negatively related to W pool when the latter decreases below a target value (Thornley & Cannell, 2000). Growth (and growth respiration) occurs only if the metabolic pool is replenished to its target value. Substrate for growth is allocated between shoot (F shoot,growth ) and root (F root,growth ) depending on root xylem water potential (Ψ xylem ), to account for the drought-sensitivity of leaf growth (e.g. Durand et al., 1995), and the fractional allocation to the shoot (f shoot,growth ) is greater in spring than summer/ autumn. New assimilates carry the isotopic signal of current assimilation, while the substrate supplied to growth and respiration carries the isotopic signal of the metabolic pool, assuming no isotopic fractionation during respiration or growth.
The 18 O : 16 O ratio of cellulose in a sample (R cellulose ) was calculated similar to Ogée et al. (2009), assuming that cellulose synthesis was proportional to shoot growth. The integration time in days was computed for a set thermal time by summing up daily mean soil temperature above a base temperature of 4°C . The biochemical fractionation (ϵ bio ) between water and substrate used in cellulose synthesis was calculated from daily mean air temperature using the empirical function found by Sternberg & Ellsworth (2011, their fig. 2). The p ex p x was calculated from midday RH of air: p ex p x = 0.016 RH − 0.393 (Fig.  2), derived by regressing the differences between Δ 18 O cellulose predicted with a constant, nonbiased p ex p x and observed Δ 18 O cellulose against RH (Fig. 3a) (Fig. 2). Alternatively, we also calculated the isotope composition of cellulose by keeping either p ex p x (0.556) or ϵ bio (27‰) or both p ex p x and ϵ bio constant.

Statistical and sensitivity analysis
Model performance was evaluated by calculating the mean bias error (MBE = P À O, with P the mean predicted value and O the mean observed value) between observed and predicted δ 18 O cellulose (or Δ 18 O cellulose ), the mean absolute error , with P j the predicted and O j the observed value of sample j, and n the number of samples (Willmott & Matsuura, 2005), and R 2 values. Linear regression and correlation analysis was performed to investigate the relationship between δ 18 O cellulose (or Δ 18 O cellulose ) and environmental parameters. All analyses were conducted in R, v.3.4.2 (R Core Team, 2017).
A sensitivity analysis was performed (1) to quantify and disentangle the effect of meteorological variables on Δ 18 O cellulose and δ 18 O cellulose , and (2) to evaluate the responsiveness of δ 18 O cellulose and Δ 18 O cellulose to the parameters of the allocation-and-growth model. Two sensitivity runs were performed for each parameter or meteorological variable, using a minimum and a maximum value, based on the range of observed or expected values (see Table 1). In each sensitivity run, one parameter (or meteorological variable) at a time was changed while all other parameters (or variables) were held the same as in the standard simulation. Regarding (1) above, the incoming short-wave radiation and wind speed, two variables that are expected to affect the leaf energy budget, were halved and doubled for each time step. Regarding (2) above, the applied parameter ranges of the allocation-and-growth model were first derived from measurements at the experimental site conducted during the study period; if no measurements were available for a parameter, relevant ranges were taken from the literature (see Table 1). Quantification of parameter effects (sensitivities) followed the procedure outlined in Hirl et al. (2019). First, systematic effects were quantified based on the 'mean sensitivity', computed as the mean difference between the sensitivity and the standard run in a sensitivity run and δ ref,j is that in the standard run. Second, the variability of the parameter effect was depicted based on the standard deviation of the sensitivity, calculated from the differences between δ sens,j and δ ref,j . To disentangle an eventual contribution of a temperature-sensitive ϵ bio from the effects of other processes on the temperature sensitivity of Δ 18 O cellulose , we also performed a range of sensitivity analyses in which air temperature was decreased (or increased) by 1, 3 or 5°C for all half-hourly values.

Metabolic pool turnover, carbon allocation and integration time
The turnover of the metabolic pool was validated by comparing the predicted 13 C-labelling of W pool (f new , the fraction of labelled C in W pool ), forced by a step-change of δ 13 C CO2 in the MuSICA input data, with the labelling kinetics of total autotrophic respiration observed at the same site in May 2007 (Gamnitzer et al., 2009;see their fig. 5). The predicted f new in W pool matched closely the observed fraction of labelled C in total autotrophic respiration (slope = 0.99; intercept = 0.06; R 2 = 0.98; P < 0.001; Fig. 4a).
Model predictions for the fraction of carbon allocated to shoot growth (f shoot,growth ) compared well with observations of f shoot,growth made by Schleip (2013)  To constrain the integration time of the 18 O-signal in the leaf cellulose samples, we compared the R 2 between observed and predicted δ 18 O cellulose (Fig. 4c) for predictions of δ 18 O cellulose based on integration times of 50-600 GDD. Those calculations were made with the model in its standard parameterization (Table 1) and assumed leaf growth started on the 15 March every year. The R 2 between observed and predicted δ 18 O cellulose increased up to an integration time of 400 GDD (R 2 = 0.57) and decreased beyond that. This integration time of 400 GDD translated to time spans varying between 22 and 59 d in different periods, depending on thermal conditions (Fig. 4d) and was subsequently used as the standard integration time for all predictions of δ 18 O cellulose and Δ 18 O cellulose .
Observed variation of δ 18 O cellulose and Δ 18 O cellulose δ 18 O cellulose displayed dynamic variation within and between years, with pronounced increases and decreases in some years (2007,2010) and lower variability in others (2008Fig. 5). The total range of observed δ 18 O cellulose was 5.1‰, and a clear seasonal pattern was not evident (right panel, Fig. 5a). Annual mean δ 18 O cellulose did not differ significantly between the individual years, except for 2010 when the mean was 1.5‰ lower than the average of the other years. The lower δ 18 O cellulose in 2010 was linked to a more negative δ 18 O rain and lower leaf water 18 O-isotopic enrichment in that year (cf. Hirl et al., 2019).
The range of observed Δ 18 O cellulose was 50% greater than that of δ 18 O cellulose (compare Fig. 5a and b). Unlike δ 18 O cellulose , average Δ 18 O cellulose declined distinctively during the growing season at a rate of 0.83‰ per month (R 2 = 0.90; P < 0.01). This was related to an increasing trend of δ 18 O source over the growing season (Hirl et al., 2019). Overall, the observed variation Δ 18 O cellulose showed statistically significant relationships with most meteorological variables averaged over the respective integration time (Table 2; Fig. 6). The correlation of Δ 18 O cellulose with midday RH (r = −0.69) was stronger than with daily mean air temperature (r = −0.43) and midday air temperature (r = −0.37). Δ 18 O cellulose was only weakly related to PAW in the rooting zone (r = −0.24; P = 0.05) and there were no statistically significant relationships with vapour pressure deficit (VPD) or annual precipitation (both P > 0.05). Correlation analysis also indicated connections between Δ 18 O cellulose and cumulative short-wave radiation (r = 0.81) and wind speed (r = 0.56). Yet, a sensitivity analysis in which the radiation input was halved or doubled demonstrated that cumulative short-wave radiation itself had only a very small effect on Δ 18 O cellulose (Fig. S1), indicating that the strong correlation with Δ 18 O cellulose was indirect and arose primarily from the relationship of short-wave radiation with midday RH (r = −0.63) and temperature (r = 0.57) (Table S1). Sensitivity analysis also indicated that wind speed itself had a negligible effect on Δ 18 O cellulose (Fig. S1).

Prediction of δ 18 O cellulose and Δ 18 O cellulose
The standard model, with RH-sensitive p ex p x and temperaturesensitive ϵ bio , was able to reproduce with great accuracy the  Red points indicate f shoot,growth as estimated by 13 C-labelling for two 2-wk periods in May and September 2007 (Schleip, 2013), and blue points denote the model predictions for the same 2-wk periods. (c) R 2 values for the comparison between predicted and observed δ 18 O cellulose . Model predictions were obtained using integration times between 50 and 600 growing-degree-days (GDD). The vertical dashed line indicates the average leaf lifespan observed by Schleip et al. (2013) for four co-dominant species of the same pasture investigated in this study. (d) Integration times (Δt) calculated from soil temperature and assuming a constant thermal age of 400 GDD with a base temperature of 4°C (see Eqn S8 in Supporting Information Methods S2). The blue points mark the days on which sampling occurred. It was assumed that leaf growth did not start before 15 March.
Ó 2020 The Authors New Phytologist Ó 2020 New Phytologist Trust New Phytologist (2020) www.newphytologist.com observed seasonal dynamics of δ 18 O cellulose (R 2 = 0.57) and Δ 18 O cellulose (R 2 = 0.74) over the entire 6-yr-long study period (Figs 5,7d). The model error and bias were small (Table 3), as the model was also able to capture the short-term variations. By contrast, predictions of Δ 18 O cellulose made with constant p ex p x and/or constant ϵ bio degraded significantly the goodness-of-fit with observed Δ 18 O cellulose (Fig. 7a-c).
Variation of midday RH at the study site, averaged over the integration time, was significant (47-73%) and implied a large range of modelled p ex p x (0.37-0.77; average: 0.59), with a generally increasing trend during the growing season (Fig. 5c). Meanwhile, daily mean air temperature varied considerably at the scale of integration times (between 7 and 22°C), leading to pronounced variation in ϵ bio (from 26 to 30‰; Fig. 5d). ϵ bio generally decreased from spring to summer and increased towards late autumn. The ability of the model to reproduce the observed relationships between Δ 18 O cellulose (or δ 18 O cellulose ) and meteorological variables or PAW also depended on the inclusion of an RHsensitive p ex p x and a temperature-sensitive ϵ bio (Table 2). In particular, the model was unable to reproduce the observed relationship between Δ 18 O cellulose and air temperature when a constant ϵ bio of 27‰ was used instead of a temperature-dependent one ( Fig. S3d-f).

Relationship between canopy conductance and Δ 18 O cellulose or δ 18 O cellulose
Observed Δ 18 O cellulose was negatively related to midday canopy conductance predicted by the MuSICA model (R 2 = 0.26; P < 0.001; Table 2; Fig. 8), with a 100 mmol m −2 s −1 increase in g canopy connected to a 1.4‰ decrease of Δ 18 O cellulose . Again, the relationship between midday canopy conductance and observed or predicted Δ 18 O cellulose was the same if the prediction was based on RH-sensitive p ex p x and temperature-sensitive ϵ bio , but not if constants for p ex p x or ϵ bio were used.
The relationship between observed δ 18 O cellulose and modelled midday canopy conductance was also negative (P < 0.05), although the slope and R 2 for that relationship were much smaller than those for Δ 18 O cellulose . The relationship between predicted δ 18 O cellulose and modelled canopy conductance was similar to the observations when predictions were based on RHsensitive p ex p x and temperature-sensitive ϵ bio , and became

Sensitivity to model parameters
Model predictions of δ 18 O cellulose or Δ 18 O cellulose were relatively insensitive to alterations of the parameters of the allocation-andgrowth model (Table 1; Fig. S1), except for integration time (Fig. 4c). That is, both the mean sensitivity (i.e. the mean difference between the sensitivity run and standard run for predictions of δ 18 O cellulose (or Δ 18 O cellulose )) as well as the standard deviation of the sensitivity was smaller than 0.2‰ for each parameter value (Fig. S1).
In the MuSICA model, a 2.5-fold increase of m gs (see Table 1) the slope of the Ball-Woodrow-Berry stomatal conductance model (Ball et al., 1987) caused an average reduction of 0.84‰ of the predicted δ 18 O cellulose and Δ 18 O cellulose . For a particular situation (time), the size of that effect depended on the conditions prevailing during the respective integration time. This conditionality was reflected in the relatively high standard deviations between individual sensitivity and standard model runs for m gs (Fig. S1). Conversely, a reduction of m gs by one-third caused an average increase of 0.3‰, with the increase again dependent on prevailing conditions, as above. On the other hand, increasing V cmax and J max by 2.2-fold or reducing it by two-thirds had only small effects on the prediction of δ 18 O cellulose and Δ 18 O cellulose (mean sensitivity ≤ 0.32‰).
In general, predicted Δ 18 O cellulose decreased (increased) with increasing (decreasing) temperature in the sensitivity analysis (Fig. S4). If ϵ bio was allowed to vary with temperature (corresponding to the default parameterization), a 1°C increase (decrease) in temperature caused a 0.28‰ decrease (increase) in the predicted Δ 18 O cellulose . Conversely, if ϵ bio was set to be constant (at 27‰), the resulting temperature-sensitivity of predicted Δ 18 O cellulose was only −0.11‰°C -1 .

Estimates of isotopic signal integration time in cellulose agree with leaf lifespan and residence time of nonstructural carbon
This work presents a new allocation-and-growth model that generates realistic estimates for the turnover of the metabolic pool used to supply autotrophic respiration, the allocation of assimilates to shoot and root growth, and the integration (or residence time) of cellulose in growing leaf material in a drought-prone, temperate grassland ecosystem. When coupled to a recently published version of an 18 O-enabled soil-vegetation-atmosphere transfer model (MuSICA) parameterized for a grassland (Hirl et al., 2019), this allocation-and-growth model provided a faithful prediction of the fortnightly, seasonal and multi-annual variation of the 18 O composition of cellulose in mixed-species leaf samples from this pasture ecosystem. The coupled model also enabled an explicit evaluation of the current theory of 18 O-enrichment of cellulose (Δ 18 O cellulose ), namely the factors that determine the relationship between δ 18 O cellulose and δ 18 O source on the one hand and those of δ 18 O leaf on the other. We found compelling evidence that: ϵ bio , the average biochemical fractionation between the organic substrate used for cellulose synthesis and water in the field, is temperature-dependent, in close agreement with studies of submerged aquatic plants as shown by Relative humidity (RH), midday air temperature (Midday T air ), vapour pressure deficit (VPD), plant-available water (PAW) and wind speed (u) all represent midday (AE 3 h around noon) values averaged over the integration time Δt (400 GDD). Daily mean T air represents daily average air temperature over Δt, and precipitation amount (Precip) and short-wave radiation (R s ↓) represent the precipitation sum and the total sum of incoming short-wave radiation within Δt. ns, not significant (P > 0.05).
Ó 2020 The Authors New Phytologist Ó 2020 New Phytologist Trust New Phytologist (2020) www.newphytologist.com Sternberg & Ellsworth (2011); and p ex p x , the proportion of oxygen in cellulose derived from source water, responded to RH in accordance with the findings of Lehmann et al. (2017) and our own investigations (unpublished data) with L. perenne in controlled environments.
To predict δ 18 O cellulose accurately, the integration time used in the model calculations was particularly important. Model-data comparison of δ 18 O cellulose indicated a best match with 400 GDD. This was 14% shorter than the average observed leaf lifespan (LLS; 463 AE 56 GDD) of the main species (L. perenne, P. pratensis, T. officinale and T. repens) at the same site . A possible explanation is that Schleip et al. considered that LLS ended when a senescing leaf had 25% chlorotic tissue while leaf sampling for the present investigation comprised only fully green leaf blades and excluded any chlorotic leaf blades. Therefore, the collected leaves were probably not older than c. 320-450 GDD (compare with fig. 2 in Schleip et al., 2013), in close agreement with the observed best match for an integration time of c. 400 GDD (Fig. 4c). This value also agrees well with the mean residence time of structural carbon in aboveground biomass of the codominant species at the same site observed by Ostler et al. (2016).

Research
New Phytologist variations of RH and temperature were not significantly correlated at the integration time-scale (P = 0.50), meaning that RH and air temperature varied largely independently from each other. This factor aided/enabled the identification and distinction of the RH effect on p ex p x and of the temperature effect on ϵ bio .
Notably, daily mean temperatures over the integration time covered most of the temperature-sensitive range of ϵ bio (Sternberg & Ellsworth, 2011), explaining why the inclusion of a temperaturedependent ϵ bio was so critical for simulating Δ 18 O cellulose . Unfortunately, it is presently unknown how temperature influences the relevant biochemical mechanisms and their associated isotope effects that combine to yield the 'average biochemical fractionation between substrate for cellulose synthesis and water' ϵ bio . In that regard, however, it is interesting and encouraging that the temperature-dependence of ϵ bio observed here matched closely the patterns observed by Sternberg & Ellsworth (2011) in aquatic and heterotrophic systems. To the best of our knowledge, convincing general, field-based evidence supporting the requirement of a temperature-sensitive ϵ bio has not been published so far. Still, a possible temperature-sensitivity of ϵ bio has been discussed in modelling studies on tree-ring δ 18 O cellulose by Keel et al. (2016) and Lavergne et al. (2017) on tree species from different functional groups growing across a range of temperate to cold environments. Partial evidence in favour of a temperature-sensitive ϵ bio was presented by both groups, although overall evidence was inconclusive because of uncertainties in either the models or the climate and site data. Strong evidence for a temperature-sensitive ϵ bio is reinforced in our comprehensive model parameterization and validation study against an extensive field dataset. In particular, the leaf water modelling, partitioning of source (δ 18 O source )

(d)
Canopy conductance (mmol m −2 s −1 ) New Phytologist (2020) www.newphytologist.com and enrichment effects (Δ 18 O leaf ) on leaf water isotope composition (δ 18 O leaf ) (Hirl et al., 2019) permitted a very robust modelling of Δ 18 O cellulose that could be contrasted with observations of Δ 18 O cellulose , across a multiyear dataset using the same set of model parameters.
Sensitivity analyses of temperature effects on predictions of Δ 18 O cellulose demonstrated that c. 60% of the temperature effect in Δ 18 O cellulose was related to ϵ bio , with the remaining 40% connected to all other temperature-dependent processes influencing Δ 18 O cellulose . Interestingly, the temperature sensitivity of Δ 18 O cellulose (−0.28‰°C -1 ) was similar in magnitude but opposite in sign to the temperature sensitivity of δ 18 O rain (+0.36‰°C -1 ; Fig.  S2b). This caused a lack of correlation between δ 18 O cellulose and temperature, and highlights the need (and benefit) of studying both Δ 18 O cellulose and δ 18 O cellulose . Previous studies found a significant correlation between δ 18 O cellulose (mainly from tree rings) and temperature. These studies used different temporal (decadal or century-scale changes in contrast to the seasonal patterns studied here), spatial (latitudinal or altitudinal, as compared to the stationary observations in the present study) or aridity gradients, and may have also differed in the temperature sensitivity of Δ 18 O cellulose relative to that of δ 18 O rain (e.g. Anderson et al., 1998;Saurer et al., 2002;Kress et al., 2010;Song et al., 2011).
Evidence for an increase of p ex p x with air humidity, possibly linked to leaf water 18 O-enrichment along leaf blades Model-data comparison also strongly supported the inclusion of an RH-dependent p ex p x in the computation of Δ 18 O cellulose . This RH dependency of p ex p x was derived from the relationship between RH and the model-data residuals when Δ 18 O cellulose was predicted with a constant p ex p x of 0.556, combined with Eqns S7 and S9 (see Methods S2 and S3). A similar agreement between predicted and observed data was obtained when using the RHsensitivity of p ex p x observed by Lehmann et al. (2017) in a controlled environment study with L. perenne and D. glomerata, two of the codominant species present in our system. The RH dependency of p ex p x required to link Δ 18 O cellulose and Δ 18 O leaf has been discussed previously for a range of C 3 and C 4 grasses (Liu et al., 2016). Failure to predict Δ 18 O cellulose with a constant (i.e. RHinsensitive) p ex p x was not related to erroneous MuSICA predictions of Δ 18 O leaf , as these exhibited a virtually identical RH-response to observed Δ 18 O leaf (Hirl et al., 2019). Predictions of Δ 18 O assim by MuSICA were based on bulk leaf water Δ 18 O leaf and the common assumption that new assimilates are in isotopic equilibrium with Δ 18 O leaf (Δ 18 O assim = Δ 18 O leaf + ϵ bio ; Barbour, 2007). Lehmann et al. (2017) showed that this assumption may not always hold in grasses where RH was affecting the δ 18 O of sucrose (δ 18 O sucrose ) more strongly than δ 18 O leaf . These results support the view that Δ 18 O sucrose is not always equal to Δ 18 O leaf + ϵ bio in grasses, and that RH affects the divergence between Δ 18 O leaf and the Δ 18 O of the water in which assimilation and sucrose synthesis takes place (Δ 18 O suc-water ). This divergence between Δ 18 O leaf and Δ 18 O suc-water may be related to the gradient of 18 O enrichment along leaf blades observed by Helliker & Ehleringer (2000, 2002a and Gan et al. (2003), which is particularly steep at low RH. Uncertainties on Δ 18 O suc-water affect the calculation of p ex p x . p x is commonly calculated from a twoend-member mixing model, the end-members being the Δ 18 O of water at the site of cellulose synthesis (Δ 18 O cel-water ) and Δ 18 O suc-water , with the latter approximated as Δ 18 O leaf : p x = 1 − Δ 18 O cel-water /Δ 18 O suc-water . A deviation between Δ 18 O suc-water and Δ 18 O leaf affects the calculation of p x . However, for grasses, this error must be small, as the water in the leaf growth and differentiation zone of grasses is close to source water (Δ 18 O cel-water of c. 0) so that p x is c. 1 (Liu et al., 2017a). Accordingly, the divergence between Δ 18 O suc-water and Δ 18 O leaf noted by Lehmann et al. (2017) mostly affects the estimate of p ex .
Other factors that have been related to variations of p ex have included 18 O-exchange during phloem transport from the leaves to the growing tissue (Gessler et al., 2007; but see Gessler et al., 2013) and futile cycling of the growth substrate within the growing tissue (Barbour, 2007), potentially linked to cellular growth rate or other cellular anatomical features, such as the lumen area of cambial cells (Szejner et al., 2020). The rate of futile cycling has been related to the turnover of the water-soluble carbohydrate pool in Ricinus communis (Song et al., 2014). However, we could not detect a significant relationship between the modelled mean residence time of substrate in the metabolic pool (τ pool ) and p ex , estimated from p ex p x with p x = 1 (Fig. S5), although statistical noise may have been a factor. Also, modelled integration time (in days)which is a function of leaf appearance, growth and senescence rates that are all closely coordinated in grasses (Lemaire et al., 2000;Schleip et al., 2013) did not correlate with estimates of p ex . Lastly, in a recent study, we did not observe an effect of RH on the hydraulic conductance of L. perenne vegetative plants in controlled environments (Baca Cabrera et al., 2020), providing no (physiological) indication of an effect of RH on xylem lumen area. Clearly, there is a great need for detailed biochemical investigations of the factors underlying the variation of p ex .
Predictions of canopy conductance from 18 O-enrichment of cellulose required a temperature-dependent ϵ bio and humidity-sensitive p ex p x Finally, this work demonstrates a significant negative relationship between midday canopy conductance (g canopy ) and both Δ 18 O cellulose and δ 18 O cellulose . Although used in multiple studies to interpret variation in Δ 13 C and intrinsic water-use efficiency, experimental or empirical evidence for the relationship between Δ 18 O cellulose (or δ 18 O cellulose ) and stomatal or canopy conductance is relatively scarce (Barbour et al., 2000;Grams et al., 2007;Sullivan & Welker, 2007;Moreno-Gutiérrez et al., 2012). In principle, the correlation between g canopy and Δ 18 O cellulose could be caused by any environmental factor that affects stomatal opening and consequently alters leaf temperature, leaf-to-air vapour pressure difference, and thus Δ 18 O leaf , Δ 18 O suc-water , ϵ bio and ultimately Δ 18 O cellulose . Correlations between Δ 18 O cellulose and stomatal conductance observed in previous studies were evoked by treatment differences in air CO 2 and O 3 concentration and by interplant competition (Grams et al., 2007), or by variation in soil temperature and soil water status (Sullivan & Welker, 2007).

New Phytologist
In a previous study, we found a significant effect of both soil water content and RH on Δ 18 O leaf , and a significant relationship between soil moisture and canopy conductance (Hirl et al., 2019). Remarkably, however, we could not identify a drought-related increase of Δ 18 O cellulose , which would indicate a drought-related decrease of stomatal conductance (Fig. S6). We believe that absence of a drought signal in Δ 18 O cellulose is caused by the low assimilation rate and deceleration (or complete cessation) of leaf growth under drought. By contrast, both Δ 18 O cellulose (Table 2; Fig. 6; Fig. S3a,c) and g canopy (R 2 = 0.13; P < 0.01) were significantly related to RH. The RH-sensitivity of Δ 18 O cellulose was even higher than that of Δ 18 O leaf (−0.15‰ % −1 ), due to the relationship between RH and p ex p x (as already discussed). If Δ 18 O suc-water is more strongly related to evaporative conditions than Δ 18 O leaf , the stomatal conductance signal in cellulose may in fact be more pronounced than previously expected solely from the climate-sensitivity of bulk leaf water. This lends further support to studies that aim to reconstruct stomatal conductance from δ 18 O cellulose .

Conclusions and perspective
This study provides strong evidence for important climate-sensitive variation of the parameters p ex p x and ϵ bio of the Barbour-Farquhar model in a terrestrial ecosystem. Elucidation of those sensitivities demanded and relied on the use of a carefully parameterized, 18 O-enabled process-based soil-plant-atmosphere transfer model and the existence of a detailed multiseasonal observational data set of 18 O signals in all relevant ecosystem water pools and leaf cellulose. As a result, our work demonstrates at the ecosystem scale the validity of earlier conclusions that were drawn from model-type, experimental studies conducted in steady-state, controlled environments, but have not been generally considered in the relevant fields of terrestrial ecosystem sciences. Our results have important implications for the interpretation of seasonal, interannual and geographical variations of Δ 18 O cellulose and δ 18 O cellulose in terms of plant physiology and climate. Existing datasets, for example on tree-ring Δ 18 O cellulose or δ 18 O cellulose from boreal to subtropical ecosystems, should benefit from a re-evaluation that considers the varying nature of biochemical fractionation and isotopic exchange with source and leaf water. In addition, there is a great need for more targeted studies investigating the mechanisms underlying the large observed effect of temperature on ϵ bio and the RH-dependent variation of p ex p x . Such knowledge will help us to better understand the physiological information that is imprinted in cellulose-18 O signals from past and present environments.

Supporting Information
Additional Supporting Information may be found online in the Supporting Information section at the end of the article.  Ball et al., 1987), to maximum rate of carboxylation and potential rate of electron transport (V cmax and J max ), and to incoming short-wave radiation and wind speed.     Methods S1 The MuSICA model. Ó 2020 The Authors New Phytologist Ó 2020 New Phytologist Trust New Phytologist (2020) www.newphytologist.com

Methods S2
The allocation-and-growth model.
Methods S3 Derivation of the p ex p x -RH function for the standard simulation.
Table S1 Correlation between meteorological variables.
Please note: Wiley Blackwell are not responsible for the content or functionality of any Supporting Information supplied by the authors. Any queries (other than missing material) should be directed to the New Phytologist Central Office.
New Phytologist is an electronic (online-only) journal owned by the New Phytologist Foundation, a not-for-profit organization dedicated to the promotion of plant science, facilitating projects from symposia to free access for our Tansley reviews and Tansley insights.
Regular papers, Letters, Research reviews, Rapid reports and both Modelling/Theory and Methods papers are encouraged. We are committed to rapid processing, from online submission through to publication 'as ready' via Early View -our average time to decision is <26 days. There are no page or colour charges and a PDF version will be provided for each article.
The journal is available online at Wiley Online Library. Visit www.newphytologist.com to search the articles and register for table of contents email alerts.
If you have any questions, do get in touch with Central Office (np-centraloffice@lancaster.ac.uk) or, if it is more convenient, our USA Office (np-usaoffice@lancaster.ac.uk) For submission instructions, subscription and all the latest information visit www.newphytologist.com New Phytologist (2020) Ó 2020 The Authors New Phytologist Ó 2020 New Phytologist Trust www.newphytologist.com