Disentangling temporal and population variability in plant root water uptake from stable isotopic analysis: a labeling study

. Isotopic labeling techniques have the potential to minimize the uncertainty of plant root water uptake (RWU) profiles estimated through multi-source (statistical) modeling, by artificially enhancing soil water isotopic gradient. 15 Furthermore, physical models can account for hydrodynamic constraints to RWU if simultaneous soil and plant water status data is available. In this study, a population of tall fescue ( Festuca arundinacae cv Soni) was grown in a macro-rhizotron setup under semi-controlled conditions to monitor such variables for a 34-hours long period following the oxygen stable isotopic ( 18 O) labeling of deep soil water. Aboveground variables included tiller and leaf water oxygen isotopic compositions as well as leaf water 20 potential ( ψ leaf ), relative humidity, and transpiration rate. Belowground profiles of root length g., temporal dynamics of leaf water potential and transpiration rate).

signature, and supported that the local increase of soil water content and formation of a peak of labelled water observed overnight were due to hydraulic lift.

Introduction
Since the seminal work of Washburn and Smith (1934) where it was first reported that willow trees did not fractionate hydrogen stable isotopes in a hydroponic water solution during root water uptake (RWU), water stable isotopologues 65 ( 1 H 2 H 16 O and 1 H 2 18 O) have been used as indicators for plant water sources in soils. In their review, Rothfuss and Javaux (2017) reported in the period 2015-2016 about no less than 40 publications in which RWU was retrieved from stable isotopic measurements. Novel measuring techniques (e.g., cavity ring-down spectroscopy -CRDS and off-axis integrated cavity output spectroscopy -ICOS) providing ways for fast and cost-effective water stable isotopic analyses certainly enable and emulate current research in that field. Water stable isotopologues are no longer powerful tracers waiting for technological 70 developments (Yakir and Sternberg, 2000) but are on the verge to be used to their full potential for addressing ecohydrological research questions and identify processes in the soil-plant-atmosphere continuum (Werner et al., 2012;Dubbert and Werner, 2019;Sprenger et al., 2016).
The isotopic determination of RWU profiles is based on the principle that the isotopic composition of xylem water at the outlet of the root system (i.e., in the first aerial and non-transpiring node of the plant) equals the mean value of the soil water 75 isotopic composition across contributing sources, weighted by their relative contribution to RWU. Results come only with reasonable precision when (i) the soil water isotopic composition depth gradient is strong and (ii) the temporal dynamics of RWU and soil water isotopic composition is relatively low. Condition (i) is fulfilled mostly at the surface of the soil, while soil water isotopic composition gradients become usually lower or null with increasing depth (due to the isotopic influence of the groundwater table). Condition (ii) is often neglected but is required due to the instantaneous nature of the sap flow 80 samples. As illustrated by Oerter and Bowen (2019), the lateral variability of the soil water isotopic composition profiles can become significant in the field and could have great implications on the representability and meaningfulness of isotopicderived estimate of RWU profiles.
To overcome these limitations, labeling pulses have been increasingly used in recent works to artificially alter the natural isotopic gradients (e.g., Beyer et al., 2016;Beyer et al., 2018;Grossiord et al., 2014;Jesch et al., 2018;Volkmann et al., 85 architectures and root segment permeabilities (Meunier et al., 2017c). Only a handful of studies coupled isotopic measurements in plant tissues and soil material with models describing RWU in a mechanistic manner. For instance, Meunier et al. (2017a) could both locate and quantify the volume of redistributed water by Lolium multiflorum by labeling of the soil with 18 O enriched water under controlled conditions.
Building on the work of Meunier et al. (2017a), the objective of the present study is to investigate on the plant and 100 environmental factors that affect the tiller xylem water isotopic composition in a population of Festuca arundinacae cv Soni.
(tall fescue), in order to unveil pitfalls and opportunities in methods used to predict RWU profiles from stable isotopic analyses following a labeling pulse.

Material and methods
Our experiment consisted in supplying labeled water from the bottom to a macro-rhizotron in which tall fescue was grown. 105 Data on soil and plant isotopic signature and hydraulic status were monitored for 34 hours. In the following, the oxygen isotopic composition of water will be expressed in per mil (‰) on the "delta" (δ 18 O) scale with respect to the international water standard V-SMOW (Gonfiantini, 1978).

Rhizotron experimental setup
The macro-rhizotron (dimensions: 1.6 m x 1.0 m x 0.2 m, see picture in Appendix A) was placed inside a glasshouse (INRA 110 Lusignan, France), where it was continuously weighed (KE1500, Mettler-Toledo, resolution: 20 g) to monitor water effluxes (i.e., bare soil evaporation or evapotranspiration). Underneath the soil compartment and in contact with it, a water reservoir (height: 0.1 m) filled with gravel acted as water table and allowed the supply of water to the rhizotron. The rhizotron was equipped with two sets of CS616 water content reflectrometer profiles (Campbell Scientific, USA) with 30 cm long probe rods positioned at six depths (-0.05, -0.10, -0.30, -0.60, -1.05 and -1.30 m) and one profile of tensiometers (SMS 2000, 115 SDEC-France) located at four depths (-0.05, -0.10, -0.30, and -0.60 m) in order to monitor the evolutions of soil water volumetric content (θ, in m 3 m -3 ) and matric potential. Finally, relative humidity (RH, %) was recorded above the vegetation with one humidity and temperature probe (HMP45D, Vaisala, Finland). The transparent polycarbonate sides (front and back) allowed the daily observations of root maximal depth. The experimental setup allowed precisely controlling the amount and δ 18 O of soil input water. Another important feature was the soil depth (i.e., 1.60 m) which minimized the influence of the 120 water table on superficial layers water content and δ 18 O.

Soil properties and installation
The soil is classified as District Cambisol (particle size distribution: sand 15%, silt 65%, clay 20%) with a dry bulk density of ρ b = 1420 kg m -3 . Prior installation in the rhizotron, the substrate was sieved at 2 mm and dried out in an air oven at 110 compacted in order to reach the above-mentioned value of dry bulk density. Soil water retention data was derived from synchronous measurements of soil water content (CS616 soil water reflectometers, Campbell Scientific, USA) and potential (Sols -PST55 micro-psychrometers, WESCOR, and SMS 2000 tensiometers, SDEC) from saturated to residual water content. The closed-form water retention curve (van Genuchten, 1980) was then parametrized by minimizing the root mean square difference between its position and measured soil water content-potential couples (Appendix B). 130

Experimental protocol
After installation, the soil was gradually flooded with local water (δ 18 O = -6.8 ‰) from the bottom reservoir up to the top of the profile for a period of three days in order to reduce as much as possible the initial lateral and vertical heterogeneities in water content and δ 18 O. The tall fescue (Festuca arundinaceae cv Soni) was sown at a seeding density of 3.6 g m -2 (which corresponds for the rhizotron surface area of 0.2 m 2 to roughly 300 plants) when soil water content reached 0.25 m 3 m -3 135 (corresponding to pF 2.3) at -0.05 m, as measured by the soil water sensors, and emerged 12 days later.
166 days after seeding (DaS 166) the following conditions were fulfilled: (i) there was a strong soil water content gradient between the soil deep [-1.5 m, -1.0 m] and superficial [-0.3 m, 0 m] layers, (ii) the tall fescue roots had reached a depth of -1.5 m (observed through polycarbonate transparent slides). That same day at 17:00, the reservoir's water was labelled and its isotopic composition measured at +470 ‰. Soil was sampled on four occasions from the surface down to -1.3 m: before 140 labeling (DaS 166 -15:45) and after labeling on DaS 167 -07:00, DaS 167 -17:00 and DaS 168 -05:00 for the determination of soil gravimetric water content (θ grav , in kg kg -1 ) and isotopic composition (δ soil , in ‰). Gravimetric water content was then converted to volumetric water content (θ = θ grav *ρ b /ρ w , in m 3 m -3 , where ρ b is the bulk soil density and ρ w is the water density). On 40 occasions during a 34-hour long period three plants were sampled from the vegetation (i.e., 120 plants were sampled in total from the cover). Each plant's tillers and leaves were pooled into two separate vials. Dead 145 material as well as the oldest living leaf around each tiller were removed in order not to contaminate tiller samples with transpiring material (Durand et al., 2007). In addition, air water vapor was collected from the atmosphere of the surrounding of the rhizotron. The air was run at a flow rate of 1.5 l min -1 through two glass cold traps in series immersed in a mixture of dry ice and pure ethanol at -80°C. Water from tillers, leaves, and soil samples were extracted by vacuum distillation and their isotopic compositions (i.e., δ tiller , δ leaf , and δ soil ) together with that of atmospheric water vapor (δ atm ) were measured with 150 an IRMS (Isoprep 18 -Optima, Fison, Great-Britain, precision accuracy of 0.15 ‰). Finally, leaf water potential (ψ leaf , in MPa) was monitored with a pressure chamber on two leaves per sampled plant, and transpiration (T, m d -1 ) was derived from the changes in mass of the rhizotron at the same temporal scale as plant sampling.
Root biomass was determined from the horizontal sampling of soil between the polycarbonate sides using a 2 cm diameter auger at -0.02, -0.08, -0.10, -0.40, -0.55, -0.70, -0.90, -1.10, and -1.30 m soil depth. Each depth was sampled once to 155 thrice. Each soil core was washed of soil particles and roots were collected over a 0.2 mm mesh filter, and dried at 60°C for 48 hours. Finally, Root Length Density (RLD, in m root (m soil) -3 ) distribution was determined from the root dry mass using the specific root length determined by Gonzalez-Dugo et al. (2005) specifically for tall fescue (95 m root (g root) -1 ). The reader is referred to Appendix C for an overview of the type and timing of the different destructive measurements during the intensive sampling period. 160

Modeling of RWU and δ tiller
60 tall fescue root systems with rooting depths ranging from -1.30 to -1.60 m depth (based on our own observations and those of the litterature, e.g., Schulze et al., 1996;Fan et al., 2016) were modeled to represent the community in the rhizotron.
The root systems were generated with the root architecture simulator CRootBox (Schnepf et al., 2018) so that the simulated RLD matched observations (Fig. 1a). Each simulated root system was considered as representative of a "class" among the 165 plant population. Each root system class was assumed to occupy one sixtieth of the total horizontal area. To simulate RWU by the 60 classes while limiting the computing time in the inverse modeling scheme, we considered each of them as a "big root" hydraulic network with equivalent radial and axial hydraulic conductances (thus neglecting architectural aspects but accounting for their respective root length density profiles).
The radial soil-root conductance between the bulk soil and each class's (i) root surface in soil layer j (K radial,j , m 3 MPa -1 d -1 ) 170 was assumed as variable in time (t): with B (dimensionless) a geometrical factor approximating the dimension of the domain between the bulk soil and root surface as radial, as given by Schroeder et al. (2009) (2) 175 where ρ (dimensionless) represents the ratio of the distance between roots and the root averaged diameter. The averaged distance between roots can be deduced from the observed root length density (RLD j , m m −3 ) and root radius (r root , m): In Eq. (1), the soil hydraulic conductivity function of Mualem (1976) and van Genuchten (1980) was used: where k sat (m 2 MPa −1 d −1 ), m (dimensionless) and λ (dimensionless) are soil hydraulic parameters (with m = 1 -2/n) and Se j , the relative water content (dimensionless), is computed from the saturated (θ sat , m 3 m -3 ) and residual (θ res , m 3 m -3 ) water contents as: The last term to define in Eq. (1) is the root length in each soil layer (l root , m). Unlike the geometrical parameter , which 185 defines a domain geometry between the bulk soil and roots of the overall population, the l root term is class specific (i) and uses the simulated root length density profiles over an area corresponding to one sixtieth of the total setup horizontal area: with ∆Z (m) and A soil (m 2 ) the soil layer thickness and horizontal surface area, respectively.
To finalize the connection between root xylem and shoot, axial conductances per root system class (K axial , m 3 MPa −1 d −1 ) 190 were derived from an equivalent "big root" specific axial conductance per root system class (k axial , m 4 MPa −1 d −1 ) as: At each time step, both the total soil-root system conductance (K soil-root , m 3 MPa −1 d −1 ) and the standard sink distribution (SSF, dimensionless, summing up to 1) were calculated from K radial and K axial , using the algorithm of Meunier et al. (2017b).
The variable K soil-root represents the water flow per unit water potential difference between the bulk soil and the leaf 195 (assuming a negligible stem hydraulic resistance), and SSF the relative distribution of water uptake in each soil layer under vertically homogeneous soil water potential conditions (Couvreur et al., 2012).
Adding soil hydraulic conductance to the one-dimensional hydraulic model of Couvreur et al. (2014) yields the following solutions of leaf water potential (ψ leaf , MPa) and water sink terms (S, d -1 ) whose formulation approaches that of Nimah and Hanks (1973): 200 Where one sixtieth of the overall transpiration rate ( ) is allocated to each class, and ψ soil,j (Mpa) is the soil water potential in soil layer j.
where axial conductances were assumed to be large enough for soil−root to control the compensatory RWU which arise 205 from a heterogeneously distributed soil water potential (Couvreur et al., 2012).
Finally, the tiller water isotopic composition (δ tiller ) was calculated as the average of local soil water isotopic compositions (δ soil ) weighted by the relative distribution of positive water uptakes (i.e., not accounting for δ soil at locations where water is exuded by the root), assuming a perfect mixture of water inside the root system (Meunier et al., 2017a): Like in the experiment, δ tiller from three plants were randomly pooled at each observation time. As the experiment included about 300 plants (value based on seeding density, see section 2.3), and we explicitly simulated 60 individual RLD profiles, we assumed that each of them was representative of a class of 5 identical plants, totaling 300 plants. A hundred pools of 3 plants (possibly including several plants of the same class) were randomly selected in order to obtain the pooled simulated δ tiller by arithmetic averaging. 215 The unknown parameters of the soil-root hydraulic model, i.e., the root radial conductivity (L pr ), the root axial conductance (k axial ), the soil saturated hydraulic conductivity (k sat ), and the soil tortuosity factor (λ) were finally determined by inverse modeling. For details on the procedure, the reader is referred to Appendix D. In order to evaluate the robustness of the hydraulic model predictions (parametrized solely based on the reproduction of shoot observations in the inverse modeling scheme) from independent perspectives, we also compared predictions and 220 measurements over 4 quantitative "soil-root domain" criteria: (i) the depth at which the transition between nighttime water uptake and exudation takes place, (ii) quantities of exuded water and overnight increase of soil water content, (iii) the enrichment of labelled water at the depth where water content increase is observed overnight, and (iv) the order of magnitude of the optimal root radial conductivity value as compared to literature data in tall fescue.
Finally, and as a comparison point, the Bayesian inference statistical model SIAR (Parnell et al., 2013) was used to 225 determine the profiles of water sink terms of ten identified potential water sources. These water sources were defined to originate from 10 distinct soil layers for which corresponding δ soil values were computed (Rothfuss and Javaux, 2017). SIAR solely bases its estimates from the comparison of δ tiller observations to the isotopic compositions of the soil water sources (δ soil ). For this, δ tiller measurements were pooled in twelve groups corresponding to different time periods, selected to best reflect the observed temporal dynamics of δ tiller . The reader is here referred to Appendix E for details on the model 230 parametrization and running procedure. 166-168 (LAI = 5.6) and the observed surface θ values lead to assume that the rhizotron water losses were due to transpiration flux solely (i.e., evapotranspiration = transpiration). The soil water isotopic exponential-shaped profiles were 245 the product of fractionating evaporation flux, and to a great extent when the soil was bare or when the tall fescue cover was not fully developed. The differences in soil water isotopic profile observed at the four different sampling dates were therefore either due to lateral heterogeneity (e.g., upper soil layers), to the soil capillary rise of labelled water from the reservoir (deep soil layers), or to the hydraulic redistribution of water through roots. The observed RLD profile (Fig. 1a) https://doi.org/10.5194/hess-2019-543 Preprint. Discussion started: 21 October 2019 c Author(s) 2019. CC BY 4.0 License.
showed a typical exponential shape, i.e., maximal at the surface (5.42 ± 0.34 cm cm -3 ) down to a minimal at -1.10 m (0.540 250 ± 0.35 cm cm -3 ), while it increased again from the latter depth up to a value of 1.660 cm cm -3 at -1.30 m. This significant trend was most probably a direct consequence of the high soil water content level in this deeper layer.

Plant water and isotopic temporal dynamics
The temporal variation of δ tiller (Fig. 3a) was found to be either (i) moderate during day and night, i.e., from DaS 167 -06:00 to 11:00 (δ tiller = -2.6 ± 1.4 ‰) and from DaS 167 -21:30 to DaS 168 -00:00 (δ tiller = -2.7 ± 0.4 ‰), or (ii) strong during the 255 day, i.e., from DaS 167 -11:00 to 18:00 (maximum value of 20.9 ‰ at DaS 167 -12:40), or else (iii) strong during the night, i.e., from DaS 167 -04:00 to 06:00 (max = 36.4 ‰ at DaS 167 -05:15) and from DaS 168 -00:00 to 06:00 (max = 14.6 ‰ at 28:00, DaS 168). Note that transpiration (Fig. 3b) occurred also at night during the sampling period, due to relatively high temperature in the glasshouse leading to a value of atmospheric relative humidity smaller than 85%, Fig. 3b). From 12:00 to 14:00 and from 16:00 to 17:00 on DaS 167 (case (ii)) high values of leaf transpiration corresponded to high values of  tiller . 260 Figure 4 shows that variables describing plant water status, i.e., T and RH (Fig. 4a) and T and ψ leaf (Fig. 4b) were well correlated: coefficient of determination R 2 was equal to 0.78 and 0.70 for the entire experimental duration, respectively. However, linear relationships between water status and isotopic variables were either inexistent, e.g., between T and δ tiller (R 2 =0.01, Fig. 4c) and between ψ leaf and δ tiller (R 2 =0.00, Fig. 4h) or characterized by a low R 2 and high p-value (e.g., between 265 T and δ leaf , R 2 =0.43, p>0.05, Fig. 4d). The partial temporal disconnection between δ leaf and T could not be attributed to problems of the isotopic methodology, during e.g., the vacuum distillation of the water from the plant tillers and leaves: water recovery rate was always greater than 99 % and Rayleigh distillation corrections were applied to standardize the observed isotopic composition values to a 100 % water recovery (based on the comparison of sample weight loss during distillation and mass of collected distillated water). The evolution of δ leaf was strongly correlated with that of δ tiller during the 270 day (R 2 = 0.90) whereas non-correlated during the night (R 2 = 0.00, Fig. 4j). These observed correlations are in agreement with the Craig and Gordon (1965) model revisited by Dongmann (1974) and extensively used in the current literature (e.g., Dubbert et al., 2017): at isotopic steady-state, δ leaf is a function of the input water isotopic composition (δ tiller ) among other variables, i.e., leaf temperature (not measured during the experiment), stomatal and boundary layer conductances, isotopic composition of atmospheric water vapor, and relative humidity. 275

