Comparing root water uptake profile estimations from an isotope-calibrated mechanistic model and a mixing model

The root water uptake profile (RWUP) reflects a plant’s survival strategy and controls evapotranspiration and carbon fluxes. Despite its importance, there is still no reliable method for reconstructing this profile. In this study, we applied and compared two possible approaches to a case study in a conifer plantation: an isotope-calibrated mechanistic model and a mixing model with a bell-shaped approximation. Our results show that, after calibrating the hydrologically-active root density profile, the mechanistic model gave a good estimation of the xylem water isotope delta (δx); even though the measured root density was greater in shallower soils, water uptake occurred throughout the entire soil profile, with more uptake in deeper soils. The RWUPs estimated by the mixing model were different from those estimated by the mechanistic model and were unrealistic. However, when we constrained the minimum thickness of the water uptake zone, there was good agreement between the RWUPs from the two approaches. We can therefore conclude that the mechanistic model calibrated with isotopes gave better results, and that sole use of the mixing model is not recommended unless appropriate constraints are applied.


INTRODUCTION
The discipline of ecohydrology seeks to determine how hydrological processes influence the distribution, structure, function, and dynamics of biological communities, and how feedback from biological communities affects the water cycle (Newman et al., 2006).At present, water uptake by plant roots may be one of the most important, but least understood, processes in ecohydrology.
We can obtain information about the depth from which root water is taken up by applying either an isotopic tracer approach (reviewed in Ehleringer and Dawson, 1992;Dawson and Ehleringer, 1998) or a mechanistic (or numerical) modelling approach (reviewed in Molz, 1981;Feddes et al., 2001).In the former approach, hydrogen and/or oxygen isotope ratios of xylem sap water are compared with those of soil or ground water at different depths to obtain a rough estimate of the water uptake zone (e.g.shallow or deep water sources).Ogle et al. (2004) developed an inverse estimation algorithm for reconstructing the root water uptake profile (RWUP) from isotope ratios, but the reconstructed RWUP has not been sufficiently validated.The numerical modelling approach allows us to simulate spatial and temporal variations in root water uptake (Campbell, 1985).However, the results rely heavily on an assumed rooting pattern and it is very hard to estimate the vertical distribution of hydrologically-active roots (Dawson et al., 2002).Therefore, there are uncertainties associated with both of these RWUP estimation methods, although unfortunately there is not enough information to permit any assessment of the reliability of the estimates.We therefore need to find ways to improve them.
In this study, we tested two approaches that combined isotopes and models, namely an isotope-calibrated mechanistic model and a mixing model with a bell-shaped approximation.The objectives of this study were (1) to obtain reliable estimates of the RWUP for a conifer plantation, and (2) to assess the usefulness and limitations of the two approaches.

Study site
The study site was a 32-year-old Japanese cypress (Chamaecyparis obtusa) plantation on the southwestern hillslope of Mt.Karasawa in the southern part of Tochigi Prefecture, Central Japan (139°36'E, 36°22'N; 198 m a.m.s.l.).Climate is temperate and humid with an annual mean temperature of 14.1°C and annual total precipitation of 1,265 mm.Snowfall sometimes occurs in winter, though annual maximum snow depths are less than 20 cm.The topsoil is an orthic brown cambisol with a silt loam texture (Sun et al., 2014a).The exact soil layer thickness has not been identified because hard rocks or roots did not allow us to excavate the soil profile >1 m below ground.
The cypress trees were, on average, 16.0 m high and had a mean diameter at breast height (DBH) of 19.1 cm.The plantation had a stand density of 2,198 stems ha -1 .The understory vegetation was sparse and included Gleichenia japonica, Cleyera japonica, and Ardisia japonica.Information about the canopy water budget and components of evapotranspiration at this site were reported by Sun et al. (2014b) and Sun et al. (2014a), respectively, and the vertical root density profile (i.e. the total length of root within a unit volume of the soil) was investigated by Kimura (2012).These investigations were performed on a southwestward slope (31°) approximately 50 m from the nearest ridge.