Partial decorrelation between water and isotopic state variables
It is generally difficult to observe a statistically significant δ leaf -δ tiller (Fig. 4j) relationship at this temporal scale under natural abundance conditions in the field since the soil water isotopic weak gradient translates into weaker δ tiller temporal dynamics.
The quality of linear fit between δ leaf and δ tiller data collected during the day (R 2 =0.90) was made possible in this specific experiment by the artificial isotopic labeling pulse that enhanced the soil water isotopic gradient, which in turn increased the range of variation of δ tiller , ultimately highlighting the δ leaf -δ tiller temporal correlation. Air relative humidity is a driving 280 variable of δ leaf in the model of Dongmann (1974)  atmospheric water vapor isotopic composition inside the glasshouse. An overall significant linear correlation was observed between RH and δ leaf during the experiment (R 2 =0.57, Fig. 4g). During the two night periods (i.e., from 04:00-06:00 and from 20:30-07:00), as relative humidity increased in the glasshouse (51 % < RH < 85 %, Fig. 3b), the influence of the isotopic labeling of the tiller water (due to the labeling of deep soil water) through term (1-RH)•δ tiller decreased to the benefit 285 of term RH•δ atm (with δ atm values ranging from -15.9 to -10.7 ‰, mean = -13.1±1.6‰, data not shown). This was especially visible between 04:50 and 06:00 on DaS 167 and between 01:00 to 06:00 on DaS 168, when δ tiller reached greater values than From a different perspective, as three plant water samples were pooled to reach a workable volume for the isotopic analysis at each observation time without replicates, the isotopic signal fluctuations may reflect both its temporal dynamics and its 290 variability within the plant population.

Rooting depth and transpiration rate control δ tiller and ψ leaf fluctuations, respectively
Despite the use of a global optimizer and 4 degrees of freedom (L pr , k axial , k sat , λ, see optimal values in Table 1) specifically aiming at matching the simulated and observed temporal dynamics of δ tiller , none of the 60 root system classes or average 295 population could reproduce the measured fluctuations in time (R 2 =0.00, Fig. 5a), regardless of the weight attributed to this criterion in the objective function. However, the predicted versus observed average δ tiller and its standard deviation for the overall dataset were not significantly different (6.6 ± 8.4 ‰ and 3.7 ± 8.4 ‰, respectively). Besides, the simulated ψ leaf fitted well the observations (R 2 =0.67, overall distributions: -0.175 ± 0.053 MPa and -0.177 ± 0.053 MPa, respectively, Fig. 5c).
When analyzing the distributions of ψ leaf and δ tiller per maximum root system depth ( Fig. 5b and d), it appears that the ψ leaf 300 signal is not sensitive to the rooting depth (Fig. 5d), while δ tiller is more sensitive to rooting depth than to the temporal evolution of the plant environment (Fig. 5b). This leaves us with two hypotheses. The "rollercoaster hypothesis": δ tiller rapidly goes up and down with all the population on board of the same car (i.e. little variability within the population, unlike predictions in Fig. 5a, but like the simulated ψ leaf in Fig. 5c). If that is correct, the physical model lacks a process that would capture the observed temporal fluctuations of 305 δ tiller . The "swarm pattern hypothesis": δ tiller is rather stable in time but its values within the plant population are dispersed like in a flying swarm, so that δ tiller values sampled at different times fluctuate, not due to temporal dynamics but to the fact that different individuals are sampled (Fig. 5a).
The model suggests that the tall fescue population ψ leaf follows a "rollercoaster" dynamics driven by transpiration rate, while the population δ tiller follows a "swarm" pattern driven by the maximum rooting depth of the sampled plants. As no correlation 310 could be expected between the drivers, transpiration rate and the maximum rooting depth of the sample plants, discrepancies between fluctuations of isotopic and hydraulic variables are not surprising any more. In future experiments and in the specific context of labeling pulses, sampling more plants at each observation time would help disentangle the spatial from temporal sources of variability of ψ leaf and δ tiller . It would however be at the cost of the temporal resolution of observations, or would necessitate a larger setup with more plants in the case of controlled conditions 315 experiments.