Isotopic survey
We collected stem core samples from three individual cypress trees using an increment borer on five separate occasions in warm weather between July and September of 2011.The stem core samples were approximately 15 cm long and had a diameter of 0.5 cm.We extracted xylem water from the stem cores by cryogenic vacuum distillation (Iizuka et al., 2004).We collected soil water samples from four depths (10, 20, 40, and 80 cm) at a distance of 1 m from trunks of each of the 3 trees using suction lysimeters (DIK-8390, DAIKI, Japan); on each sampling day we collected a total of 12 soil water samples.We also collected precipitation, mostly weekly but sometimes biweekly or monthly, over a two-year period from July 2010 to July 2012 using two precipitation collectors specially developed to avoid evaporation of stored water (Yamanaka et al., 2005).One collector was installed in an open site near the plantation stand and the other was within the stand and collected throughfall.Some results on isotopic characteristics of open rainfall and throughfall were given by Kato et al. (2013).
Oxygen and hydrogen stable isotope ratios ( 18 O/ 16 O and 2 H/ 1 H) in the collected water samples were measured by tunable-diode laser spectroscopy with a liquid water isotope analyzer (L1102-i, Picarro, CA, USA) at the Center for Research in Isotopes and Environmental Dynamics at the University of Tsukuba.The results are expressed using the δ notation (i.e.δ 18 O and δ 2 H) relative to the Vienna Standard Mean Ocean Water (V-SMOW).The analyzer measurement errors for δ 18 O and δ 2 H were 0.1‰ and 1.0‰, respectively (Yamanaka and Onda, 2011).Since, for unknown reasons, analysis methods tend to underestimate the hydrogen stable isotope ratio of xylem water (Iizuka et al., 2004), we only used the oxygen isotope data in the present study.

Mechanistic model
We developed a mechanistic model called the Root-zone Hydrology coupled with the Ecosystem and Atmosphere (RHEA) model to simulate root water uptake dynamics (Text S1, Tables SI-SIII, and Figures S1-S5).The model computes the energy balance and water partitioning in three layers (canopy, understory, and ground layers), including transpiration from trees, and then calculates the water flux between the soil layers and the root water taken up from them.Data inputs to RHEA include information about daily atmospheric conditions (shortwave/longwave downward radiation fluxes, air temperature, relative humidity, wind speed, air pressure, and precipitation), above-ground vegetation conditions (vegetation height, leaf area index, stomatal properties, and rooting profile), and soil conditions (porosity, heat capacity, thermal conductivity, saturated hydraulic conductivity, and water retention curves) that can be determined from the soil texture.In this study, the soil domain was 1.5 m thick.It was subdivided into 30 layers at 5-cm increments, and the bottom boundary condition was assigned a constant potential, namely the air entry potential.The root density profile (RDP) was expressed by the following exponential function (Gerwitz and Page, 1974): where L r is the root density for a layer (m m -3 ), L r0 is the L r at the top layer (m m -3 ), c Lr is an experimental coefficient, and z is the depth (m).Regression analysis of RDP measured by Kimura (2012) gave L r0 = 5305.5 and c Lr = 6.
To validate the simulation results, isotopic processes (i.e.fractionation, transport, and mixing) were incorporated into a version of the model called isoRHEA.As mentioned earlier, the hydrologically active roots and those physically measured are not always the same.By incorporating isotopes into the calibration, it may be easier to constrain the active rooting pattern, because only the roots that actively take up water can be identified by isotopic tracers (Dawson et al., 2002).Details of formulations in isoRHEA are provided in Text S2.
Numerical simulations were run for a 3-year period from 1 January 2010.We used meteorological data from the Japan Meteorological Agency for Sano (station ID = 41361; air temperature, precipitation, wind speed, and sunshine duration; about 5 km apart from the site) and Utsunomiya (station ID = 41915; relative humidity and air pressure; about 35 km apart from the site).We estimated the downward shortwave/longwave radiation fluxes from these observed values using the method of Allen et al. (1998).

Mixing model
In general, mechanistic models can represent complex systems with many parameters, but the parameters are often highly uncertain and difficult to constrain.Romero-Saltos et al. (2005) proposed a simple approach based on only isotopic data, which was also applied by Li et al. (2007).It is a kind of mixing model that gives rough estimates of the depths at which root water uptake mainly occurs.
In this model, the xylem water is considered a mixture of T. YAMANAKA ET AL.
soil water at different depths, and the xylem water isotope delta (δ x ) can be given as: where i is the layer number, N is the total number of layers, δ sw (i) is the isotope delta (‰) of soil water at layer i, and f u (i) is the fraction of root water uptake at layer i relative to the total uptake (i.e.transpiration).
The RWUP or vertical profile of f u (i) is assumed to be bell-shaped, that is, approximated by a normal distribution function defined by only two parameters, the mean (μ) and the standard deviation (σ).Although the normal distribution function has no limit, the soil segment is limited in this study, so that f u (i) is calculated as follows: where z(i) is the depth (m) of the i-th layer, and z m is the maximum depth of the soil segment, which was 1.5 m in this study.For this reason, the μ value is not always the mean depth of the water uptake zone and σ is not strictly the standard deviation.Interpolation of the observed δ sw into each layer allows us to calculate δ x with Equation (1) using arbitrary sets of μ and σ.Out of all possible combinations of μ and σ, a combination that gives the minimum difference between the estimated and measured δ x can be considered the most likely.In this study, we scanned values of μ and σ between 0 and 1.5 m at a resolution of 0.05 m and then obtained the optimum values using an optimization procedure known as systematic exploration of the inverse problem (Tarantola, 1987).
This mixing model does not involve any of the physical root water uptake processes, so it is difficult to precisely explain why the RWUP is bell-shaped.However, it is easier to apply than the mechanistic model and it reflects isotopic evidence about the root water uptake zone.Also, the bellshape approximation can help us identify the peak depth and vertical range of the uptake zone.