Independent observations support the validity of the hydraulic model predictions
In the last 12 hours of the experiment (DaS 167 -17:00 to DaS 168 -05:00), the measured soil water content increased by 2.9% at -0.9 m depth, which could be a sign of nighttime hydraulic redistribution. During the same period, the physical model predicted a cumulative water exudation sufficient to increase soil water content by 0.38 %, as soil water potential was 320 sufficiently low to generate reverse flow, but high enough not to disrupt the hydraulic continuity between soil and roots (Carminati and Vetterlein, 2013;Meunier et al., 2017a). While this increase is smaller than the observed water content change, it is only a component in the soil water mass balance. Given the soil water potential vertical gradient, upward soil capillary water flow may have accounted for another part of the observed moisture change. Experimental observations also show that δ soil increased by 1.0 ‰ at 0.9 m depth during that time (-6.2 ‰, a value significantly higher than -7.1 ‰ ± 0.1 ‰ 325 at earlier times), while our simulations of hydraulic redistribution generated an increase of δ soil by 0.41 ‰. As soil capillary flow may not generate local maxima of δ soil (no enrichment observed at surrounding heights, see Fig. 2b), and soil evaporation is assumed negligible at that depth, it is likely that the observed local enrichment was entirely due to hydraulic redistribution, which would then be underestimated by a factor of about 2.5 in our simulations. Increasing water exudation by a factor 2.5 would imply a simulated water content change due to exudation of 0.95% absolute water content, which 330 remains compatible with the experimental observation.
Between -1.1 m and -0.9 depth, the nighttime water flow pattern transitioned from exudation to uptake in both measurements and predictions. At -1.1 m, the model predicted a cumulative water uptake sufficient to decrease soil water content by 0.89 %, as compared to the observed 1.41 % total soil water content decrease. The remaining 0.5% water content decrease may have contributed to the recharge to the soil layers above through capillary flow, which was not simulated. 335 The fitted L pr value (2.3 10 -7 m MPa -1 s -1 ) was in the range found by Martre et al. (2001) in tall fescue (2.2 10 -7 ± 0.1 m MPa -1 s -1 ) and falls in the range obtained by Meunier et al. (2017a) for another grass (Lolium multiflorum Lam., 6.8 10 -8 to 6.8 10 -7 m MPa -1 s -1 ). Our k axial value cannot be compared to values of axial root conductance from the literature as it transfers the water absorbed by roots from 5 plants in a single "big root" per class of root system. The optimal value of k sat was quite high (Table 1) but reportedly very correlated to λ, so that the low value of the latter compensated the high value of 340 the former, thus they should be considered as effective rather than physical parameters.