Mechanistic model
Figure 1 shows the day-to-day variations in air temperature, precipitation, and model outputs, such as the net radiation and transpiration fluxes from the canopy and understory.The modeled transpiration from the canopy was up to ten times greater than that from understory, and the deviations were more noticeable during the summer season.Estimates of seasonal variation in canopy transpiration were intermediate and fell between those of air temperature and net radiation.Both transpiration fluxes dropped slightly in response to rainfall events, when incoming radiative energy was mainly consumed by evaporation of the water intercepted in the canopy.
In July, August, and September of 2011, there was good agreement between the observed and modelled mean monthly ratios of interception/precipitation (I/P), transpiration/precipitation (T/P), soil evaporation/precipitation (E s /P), and total evapotranspiration/precipitation (ET/P) (Figure 2).Modeled transpiration (including transpiration from the canopy and understory) was between 0.04 and 0.24 mm d -1 higher than the observed transpiration, and corresponded to between 10% and 40% of the monthly precipitation, respectively, with the highest percentages observed in August.

ISOTOPE-CALIBRATED ROOT WATER UPTAKE PROFILE
When calculated on an annual basis, modeled transpiration for 2011 and 2012 accounted for 32.5% and 36.6% of annual total precipitation, and 67.1% and 82.2% of evapotranspiration, respectively.The interception loss was almost the same as or greater than transpiration, and evaporation from soil accounted for less than 0.3% and 0.4% of precipitation in 2011 and 2012, respectively.
Using the measured RDP (i.e.c Lr = 6), the model values of δ x in July and August were overestimated relative to the measured values (Figure 3).For a gentler RDP gradient and smaller values of c Lr , there was less variation in δ x .As shown in Figure 3, δ x was best reproduced when c Lr = 1.These results indicate that the RDP measured by direct sampling at a distance from the tree trunk (i.e. 1 m) cannot be used to represent the hydrologically active RDP.
The simulated RWUP depended heavily on the RDP (Figure 4).Estimates based on the measured RDP showed that there was water uptake in shallower soils.However, when c Lr was 1, there was water uptake through the entire soil profile, with more uptake from deeper soils.In either case, there was little temporal variation in the RWUP.Undoubtedly, RWUP estimated with calibrated RDP must be more valid because of better reproducibility of δ x .
Because some vertical profiles of δ sw showed neither monotonic increases nor decreases, and spatial heterogeneity was relatively large especially in shallow zones (Figure 5), it is difficult to directly determine the depth of water uptake based on these profiles.However, mean values and small standard deviations for δ x are similar to those for δ sw in a deeper zone, suggesting that deep water uptake dominates.

Mixing model
Figure 6 shows the RWUP estimated by the mixing model.Without any constraints, the model predicted considerable changes in the RWUP; deep uptake dominated on 26 July and 11 Aug., but did not dominate on 30 Aug. and 29 Sep.These large changes in the RWUP are problematic because the RWUP derived from the mechanistic model showed no such tendencies to change (Figure 4).As shown or not monotonic (26 July, 30 Aug., and 29 Sep.).In such cases, when the σ is small, the RWUP is highly sensitive to noise errors in the measured δ sw .
In our trials, we constrained σ to 0.5 m or 1.0 m.When σ was greater than 0.5 m, the day-to-day variations in RWUP were not so drastic and deep water uptake dominated on most days.When σ was greater than 1.0 m, the patterns were clearer and agreed well with those from the mechanistic model.These results suggest that a bell-shaped RWUP can be a good approximation of the actual RWUP, although the mixing model cannot provide reliable estimates of the RWUP unless appropriate constraints on σ are given.