Do root water uptake profiles predicted by hydraulic and Bayesian models differ?
The root water uptake dynamics predicted by the mechanistic model are shown in Fig. 6a. Appendix E). The main differences were the following: (i) in the upper soil layers where the soil water potential was lower -1.5 MPa, the statistical model predicted water uptake, which is theoretically impossible given the leaf water potential above -0.4 MPa (van Den Honert, 1948); (ii) In the top half of the profile, the physical model predicted exudation at a rate limited by the low hydraulic conductivity between root surface and bulk soil, with a peak at night, at -0.9 m depth (quantitative analysis in previous section); (iii) Below -1.0 m depth, the water uptake rate predicted by the statistical model steadily 350 increased with depth while that of the physical model was more uniform, likely due to axial hydraulic limitation (e.g., Bouda et al., 2018) counteracting the increasing soil water potential with depth.

Conclusion
In the present study, light could be shed on RWU of Festuca arundinacae by specifically manipulating the lower boundary conditions for water content and oxygen isotopic composition. The new version of the one-dimensional model of Couvreur 355 et al. (2014) implemented here accounted for both root and soil hydraulics in a population of "big" root systems of known root length density profile. This approach underlined the high sensitivity of δ tiller to rooting depth and suggested that if δ tiller is measured on a limited number of individuals, its variations in time may reflect the heterogeneity of rooting depth within the population, rather than temporal dynamics which was minor in our simulations. The model avoided the prediction of water uptake at locations where it was physically unavailable (e.g., in the first half of the soil profile), by accounting for water 360 potential differences observed between the leaves and the soil, and explained quantitatively the local isotopic enrichment of soil water as the occurrence of nighttime Hydraulic Lift at -0.9 m depth. On the other hand, the Bayesian statistical approach tested for comparison, which was driven by isotopic information solely, naturally translated the observed changes of δ tiller into profound temporal dynamics of RWU, at the expense of eco-physiological consideration (e. g., temporal dynamics of leaf water potential and transpiration rate). 365 This case study highlights (i) the potential limitations water isotopic labeling techniques for studying RWU, especially when water addition is localized and not broadcasted in the soil therefore creates strong isotopic depth-gradients. As already pointed out in the review of Rothfuss and Javaux (2017), the study also (ii) underlines the interest of complementing in-situ isotopic observations in soil and plant water with information on soil water status and plant eco-physiology; it finally (iii) calls for the use of simple soil-root models for inversing isotopic data and gain insights into the RWU process.
The experimental part of this study was financed by the French national scientific program EC2CO/CITRIX.

Data sets
Upon acceptance, all research data needed for creating plots will be available in reliable FAIR-aligned data repositories with 375 assigned DOIs.

Author contribution
TB, JLD, and PB designed the experiments and TB, JLD, PB, and YR carried them out. VC and FM developed the physically-based root water uptake model code and performed the simulations. YR performed the statistical simulations. VC and YR prepared the manuscript with contributions from all co-authors. 380
For each of the twelve time periods: (i) the function siarmcmcdirichletv4 of the SIAR R package (https://cran.r-30 project.org/web/packages/siar/index.html) was run 500,000 times with prescribed burnin and thinby equal to 50000 and 15, respectively. The output of the model (i.e., the a posteriori rRWU distribution across the ten soil water sources J) was obtained from a flat Dirichlet a priori rRWU distribution (i.e., rRWU J =1/10); (ii) the 'best run' (br, dimensionless) was selected from SIAR's output. It was defined as the closest solution of relative contributions across sources from the set of most frequent values (mfv, dimensionless), i.e., the relative 35 contribution with the greatest probability of occurrence. The best run was identified as minimizing the objective function below, i.e., the RMSE (root mean square error) with respect to the set of mfv J :