DISCUSSION AND CONCLUSIONS
Using seasonal variations in deuterium excess for xylem water, soil water and groundwater, Kubota et al. (2012) showed that 20-year-old Japanese cypress took up water from shallow depths (< 1.0 m below ground) rather than deeper depths (> 1.5 m).This is not always contradictory to our results, because we focused on the relative proportion of shallow/deep water uptake within a 1.5 m-thick soil layer.Unfortunately, no detailed information on RDP and hydraulic properties were available for the soil (or bedrock) deeper than 1.5 m, so that we could not evaluate water uptake from such a very deep zone.In addition, Kubota et al. (2012) could not clarify RWUP for the soil down to 1.0 m.Therefore, we cannot make a meaningful comparison of their results with our results, though the effects of soil thickness and/or stand age on RWUP will be an interesting issue.
Our results indicate that the mechanistic model calibrated with isotopic data provides the best method for estimating RWUPs.The isotope-calibrated model predicted water uptake throughout the entire soil profile and, even though the measured and calibrated root densities were both higher in shallower soils, uptake was greatest from deeper soils.Various previous studies reported such a dissimilarity between RDP and RWUP (Campbell, 1985;Hamblin and Tennant, 1987).The reduction in the water uptake in shallower soils can be explained by increased resistance in the soil matrix that impeded the flow of water from the soil matrix to the root surface (Equation (S1-5) of Text S1) and the root surface that impeded the flow of water across the root surface (Equation (S1-6) of Text S1) when the soil dried out.
Water isotopes are therefore powerful and essential calibration tools for mechanistic models.To improve the method, we need to address the shape of the hydrologically active RDP.In this study, we expressed the RDP with an exponential function as a first approximation, as has been done in previous studies (Campbell, 1985;Li et al., 1999Li et al., , 2001)).However, if the shape of the RDP was more complicated, we could obtain different RWUPs and similar δ x .In recent decades, root growth and root system architecture have been simulated by numerical models (Somma et al.,Figure 6. RWUPs estimated by the mixing model with different constraints.Estimated results from the mechanistic model after isotope calibration are also shown 1998; Pagès et al., 2004), and RWUPs have been estimated from the simulated results (Javaux et al., 2008;Schneider et al., 2010).We may be able to combine this approach with our method for complicated RDPs.
It should be noted that, when the RDP was calibrated with δ x , there were no significant changes in either the energy balance or water partitioning at the land-atmosphere interface, because stomatal conductance was assumed to be independent of water uptake dynamics in the present model.Because there was good agreement between the estimated and observed transpiration fluxes (Figure 2), this seems to be a minor issue for temperate forests with humid climates.Although this should be examined further and the model may need to be modified if it is applied to arid and semi-arid regions, the isotope-calibrated root-zonehydrology model has a great potential in tracing the terrestrial water cycle by coupling with isotope-incorporated runoff models (Ma andYamanaka, 2013, 2016;Yamanaka and Ma, 2017) or groundwater flow models (Liu et al., 2014).
Finally, we need to admit that the mixing model with the bell-shaped approximation is not always capable of reproducing the actual RWUP.In particular, when vertical changes in δ sw are small or nonmonotonic, the estimated RWUP should be checked with those values derived from other approaches.We do not recommend the sole use of the mixing model approach unless appropriate constraints are applied.

Figure 1 .
Figure 1.Day-do-day variations in air temperature and precipitation as model input, and estimated net radiation and transpiration fluxes from the canopy and understory

Figure 2 .
Figure 2. Comparison of observed and modeled monthly mean flux ratios; I is the interception, P is the precipitation,T is the transpiration, E s is the soil evaporation, and ET is the total evapotranspiration.Observed values are adopted fromSun et al. (2014a, b)

Figure 3 .
Figure 3. Day-to-day variations in measured and modeled δ 18 O of xylem water.Model estimates were given with different c Lr values, which is a parameter relating to curvature of exponential root density profile (larger value represents steeper slope at shallow soils)

Figure 4 .
Figure 4. Root water uptake profiles (RWUPs) estimated by the mechanistic model with different c Lr Text S1.Detailed description of the RHEA model Text S2.Detailed description of the isoRHEA model Figure S1.Schematic illustration of the soil sub-model Figure S2.Short-wave radiation transfer in the canopy submodel Figure S3.Long-wave radiation transfer in the canopy submodel Figure S4.Resistance network for sensible heat in the canopy sub-model Figure S5.Resistance network for latent heat in the canopy sub-model

Table SI .
List of symbolsTable SII.Formulation and parameter setting in the soil sub-model Table SIII.Formulation and parameter setting in the canopy sub-model