Coordination of rooting, xylem, and stomatal strategies explains the response of conifer forest stands to multi-year drought in the southern Sierra Nevada of California

Extreme droughts are a major determinant of ecosystem disturbance that impacts plant communities and feeds back into climate change through changes in plant functioning. However, the complex relationships between aboveground and belowground plant hydraulic traits and their role in governing plant responses to drought are not fully understood. In this study, we use a model, the Functionally Assembled Terrestrial Ecosystem Simulator in a configuration that includes plant hydraulics (FATES-Hydro), to investigate ecosystem responses to the 2012–2015 California drought in comparison with observations at a site in the southern Sierra Nevada that experienced widespread tree mortality during this drought.

We conduct a sensitivity analysis to explore how different plant water sourcing and hydraulic strategies lead to differential responses during normal and drought conditions.

The analysis shows the following.

1. Deep roots that sustain productivity through the dry season are needed for the model to capture observed seasonal cycles of evapotranspiration (ET) and gross primary productivity (GPP) in normal years, and deep-rooted strategies are nonetheless subject to large reductions in ET and GPP when the deep soil reservoir is depleted during extreme droughts, in agreement with observations.

2. Risky stomatal strategies lead to greater productivity during normal years as compared to safer stomatal control, but they also lead to a high risk of xylem embolism during the 2012–2015 drought.

3. For a given stand density, stomatal and xylem traits have a stronger impact on plant water status than on ecosystem-level fluxes.

Our study highlights the significance of resolving plant water sourcing strategies to represent drought impacts on plants and consequent feedbacks in models.

1 Introduction

Understanding plant water use strategies and the resulting ecohydrologic processes in forests is critical for predicting surface water and energy exchange, carbon dynamics, and vegetation dynamics of water-constrained ecosystems in a changing climate. Mediterranean-type climates, as in California, are characterized by hot, dry summers and cool, wet winters, resulting in asynchronous supplies of energy and water. In addition to these climatic stresses, plants in California are further subject to high interannual variability in precipitation and periodic severe drought events, such as the recent 2012–2015 drought, which led to widespread tree mortality (Fettig et al., 2019). Together, these two climatic constraints present a unique challenge to the success of forests in California, which is likely to be exacerbated by a warming climate.

On evolutionary timescales, natural selection has led to a wide array of strategies and functional traits that allow plants to both grow and survive under a range of environmental conditions (Grime, 1977, 1979; Coley et al., 1985; Westoby et al., 2002; Craine, 2002; Reich et al., 2003). Given the centrality of water sourcing to plant physiology, plant hydraulic traits play an important role in water-constrained ecosystems. Once absorbed by fine roots, water flows through the vascular system via coarse roots, stems, and branches to leaves, where it evaporates through stomata. The rate of water flow through stems, and thus the supply to leaves, is determined by the hydraulic conductivity along this pathway. If the water potential of xylem tissue becomes too low, cavitation can occur and cause a loss of conductivity. Because this cavitation can damage the xylem network, trees have developed different strategies to mitigate this effect, all of which come at some cost. These strategies include (1) early stomatal closure or leaf deciduousness to reduce the flow of water at the cost of reduced carbon intake, (2) building cavitation-resistant xylem at the cost of increased hydraulic resistance, and (3) growing deep roots to access more moisture at the cost of higher carbon investment. In this study, we focus on the potential hydraulic strategies that trees in Californian ecosystems use, with a particular emphasis on how the coordination of hydraulic functional traits at the leaf, stem, and root levels is critical to carbon assimilation, transpiration, and consequently the productivity and response of trees to drought (Matheny et al., 2017a, b; Mursinna et al., 2018a).

The traits that regulate stomatal conductivity are the most important hydraulic traits of leaves and the primary ones through which photosynthesis and transpiration are coupled. Stomatal behavior falls along a gradient between two extremes: stomata may close early during water stress to avoid the risk of hydraulic failure or remain open to maximize carbon uptake while exposing xylem to a higher risk of embolism (Martínez-Vilalta et al., 2004; McDowell et al., 2008; Skelton et al., 2015; Matheny et al., 2017). The sensitivity of stomata to water stress determines where the stomata operate along the safety–risk gradient and thus the degree to which carbon intake is traded to prevent the cavitation of xylem. Where the best stomatal strategy sits along the safety–risk gradient would depend on the physical environment.

Maximum hydraulic conductivity and vulnerability to cavitation are the two key xylem hydraulic traits. Differences in the anatomy and morphology of the conductive xylem cell structure and anatomy (Hacke et al., 2017) lead to differences in maximum conductivity and the water potential at which cavitation starts to occur (Pockman and Sperry, 2000; Sperry, 2003). Within the conifers, there are at least three mechanisms that lead to a tradeoff between xylem safety and efficiency. The first is the morphology of the xylem conduit. It is widely acknowledged that narrow (or short) tracheids are safer than wider (or longer) tracheids but have a lower conductance per sap area (Choat and Pittermann, 2009). Second are the inter-vessel pit membranes. Thicker and less porous membranes prevent the spread of air but increase the hydraulic resistance of xylem (e.g., Li et al., 2016; Pratt and Jacobsen, 2017). The third mechanism comes from the division of limited space (Pratt and Jacobsen, 2017). With the same cross-sectional area of conduits, vessels with a thicker cell wall provide stronger mechanical support, so that the conduits are less likely to collapse when xylem water potential becomes more negative. However, this reduces the area that can be used for conduits transporting water. While these physiological constraints require that the tradeoff exists to some extent, in many studies, this tradeoff appears to be weak, and there are certainly species that have both safe and efficient xylem. Further, there are many other plant traits that can affect safety, such as wood density (Pratt and Jacobsen, 2017), pit anatomy (Sperry and Hacke, 2004; Lens et al., 2011), and biochemistry (Gortan et al., 2011). These traits can have large variations among different plant types. The tradeoff will be weakened when grouping plants at a coarse scale, e.g., by biomass, families, and/or a range of geological and climatic regions. However, when focusing on certain species in a particular region, the tradeoff becomes stronger, as demonstrated by many local studies (e.g., Barnard et al., 2011; Corcuera et al., 2011; Baker et al., 2019). For example, Kilgore et al. (2021) show that there is a clear safety–efficiency tradeoff across pine trees in a specific location in the western US. Thus, while we acknowledge that there are many exceptions to the xylem safety–efficiency tradeoff, it is a useful framework for examining plant strategies for dealing with drought.

The traits that govern the hydraulic function of plant root systems are also critically important, but they are the least understood, studied, and quantified. These traits include the rooting depth, the root-to-shoot ratio, the vertical and lateral distributions of roots, and the fine root density and diameter, all of which are related to water uptake (Canadell et al., 2007; Allen, 2009; Reichstein et al., 2014; Wullschleger et al., 2014). In general, species with deeper roots can access water at greater depths, which is unavailable to more shallow-rooted species (Jackson et al., 1996; Canadell et al., 1996). The vertical root distribution can affect the water uptake and thus the evapotranspiration (ET) pattern during the dry-down period (Teuling et al., 2006). This in turn affects the seasonal distribution of water over the soil depth and thereby the resilience of plants to seasonal droughts (Yu et al., 2007). The vertical root distribution is also a means of belowground niche differentiation (Ivanov et al., 2012; Kulmatiski and Beard, 2013), whereas the extent of the lateral root distribution dictates the competition for water (Agee et al., 2021). Whether a plant can benefit from having deep roots is related to the plant's leaf and xylem hydraulic traits (e.g., Johnson et al., 2018; Mackay et al., 2020), thus requiring coordination of rooting and hydraulic traits.

Given the strength of the Mediterranean-type climate of California, the coordination of rooting and hydraulic strategies will play a critical role in forest dynamics. However, the interplay of rooting and hydraulic strategies and their impact on ecosystem processes have not been well understood. In this study, we address this question at the Soaproot site (CZ2) of the southern Sierra Nevada of California as the study area. The CZ2 site was strongly affected by the 2012–2015 drought, with extremely high tree mortality rates (90 % of the pine died) (Fettig et al., 2019). While the 2012–2015 drought was widespread across California, the highest rates of tree mortality occurred in the southern Sierra Nevada, centered around an elevation similar to this site (1160 to 2015 m; Asner et al., 2016; Goulden and Bales, 2019). This mid-elevation region is also characterized by the highest forest productivity along an elevation gradient from foothill woodlands to subalpine forest (Kelly and Goulden, 2016). This leads us to ask whether strategies associated with high productivity have exposed trees to a high mortality risk under prolonged drought.

Specifically, here we use the Functionally Assembled Terrestrial Ecosystem Simulator in a configuration that includes plant hydraulics (FATES-Hydro) to explore the tradeoffs associated with differing hydraulic strategies and, in particular, their implications for plant productivity and the risk of drought-induced mortality. We conduct a sensitivity analysis using FATES-Hydro in comparison with observations from the CZ2 eddy-covariance site to investigate how stomatal, xylem, and rooting strategies affect the ecosystem and physiological processes of the forest and whether that may explain the high rates of both productivity and drought-associated mortality of conifers at CZ2. We note that this is not an exhaustive model parameter sensitivity study. The main purpose is to use a sensitivity analysis to explore scientific questions around hydraulic trait tradeoffs.

2 Methods

2.1 Study site

The Soaproot site is a 543 ha headwater catchment at 1100 m elevation (37^∘2.4^′ N, 119^∘15.42^′ W), which is at the lower boundary of the rain–snow transition line with warm, dry summers and cool, wet winters (Geen et al., 2018). The mean annual temperature is about 13.8 ^∘C (Goulden et al., 2012). Under normal conditions, the annual precipitation is about 1300 mm, but during a dry year, the precipitation can drop to 300–600 mm (Bales et al., 2018a). The site is a ponderosa pine (Pinus ponderosa)-dominated conifer ecosystem exhibiting high productivity. Kelly and Goulden (2016) reported 2.1 t C ha⁻¹ yr⁻¹ average annual gross stem wood production. Other species include California black oak (Quercus kelloggii Newberry) and incense cedar (Calocedrus decurrens).

Soils at the Soaproot site are mainly of the Holland (fine-loamy, mesic Ultic Haploxeralfs) and Chaix (coarse-loamy, mesic Typic Dystroxerepts) series, which are representative of soils across a similar elevation band of the western Sierra Nevada (Mooney and Zavaleta, 2003). Soils of the Holland series have sandy loam surface textures and underlying Bt horizons with sandy clay loam textures, while soils of the Chaix series have sandy loam textures throughout the profile. The regolith depth is estimated to be 15 m (Holbrook et al., 2014). The total porosity over the whole regolith depth of the site is estimated to be 1620 mm, and the total available storage porosity (plant-accessible water storage capacity), which is the difference in the volumetric water content between the field capacity and the permanent wilting point ($\sim -\mathrm{6}$ Mpa), is estimated to be 1400 mm (Klos et al., 2017). The available water storage capacity is approximately 0.20 cm³ cm⁻³ in the upper regolith (0–5 m depth), which decreases to 0.05 cm³ cm⁻³ or less in the lower regolith (below 5 m depth) (Holbrook et al., 2014).

An eddy-covariance flux tower was installed at this site in September 2010. The elevation of the tower is 1160 m a.s.l. Instruments on the flux tower track changes in carbon dioxide, water vapor, air temperature, relative humidity, and other atmospheric properties. We compare the simulated gross primary productivity (GPP) and latent heat flux with the flux tower measurements over the period from 2011 to 2015 (Goulden and Bales, 2019). We computed the root mean square error (RMSE) of the hourly mean diurnal cycle of each month. This allows us to examine the capacity of FATES-Hydro to predict the carbon and water fluxes. The transpiration at the site contributed to the majority of the ET, as indicated by the measurements from an adjacent catchment and the fact that the site is fully vegetated with an annual leaf area index (LAI) of around 3 to 4.

2.2 FATES-Hydro model and parameterization

2.2.1 The FATES-Hydro model

FATES is a cohort-based and size- and age-structured dynamic vegetation model, where long-term plant growth and mortality rates and plant competition emerge as a consequence of physiological processes. In the model, multiple cohorts grow on the same land unit, share the soil water, and interact with each other through light competition. FATES is coupled within both the Community Land Model V5 (CLM5) (Lawrence et al., 2019) and Energy Exascale Earth System Land Model (ELM) (Golaz et al., 2020) land surface models (LSMs). In this study, FATES is coupled with the CLM5 model. FATES-Hydro is a recent development of the FATES model (Fisher et al., 2015; Koven et al., 2020), in which a plant hydrodynamic module, originally developed by Christoffersen et al. (2016), is coupled to the existing photosynthesis and soil hydraulic modules. FATES-Hydro is described in more detail by Xu et al. (2023) and its Supplement.

Conceptually, plant hydraulic models can be broadly grouped into two types. The first group represents the plant hydraulic system as analogous to an electrical circuit (e.g., Mackay et al., 2011; Huang et al., 2017; Eller et al. 2018; Kennedy et al., 2019). The total resistance of the plant is calculated from the resistance of each compartment using Ohm's law. There is no storage of water in the plants, and the transpiration from plants at any given time step is considered to come directly from soil storage. The second group represents plant hydraulics by a series of connected porous media corresponding to each plant compartment (e.g., Bohrer et al., 2005; Janott et al., 2011; Xu et al., 2016; Christoffersen et al., 2016). The porous media model takes into account the water storage in the plant. The flow between two adjacent compartments is driven by the difference in the water potential mediated by the hydraulic conductivity. FATES-Hydro falls into the second group. The various models in the second group differ in the exact formulas used to describe the pressure–volume and pressure–conductivity relations as well as the different numbers and arrangements of nodes within the soil–plant–atmosphere system.

In FATES-Hydro, for each plant cohort, the hydraulic module tracks water flow along a soil–plant–atmosphere continuum of a representative individual tree based on hydraulic laws and updates the water content and potential of leaves, stems, and roots with a 30 min model time step. Water flow from each soil layer within the root zone into the plant root system is calculated as a function of the hydraulic conductivity as determined by root biomass and root traits such as specific root length and the difference in water potential between the absorbing roots and the rhizosphere. The vertical root distribution is based on the Zeng (2001) two-parameter power-law function, which takes into account the regolith depth:

$$\begin{array}{}\text{(1)}& Y_i={\displaystyle \frac{\mathrm{0.5}(e^{-r_{\mathrm{a}}{z}_{\text{l}i}}+e^{-r_{\mathrm{b}}{z}_{\text{l}i}})-\mathrm{0.5}(e^{-r_{\mathrm{a}}{z}_{\text{u}i}}+e^{-r_{\mathrm{b}}{z}_{\text{u}i}})}{\mathrm{1}-\mathrm{0.5}(e^{-r_{\mathrm{a}}z}+e^{-r_{\mathrm{b}}z})}},\end{array}$$

where Y_i is the fraction of fine or coarse roots in the ith soil layer, r_a and r_b are the two parameters that determine the vertical root distribution, Z_li is the depth of the lower boundary of the ith soil layer, Z_ui is the depth of the upper boundary of the ith soil layer, and Z is the total regolith depth. The vertical root distribution affects water uptake by the hydrodynamic model by distributing the total amount of root, and thus root resistance, through the soils.

The total transpiration of a tree is the product of the total leaf area (LA) and the transpiration rate per unit leaf area (J). In this version of FATES-Hydro, we adopt the model developed by Vesala et al. (2017) to take into account the effect of leaf water potential on the within-leaf relative humidity and transpiration rate:

$$\begin{array}{}\text{(2a)}& {\displaystyle }E=\text{LA}\cdot J,\text{(2b)}& {\displaystyle }J={\mathit{\rho }}_{\text{atm}}{\displaystyle \frac{(q_{\text{l}}-q_{\text{s}})}{\mathrm{1}/g_{\text{s}}+r_{\mathrm{b}}}},\text{(2c)}& {\displaystyle }{q}_{\text{l}}=\exp \left({\displaystyle \frac{k_{\text{LWP}}\cdot \text{LWP}\cdot V_{{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}}{R\cdot T}}\right)\cdot q_{\text{sat}},\end{array}$$

where E is the total transpiration of a tree; LA is the total leaf area (m²); J is the transpiration per unit leaf area (kg s⁻¹ m⁻²); ρ_atm is the density of atmospheric air (kg m⁻³); q_l is the within-leaf specific humidity (kg kg⁻¹); q_s is the atmosphere specific humidity (kg kg⁻¹); g_s is the stomatal conductance per leaf area, r_b is the leaf boundary layer resistance (s m⁻¹); w is a scaling coefficient (unitless), which can vary between 1 and 7, and here we use a value of 3; LWP is the leaf water potential (Mpa); $V_{{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}$ is the molar volume of water ($\mathrm{18}\times {\mathrm{10}}^{-\mathrm{6}}$ m³ mol⁻¹); R is the universal gas constant; and T is the leaf temperature (K).

The sap flow from absorbing roots to the canopy through each compartment of the tree along the flow pathway (absorbing roots, transport roots, stems, and leaves) is computed according to Darcy's law in terms of the plant sapwood water conductance and the water potential gradient:

$$\begin{array}{}\text{(3)}& Q_i=-K_i\left[{\mathit{\rho }}_{\text{w}}g\right(z_i-z_{i+\mathrm{1}})+({\mathrm{\Psi }}_i-{\mathrm{\Psi }}_{i+\mathrm{1}}\left)\right],\end{array}$$

where ρ_w is the density of water; z_i is the height of the compartment (m); z_i+1 is the height of the next compartment down the flow path (m); Ψ_i is the water potential of the compartment (Mpa); Ψ_i+1 is the water potential of the next compartment down the flow path (Mpa); and K_i is the hydraulic conductivity of the compartment (kg Mpa⁻¹ m⁻¹ s⁻¹). The hydraulic conductivity of the compartments is by the water potential and maximum hydraulic conductivity of the compartment through the pressure–volume (P–V) curve and the vulnerability curve (Manzoni et al., 2013; Christoffersen et al., 2016).

The plant hydrodynamic representation and the numerical solver scheme within FATES-Hydro follow Christoffersen et al. (2016). We made a few modifications to accommodate the multiple soil layers and to improve the numerical stability. First, to accommodate the multiple soil layers, we have sequentially solved the Richards equation for each individual soil layer, with each layer-specific solution proportional to each layer's contribution to the total root–soil conductance. Second, to improve the numerical stability, we have the option to linearly extrapolate the P–V curve beyond the residual and saturated tissue water content to avoid the rare cases of overshooting in the numerical scheme under very dry or wet conditions. Third, Christoffersen et al. (2016) use three phases to describe the P–V curves: (1) dehydration phases representing capillary water (sapwood only), (2) elastic cell drainage (positive turgor), and (3) continued drainage after cells have lost turgor. Due to the possible discontinuity of the curve between these three phases, this leads to the potential for numerical instability. To resolve this instability, FATES-Hydro added the van Genuchten model (Van Genuchten, 1980; July and Horton, 2004) and the Campbell model (Campbell, 1974) as alternatives to describe the P–V curves.

In this study, we use the van Genuchten model because of two advantages: (1) it is simple, with only three parameters needed for both curves; and (2) it is mechanistically based, with both the P–V curve and the vulnerability curve derived from a pipe model and thus connected through three shared parameters:

$$\begin{array}{}\text{(4a)}& {\displaystyle }\mathrm{\Psi }={\displaystyle \frac{\mathrm{1}}{-\mathit{\alpha }}}\cdot {\left({\displaystyle \frac{\mathrm{1}}{S_{\text{e}}^{\mathrm{1}/m}}}-\mathrm{1}\right)}^{\mathrm{1}/n},\text{(4b)}& {\displaystyle }\text{FMC}={\left(\mathrm{1}-{\left({\displaystyle \frac{{\left(-\mathit{\alpha }\cdot \mathrm{\Psi }\right)}^n}{\mathrm{1}+{\left(-\mathit{\alpha }\cdot \mathrm{\Psi }\right)}^n}}\right)}^m\right)}^{\mathrm{2}},\end{array}$$

where Ψ is the water potential of the media (xylem in this case) (Mpa); FMC is the fraction of the xylem conductivity, $K/K_{\text{max}}$ (unitless); α is a scaling parameter for the air entry point (Mpa⁻¹); S_e is the dimensionless standardized relative water content expressed as $S_{\text{e}}=(\mathit{\theta }-{\mathit{\theta }}_{\text{r}})/({\mathit{\theta }}_{\text{sat}}-{\mathit{\theta }}_r)$, with θ, θ_r, and θ_sat the volumetric water content (m³ m⁻³), residual volumetric water content, and saturated volumetric water content, respectively; and m and n are dimensionless (xylem conduit) size distribution parameters. The model assumes that xylem conductance can be restored as xylem water content increases due to increased water availability after a dry period without any hysteresis in the FMC curve.

The stomatal conductance is modeled in the form of the Ball–Berry conductance model (Ball et al., 1987; Oleson et al., 2013; Fisher et al., 2015):

$$\begin{array}{}\text{(5)}& g_{\text{s}}=b_{\text{slp}}{\displaystyle \frac{A_{\text{n}}}{c_{\text{s}}/P_{\text{atm}}}}{\displaystyle \frac{e_{\text{s}}}{e_i}}+b_{\text{opt}}{\mathit{\beta }}_{\text{t}},\end{array}$$

where b_slp and b_opt are parameters that represent the slope and intercept in the Ball–Berry model, respectively. These terms are plant-strategy-dependent and can vary widely with plant functional types (Medlyn et al., 2011). The parameter b_opt is also scaled by the water stress index β_t. A_n is the net carbon assimilation rate (µmol CO₂ m⁻² s⁻¹) based on the Farquhar (1980) formula. This term is also constrained by the water stress index β_t in the way that V_cmax,25 is scaled by β_t as V_cmax,25β_t (Fisher et al., 2018). c_s is the CO₂ partial pressure at the leaf surface (P_a), e_s is the vapor pressure at the leaf surface (P_a), and e_i is the saturation vapor pressure (P_a) inside the leaf at a given vegetation temperature when A_n=0.

The water stress index β_t, a proxy for stomatal closure in response to desiccation, is determined by the leaf water potential adopted from the FMC_gs term from Christoffersen et al. (2016):

$$\begin{array}{}\text{(6)}& {\mathit{\beta }}_{\text{t}}={\left[\mathrm{1}+{\left({\displaystyle \frac{{\mathrm{\Psi }}_{\text{l}}}{P_{\mathrm{50}\mathrm{gs}}}}\right)}^{\text{ags}}\right]}^{-\mathrm{1}},\end{array}$$

where Ψ_l is the leaf water potential (MPa), P_50 gs is the leaf water potential of 50 % stomatal closure, and a_gs governs the steepness of the function. For a given value of a_gs, the P_50 gs controls the degree of risk of xylem embolism (Christoffersen et al., 2016; Powell et al., 2017). A more negative P_50 gs means that, during leaf dry-down from full turgor, the stomatal aperture stays open and thus allows the transpiration rate to remain high and xylem to dry out, which can thus maintain high photosynthetic rates at the risk of exposing xylem to embolism and thus plant mortality. Conversely, a plant with a less negative P_50 gs will close its stomata quickly during leaf dry-down, thus limiting transpiration and the risk of xylem embolism and mortality associated with it, at the cost of reduced photosynthesis.

Figure 1 Sensitivity analysis setup for (a) root parameters that give five root distribution scenarios with effective rooting depths of 1, 3, 5, 6.5, and 8 m and (b) two xylem scenarios for safe xylem ($P_{\mathrm{50}\mathrm{x}}=-\mathrm{4.8}$, K_max=0.64) and efficient xylem ($P_{\mathrm{50}\mathrm{x}}=-\mathrm{2.5}$, K_max=0.88).

Table 1 Parameters used in the FATES-Hydro sensitivity analysis.

^a Values for efficient or unsafe xylem. ^b Values for inefficient or safe xylem.

2.2.2 Sensitivity analysis and parameterization

The goal of this analysis is to better understand how coordinated aboveground and belowground hydraulic traits determine plant physiological dynamics and the interplay between ecosystem fluxes and tissue moisture during the extreme 2012–2015 drought at the Soaproot site. We thus conduct a global sensitivity analysis of selected hydraulic parameters to explore the linkages of aboveground and belowground hydraulic strategies. We use a full-factorial design for the parameter sensitivity analysis in order to best investigate the relationships between the parameters. Because this design requires a relatively small set of parameters or groups of parameters to vary, we chose parameters that represent the major axes of relatively well-understood stomatal, xylem, and rooting mechanisms or strategies that control the hydraulic functioning of trees. We set the values of these parameters within the realistic (allowable biological) range based on an online database and the literature where the species and physical environment are as close to our system as possible. We list other major parameters and their estimates that do not vary in the sensitivity analysis (Table 2). We acknowledge that the biggest disadvantage of this study is the lack of sufficient field data to constrain the model. This is a result of using a natural drought as an experiment of opportunity, which, because it was not anticipated, did not allow for coordinated planning, as would be the case in an experimentally manipulated drought. The trees at that site had all died by the time we started this study.

The parameters that we vary here are (1) the pair r_a and r_b, which control the vertical root distribution as deep vs. shallow roots; (2) two sets of xylem parameters (P_50x,K_max, m, n, and α) that jointly represent two distinct xylem strategies: efficient or unsafe xylem and inefficient or safe xylem within the range observed for temperate conifer trees; and (3) the stomatal parameter P_50 gs, which represents the stomatal strategy along a risky-to-safe gradient (Table 1). The ranges of the root parameters are chosen so that the effective rooting depth, above which 95 % of the root biomass stays, varies from 1 to 8 m, which is the possible range at the Soaproot site, as indicated by current knowledge of the subsurface structure (Klos et al., 2017). Note that here we refer to a higher proportion of roots in the deep subsurface layers as “deep rooting” (e.g., effective rooting depth = 8 m; r_a=0.1, r_b=0.1) as compared to “shallow rooting” (e.g., effective rooting depth = 2; r_a=1, r_b=5), which represents a larger proportion of fine roots in the upper layers (Fig. 1a).

Table 2 List of major parameters.

The safety–efficiency tradeoff of xylem has been widely discussed in the literature (e.g., Gleason et al., 2016; Hacke et al., 2006, 2017; Martnez-Vilalta et al., 2004). Given that we do not have any measurements that can be used to generate a vulnerability curve at our study site, we consulted the literature (Domec et al., 2004; Barnard et al., 2011; Corcuera et al., 2011; Anderegg and Hillerislambers, 2016; Baker et al., 2019; Kilgore et al., 2021) for observed curves from sites that are as similar in both climate (e.g., mean annual precipitation and temperature) and the set of conifer species (P. ponderosa) to our study site as possible, together with values of xylem traits (K_max and P_50x) of P. ponderosa in temperate regions of the Plant Trait Database (Kattge et al., 2020), to determine the two hypothetical vulnerability curves representing the safe–inefficient and unsafe–efficient xylem strategies. We set the parameters of the van Genuchten model to represent these two sets of the P–V and vulnerability curves as shown in Fig. 1b and c. It is worth noting that, with the same K_max and P₅₀, the exact shape of the vulnerability can differ depending on the formula used and the parameter values. However, this should not be an issue in our study because the vulnerability curve is mainly constrained by P₅₀ and K_max. Second, given that there is a large range of variation in the measured values, the effect caused by the exact shape of the curves is minor. Third, since the objective of our study is not to accurately predict mortality but rather to examine the effect of different combinations of stoma, xylem, and root strategies, even if the shape of our vulnerability curve is not the most accurate, as long as the curve captures the overall pattern of the pressure–conductivity relation, it will not affect the relative outcome of this study.

We follow the theory of Skelton et al. (2015) to define safe vs. efficient stomatal strategies. In FATES-Hydro, there are two key stomatal parameters: P_50 gs and a_gs. Here, we only vary P_50 gs while keeping a_gs constant because the objective here is to choose the parameters that are relatively well understood and to catch the safe vs. risky strategies as described by Skelton et al. (2015) rather than exhaust the parameter space within the model. In essence, the different combinations of P_50 gs and the shape parameter (a_gs) can generate similar stomatal response curves. For example, a small negative P_50 gs with a small a_gs would result in a flat stomatal response curve, which is similar to a large negative P_50 gs combined with a large a_gs. Further, P_50 gs is well understood and has more observed data, while a_gs is less studied and barely has any observed data. With a given a_gs, the variance of P_50 gs for a given P_xylem value controls the degree of embolism risk from a “risky” strategy, where P_50 gs is equal to or lower than P_xylem, to a “conservative” strategy, where P_50 gs is higher than P_xylem. The P_xylem in Skelton et al. (2015) is for a fynbos species and therefore is not appropriate for our study because our species are pine trees, which are woody plants. Trees have woody tissue, which contributes to strengthening the conduits and which makes them less likely to collapse when embolized, hence allowing their stomata to be riskier than those of herbaceous plants. From the observed P_50 gs and xylem traits of the closely related pine species in the TRY database (Kattge et al., 2020) and elsewhere in the literature (Bartlett et al., 2016) as well as the observed soil water potential at the study site, we choose to vary P_50 gs between P_50 xylem and P_20 xylem (corresponding to the point at which xylem have lost 50 % and 20 %, respectively, of their maximum conductivity).

Figure 2 Impact of hydraulic strategies on ecosystem water and energy fluxes: (a) monthly mean gross primary productivity and (b) monthly mean latent heat flux of the endmember cases.

The emergent behavior of FATES or any model with a dynamic ecosystem structure can make analysis of physiological rate variation difficult, as the stand structure will respond and thus also vary when parameters are changed. Here, we want to first understand the direct trait control in the absence of structural differences. To overcome complications of the dynamic structure, we use a reduced-complexity configuration for running the model which we refer to as the “static stand structure” mode. To investigate dynamic competitive effects when growth and mortality are the next step, in this mode, the stand structure is initialized from observed forest census data and subsequently is fixed; i.e., the model does not permit plant growth or death to change the vegetation structure. This allows the direct assessment of hydraulic and physiological parameter variation in the model without the consequent feedback loops associated with a varying ecosystem structure. The stand structure is initialized with census data from the CZ2 site (Table S1 in the Supplement) and thus includes multiple cohorts of differently sized trees. Because this type of model configuration ignores prognostic plant mortality, in the interest of being able to compare across simulations where mortality rates might otherwise be very high, we use the loss of xylem conductivity as a measure of the mortality risk of conifer trees at CZ2, which has widely been used as an indicator of the drought mortality of forest (e.g., Hammond et al., 2019).

To force the model with an atmospheric upper boundary, we use the Multivariate Adaptive Constructed Analogs (MACA) climate data (Abatzoglou and Brown, 2012) from 2008 to 2015 of a 4 km × 4 km grid that covers the study area. The daily average MACA data are disaggregated to 3-hourly climate data (see Appendix S2 in Buotte et al., 2018, for details). We set the initial soil water content to 75 % of the saturated water content, close to field capacity. We believe that this is a realistic value because the model is initialized in January, when the study area has high precipitation and trees are all in a dormant state, and in a year when there is no drought. We also tried to initialize the soil with a higher water content (e.g., saturation) but did not find any differences, as the extra water drained quickly in the winter when transpiration is low.

3 Results

3.1 Sensitivity of GPP and ET to parameter perturbations

The parameter sensitivity analysis revealed that, in a monthly mean flux comparison, the simulations with deep roots provided a better match to the overall observed pattern of the GPP and ET (Fig. 2). In general, the simulated transpiration contributed 90 % of the ET. The deep-rooted cases more accurately captured the seasonality (e.g., the peak time) and the declining trend of observed GPP from 2011 to 2015. The deep-rooted cases also matched the observed ET fairly well. The simulated GPP of the shallow-rooted cases was higher than the observed values during the wet seasons (December to March) but much lower than those during the dry season of the pre-drought period. Overall, the simulated ET of the shallow-rooted cases was lower than the observed values. To quantify this assessment, we computed the RMSE from the hourly mean GPP and ET of each month and each year of all 40 cases (Fig. S2). We chose the RMSE as it is a common and compact metric for assessing model performance, though we note that other metrics could in principle be used, each of which has different advantages and disadvantages (e.g., Collier et al., 2018). The RMSE of the GPP and ET decreased with both effective rooting depth and P_50 gs for both xylem strategies (Fig. 3). The P_50 gs had less of an impact on the RMSE of the GPP for the case with safe xylem than on that of the GPP for the case with efficient xylem. In terms of the GPP, the effective rooting depth of 6.5 m provided the best fit, as indicated by the darkest color (GPP RMSE = 1.12 gC m⁻² s⁻¹, ET RMSE = 250 W m⁻²). This underscores the importance of deep roots for maintaining transpiration and photosynthesis during the dry season as well as the role of deep roots in increasing the relative decline in these fluxes during droughts.

Among the parameters that we varied in the sensitivity analysis, the vertical root distribution had the largest impact on the GPP and ET at CZ2. Figure 2a–b show the monthly mean GPP and ET of the endmembers of the sensitivity analysis (see Fig. S1 in the Supplement for the complete set of outcomes). We acknowledge that the variation in rooting depth across the ensemble was large; however, we also highlight that so too was the uncertainty in plant rooting depth and also that the uncertainty in rooting depth was less well quantified than other plant traits (e.g., P₅₀), such that this wide variation reflects a real and deep uncertainty in plant rooting profiles. Deep roots resulted in substantially higher GPP and transpiration during normal years (2011 and 2012). During long-term droughts, when deep soil moisture was depleted, the relative advantage of deep roots over shallow roots was reduced. Shallow roots resulted in substantially lower GPP and transpiration during the dry season (August to October), with the seasonal maximum occurring earlier – in May – as opposed to July with the deep-rooted cases. The shallow-rooted cases also had much lower GPP and ET during the dry seasons of the pre-drought period. During the late stages of drought (2014 and 2015), the GPP and ET of the different cases became more similar between the shallow- and deep-rooted cases.

Figure 3 Root mean square error of GPP (a–b) and latent heat flux (c–d) with respect to variation in input parameters.

The second most important set of parameters for controlling carbon and water fluxes contains those that govern the stomatal strategy. The simulations with a riskier strategy (P_50 gs=P_50 x) provided higher GPP and ET than those with a safer strategy (P_50 gs=P_20 x) during the pre-drought periods and the early stage of the drought (2011 to 2013), but they provided slightly lower GPP and ET in the late stage (2014 and 2015) for the deep-rooted cases. However, risky stomata provided slightly higher GPP and ET at all times for the shallow-rooted cases. The xylem strategy had the smallest effect on the GPP and ET of the parameters that we varied (e.g., the RMSEs of ET were both approximately 260 W m⁻² for safe and efficient xylem, respectively, with P_50 gs=P_20 x and an effective rooting depth of 8 m). In the deep-rooted cases, the safe xylem and efficient xylem strategies resulted in almost identical GPP and ET, which can be seen in the widest overlap between the dashed and solid lines in Fig. 1. In the shallow-rooted cases with safe stomata, safe xylem generated slightly higher GPP and ET than efficient xylem. In addition, the strength of the effects of the stomatal and xylem strategies also depends on the rooting depth. The deeper the effective rooting depth, the less significant the impacts of the stomatal strategy (Fig. S1).

Figure 4 Seasonal and interannual variation of plant physiological characteristics: (a) monthly mean stem fraction of conductance $K/K_{\text{max}}$ (b), monthly mean leaf water potential, and (c) monthly mean overall absorbing root water potential of the 55 cm DBH cohort throughout the 2011–2015 period.

3.2 Sensitivity of plant water status to parameter perturbations

We examined the impact of vertical root distributions and stomatal and xylem strategies on the seasonal variation of the following three plant physiological variables, which served as indices of plant water stress (Fig. 4): the fraction loss of the xylem conductivity of the stem (SFL), LWP, and an overall absorbing root water potential (AWP). In the model, absorbing roots in different soil layers had different water potentials associated with the soil water potential of that layer. We calculated a cohort-level effective AWP as the root-fraction-weighted average of water potential in absorbing roots across all the soil layers. Thus, the AWP represents the overall rhizosphere soil moisture condition that is sensed by the tree. These physiological variables were tracked for each cohort. In any given case, the differences in these variables among differently sized cohorts were negligible (Fig. S3). Therefore, we present the outcome of all cohorts with a diameter at breast height (DBH) between 50 and 60 cm, the size class that was most abundant at CZ2.

Stomatal and rooting strategies together controlled the loss of xylem conductivity during the dry season of the pre-drought period and the whole period of the long-term drought (Fig. 4a). In all the cases, the xylem conductivity reached a maximum during the wet season (December to January), started to decline during the growing season (April to June), and then reached its minimum in the dry season. With the same stomatal strategy, deep roots led to a less extreme loss of xylem conductivity than shallow roots. A deep rooting strategy was also able to maintain xylem conductivity with very little seasonal loss during the pre-drought period; however, as deep soil moisture was depleted, this effect was reduced. With a shallow rooting profile, the xylem conductivity started to decline earlier, and the minimum was much lower than that of a deep rooting profile. For example, with risky stomata, the minimum fraction of xylem conductivity in the deep-rooted cases in 2012 was 0.4, but it was lower than 0.2 with shallow roots. Unlike the deep-rooted cases, the seasonal variation of the loss of xylem conductivity did not differ too much during the pre-drought and drought periods in the shallow-rooted cases. Furthermore, during the very late stage of the drought, the deep-rooted cases had a lower fraction of xylem conductivity than the shallow-rooted cases (e.g., in January 2015).

In general, risky stomata allowed a greater loss of xylem conductivity ($K/K_{\text{max}}$) than safe stomata, but the extent depended on the vertical root distribution. The effect of the stomatal strategy was more obvious in the shallow-rooted cases. Risky stomata combined with shallow roots resulted in an increase in the duration of the 50 % loss of xylem conductivity as well as a maximum loss of xylem conductivity during the dry season. With a deep rooting strategy, the difference in the percentage loss of xylem conductivity between the safe and risky stomatal cases increased with the progression of the drought; however, with a shallow rooting strategy, this difference remained approximately the same over time. In addition, in 2011, a very wet year with deep roots, a safe xylem strategy was able to maintain the maximum xylem conductivity even during the dry season (Fig. 4a). The impact of the xylem strategy on the percentage loss of xylem conductivity was relatively weak. In both the deep- and shallow-rooted cases, trees with safe xylem lost less xylem conductivity during the wet season but more conductivity during the dry season.

Figure 5 Impact of different combinations of rooting depth, xylem, and stomatal traits on soil moisture; the left column shows deep-rooted cases with (a) efficient xylem and risky stoma, (c) efficient xylem and safe stoma, (e) safe xylem and risky stoma, and (g) safe xylem and safe stoma. The right column shows shallow-rooted cases with (b) efficient xylem and risky stoma, (d) efficient xylem and safe stoma, (f) safe xylem and risky stoma, and (h) safe xylem and safe stoma. The red circles highlight the effects of stomatal traits on deep water storage during the wet season of the pre-drought period.

The safe stomata and safe xylem cases for both deep- and shallow-rooted trees experienced greater declines in stem conductivity compared with the safe stomata and efficient xylem for the corresponding rooting depths (Fig. 4a). This is because, with safe stomata, trees operate at the right end of the vulnerability curve displayed in Fig. 1b, where the hydraulic conductivity of efficient xylem is much higher than that of safe xylem. Thus, when the same amount of water is transpired, efficient xylem will lose less water potential than safe xylem. This keeps the xylem water potential of a plant with efficient xylem higher than that of one with safe xylem, and consequently it also keeps the xylem conductivity, K, higher. This is also because we set P_50 gs based on P_xylem; thus, the P_50 gs of safe stomata for plants with efficient xylem was higher (less negative) than that of plants with safe xylem, resulting in lower transpiration rates, which in turn reduced the loss of xylem water potential. As a result, plants with both safe stomata and efficient xylem not only transpired less water, but also lost less water potential per volume of transpired water. Together, these two mechanisms contributed to keeping the xylem conductivity of the efficient xylem cases higher.

In addition, stomatal, rooting, and xylem strategies had similar impacts on the seasonal variation of both leaf and fine root water potentials (Fig. 4c and d). Leaf and fine root water potentials peaked during the winter, started to decline in early spring, and reached their lowest point in the dry season. Deep roots, safe stomata, and safe xylem traits all contributed to the maintenance of higher leaf and fine root water potentials during the growing and dry seasons. With deep roots, there was less of a difference in leaf and fine root water potentials between stomatal and xylem strategies in the very wet year of 2011. Plants that combine safe stomata and/or safe xylem with deep roots were able to keep the leaf and fine root water potentials relatively high (less than −5 Mpa) during the dry season of the drought period. However, while plants that combine risky stomata or efficient xylem with deep roots could keep the dry season leaf water potential above −5 Mpa during the pre-drought period, their traits led to the dry season leaf water potential dropping below −8 Mpa or even −10 Mpa during the drought period. In both the deep- and shallow-rooted cases, safe xylem led to much lower leaf and fine root water potentials during the dry season. The seasonal and interannual variation of fine root water potentials was almost identical to the leaf water potential, except that the water potential of fine roots was slightly higher (∼0.5 Mpa) than the leaf water potential.

3.3 Sensitivity of subsurface hydrology to parameter perturbations

In the simulation outcomes, the vertical root distributions again had the largest impact on hydrological processes, subsurface water content, and how they changed over the drought. With deep roots, there was less drainage loss from surface and subsurface runoff compared with shallow roots, especially during the growing season (Fig. 5a, c, e, and g). The subsurface water content exhibited different vertical and temporal patterns between the cases with different vertical root distributions. In the deep-rooted cases, during the pre-drought period, the water content in the deepest layers fluctuated between the wet and dry seasons. During the first year of the drought, the water content of the deepest layers (6–8 m) slightly increased during the wet season, but as the drought progressed, the soil water content became consistently depleted in the middle and deep layers (5–8 m), and only the shallow layer's (<0.16 m) water content increased during the wet season. In the shallow-rooted cases (Fig. 5b, d, f, and h), soil moisture in the surface layers (top 2 m) exhibited seasonal variation, but this became weaker with depth. Moreover, the soil moisture at 6–8 m depth stayed consistently high throughout the year during the pre-drought period and remained slightly low throughout the entire drought period, while the water content of the middle and upper layers of the shallow-rooted case exhibited a similar pattern of seasonal variation before and during the drought.

The stomatal strategy, as quantified by P_50 gs, had a weak impact on hydrological processes and soil moisture. In both the deep- and shallow-rooted cases, riskier stomata led to a slightly lower total subsurface water content (Fig. 6a). The effect of P_50 gs was less significant during the pre-drought period for both the deep- and shallow-rooted cases and became more significant as the drought progressed. The effect of P_50 gs on the total subsurface water content was less significant in the shallow-rooted cases. Figure 5c presents the effect of P_50 gs on the water content of the shallow and deep soil layers. In both the shallow- and deep-rooted cases, increasing P_50 gs had a negligible impact on the water content of the shallow layers during both the pre-drought and drought periods (Fig. 5c, left). For the deeper layers, in the shallow-rooted case, P_50 gs had no impact on the water content at any time; in the deep-rooted cases, a risky P_50 gs resulted in a lower dry season water content of the deep layers (7–8 m) during the pre-drought period (as indicated by the red circles in Fig. 5a and c), but it decreased the water content of those layers year-round during the drought period (Fig. 5a and e). In the deep-rooted cases, safe stomata with efficient xylem led to a slightly higher water content in the deep layers (5–8 m) during the pre-drought period and in the shallow layers (0–3 m) during the drought period (Fig. 6a). Risky stomata with safe xylem in the deep-rooted cases were most effective at accessing soil water. Although the soil water contents were generally high in the shallow-rooted cases, stomatal and xylem strategies exhibited a similar impact on soil water storage to those in the deep-rooted cases (Fig. S4).

Figure 6 Impact on hydrologic processes: (a) monthly mean total runoff and (b) monthly mean total soil water content of the entire soil column.

Furthermore, simulations with deep roots resulted in almost no loss of soil water to drainage during the dry season in normal years or during the whole drought period; by contrast, in simulations with shallow roots, the drainage loss was high during the pre-drought period and decreased through the drought period, but there was still some runoff even at the end of the period (Fig. 6a). The observed total annual runoff from the 2008–2011 pre-drought period was approximately 250 mm yr⁻¹, but it was zero during the 2012–2015 drought period (from Fig. 4; Bales et al., 2018b). This observed difference in runoff between the pre-drought (∼290 mm yr⁻¹, 2011–2012) and drought periods (∼0 mm yr⁻¹) in the deep-rooted case was consistent with the predicted pattern. During the pre-drought period, the wet-season total subsurface water contents from December to February were similar between the cases with deep and shallow roots; however, during the dry season (from June to September), the total subsurface water content with shallow roots was substantially higher than the case with deep roots (Fig. 6b).

4 Discussion

4.1 Vertical root distribution as the first-order control

The outcomes of our simulations indicated that the vertical root distribution exerted first-order control on both ecosystem-level fluxes and plant physiology at CZ2. This dominance of the rooting strategy over other hydraulic traits is related to the nature of the rainfall pattern in the Mediterranean-type climate of that region. The CZ2 site receives effectively all of its rain during the winter. This water is stored in the soil column and is slowly released through the growing season. The root zone soil moisture exhibits strong seasonal variation, which constrains plant water use and gas exchange as a function of the gradual drying of the soil column (Bales et al., 2018b). In the model, the stomatal behavior was controlled by the leaf water potential, which itself was strongly affected by the root zone soil moisture. In our simulations, the daytime average leaf water potential of a 55 cm DBH cohort was correlated well with the fine root water potential and was always approximately 0.5 Mpa lower (Fig. S5). This offset is consistent with the relationship between the midday leaf water potential and the pre-dawn leaf water potential found by Martínez-Vilalta et al. (2014) at the global scale.

With deep roots, trees used more subsurface storage capacity at the CZ2 site. In wet years such as 2011, the root zone water potential of deep-rooted trees was kept relatively high (Fig. 4b), and the trees operated at the upper end of their vulnerability curve throughout the year, with a typical loss of conductivity of 10 % (Fig. 7). Therefore, we did not observe much of an effect of stomatal strategy on the GPP and transpiration in wet years. At the upper end of the vulnerability curve, stomata were fully open regardless of stomatal strategy (either safe or risky). When the drought began in late 2012, annual rainfall fell below the total root zone storage, and therefore the deep storage remained depleted throughout the year. During the drought, the deep-rooted trees were able to operate at the high end of the vulnerability curve in the wet season, when the rainfall recharged the surface layer. As the surface layers dried, water potential then gradually fell to the lower end of the vulnerability curve; consequently, photosynthesis and transpiration started to drop as the dry season progressed. With risky stomata, trees can drive soil moisture to lower values. This is why we observed the difference in the effect on the GPP and transpiration between different stomatal strategies during the dry season as the drought progressed.

With shallow roots, trees can only use surface soil moisture storage. As a result, the surface water storage was quickly used up after the wet season, and the root zone water potential dropped to near the low end of the vulnerability curve during the dry season. Thus, shallow-rooted trees operated along the full extent of the vulnerability curve year-round during both the pre-drought and drought periods. Therefore, for the deep-rooted cases, we did not observe a strong effect of stomatal strategy on the GPP and transpiration during the wet season; however, unlike the deep-rooted cases, the effect of stomatal strategy on the GPP and transpiration during the dry season could be observed throughout the whole simulation period.

Figure 7 Simulated leaf water potential and fraction loss of conductivity ($K/K_{\text{max}}$) of all the cases, which follow the two vulnerability curves.

Rooting strategies greatly control the spatial pattern of vertical soil water content (Fig. 5). With deep roots, the vertical soil moisture variation was more homogeneous due to the extensive root distribution. With shallow roots, the soil became extremely dry at the surface (<1 m) and extremely wet in deep layers (>5 m) as a result of the aggregated root distribution in the upper layers. This finding is similar to that of Agee et al. (2021), who found that extensive lateral root spreading results in homogeneous soil moisture distribution. The homogeneous soil moisture pattern may contribute to a more energy-efficient system that reduces plant water stress (Agee et al., 2021) because it minimizes the loss of energy dissipation through water transport (Hildebrandt et al., 2016). Both our study and that of Agee et al. (2021) emphasize the importance of the means by which the root distributions determine how the subsurface storage is used.

Given the shape of the vulnerability curves, plants stopped transpiring in all simulations when their leaf water potential reached approximately −10 MPa with efficient xylem or −15 MPa with safe xylem, depending on their stomatal strategy (Fig. 7). Because we held the stand structure and leaf area constant to allow comparisons between cases, the simulated leaf water potential of the shallow-rooted, risky stomata combination could reach as low as −15 Mpa (Fig. 4b) during dry seasons, even during the pre-drought period, which is well below the lowest possible leaf water potential observed (−10 Mpa; Vesala et al., 2017). Leaves would likely be wilted before the water potential drops below −10 MPa, and the tree would already have shed the leaves due to canopy desiccation. However, we specifically did not permit that in these simulations to keep the different cases comparable. Although it might be unrealistic, the leaf water potential can serve as an indicator of the degree of canopy desiccation. With no or very few leaves, trees would rely on stored carbon to support respiratory demand until the wet season arrives to regrow leaves. Depending on the duration of the dry season, trees may exhaust their stored carbon and die from carbon starvation. While risky stomata can generate a higher GPP (Fig. 1a), they also result in a longer duration of more negative leaf water potential (Fig. 4b). This suggests that shallow-rooted pines at CZ2 with risky stomata would benefit from allocating more net primary productivity to their storage pools rather than growth in order to reduce carbon-starvation mortality. Therefore, even though the model generated unrealistically low leaf water potentials, the extent and duration of these simulated potentials allowed us to gain some insights into the interaction of plants' hydraulic strategies and the life history strategies of conifer trees under a Mediterranean-type climate. Furthermore, the unrealistic leaf water potential from the shallow root simulations indicated that the trees at that site must have very deep roots to exist there, which is in agreement with the conclusions of Goulden and Bales (2019).

In this simulation, the impacts of xylem traits on the GPP and ET were weak and subtle. This was the result of the relative position of the two vulnerability curves, particularly the intersection of the two vulnerability curves in absolute conductivity. When the absolute conductivity was plotted as a function of pressure (Fig. 1b and solid lines in Fig. S6), we observed, on the left side of the intersection, that the safe xylem was not only safe but also efficient. Thus, a safety–efficiency tradeoff of xylem only occurred on the right side of the intersection point. Therefore, in the shallow-rooted cases, when the root zone water content and hence the plant water status are low, safe xylem can generate slightly higher GPP and ET than unsafe xylem. Furthermore, the two pressure–conductivity curves diverged mainly at the wet end (corresponding to the wet season). This is likely to be because the xylem structures of conifers are very similar and the range of variation of xylem traits in the sensitivity analysis was limited to the dominant species at the site. Therefore, the difference in the xylem traits of conifers does not cause significant impacts on the ecosystem-level fluxes under the Mediterranean-type climate of CZ2, which are constrained by energy during the wet season (Goulden et al., 2012). In addition, the maximum rates of GPP and ET are co-constrained by the stand density, total leaf area, maximum stomatal conductance, and vapor pressure deficit (VPD). In this study, we used the static stand structure mode of FATES-Hydro, where the stand density and the total leaf biomass (i.e., total leaf area) of the trees were held constant. This further limited the effect of xylem traits on the GPP and ET.

Figure 8 Simulated average annual GPP and fraction of conductance of a 55 cm DBH cohort with (a) deep roots (effective rooting depth = 8 m) and (b) shallow roots (effective rooting depth = 1 m).

4.2 Balancing productivity and mortality risk

The hydraulic traits that contribute to high carbon fixation rates often make trees more susceptible to drought. Stomatal strategy (P_50 gs) can have both positive and negative impacts on trees, creating a tradeoff in the balance between productivity and physiological stress. The risky stomata (P_50 gs P_50x) can generate a higher GPP but also result in a greater loss of xylem conductivity and lower leaf water potential. The tradeoff varies depending on the plant's root strategy (i.e., having a deep-rooted distribution vs. a shallow-rooted distribution) and moisture state.

To better understand the tradeoff between productivity and mortality risk, we plotted the simulated annual average GPP for each year against the fraction of conductivity ($K/K_{\text{max}}$) of a 55 cm DBH cohort for two scenarios – deep roots (Fig. 8a) and shallow roots (Fig. 8b) – with different combinations of xylem and stomatal strategies. In both scenarios, for each pair of xylem and stomatal strategies, the GPP per tree increased almost linearly with $K/K_{\text{max}}$. However, as the safety of the stomata increased, the GPP declined more quickly with the loss of conductivity. This response was stronger in the deep-rooted scenarios. Having efficient xylem only slightly increased the steepness of the lines. The stomatal strategies thus represented points along a gradient of the tradeoff between growth and mortality risk: the safer the stomata, the more the GPP was traded to reduce the mortality risk.

Along this tradeoff space, the point at which trees can maximize their net carbon gains likely depends on the xylem traits. Studies have demonstrated that trees may temporarily lose xylem conductivity during mild droughts, which can recover once soil water becomes available. However, under extreme drought, their xylem could collapse and permanently damage the xylem conduits. In this case, trees rely on new sapwood growth to support transpiration (Brodribb et al., 2010; Anderegg et al., 2013). At one extreme, if the stomatal behavior is too safe, it will produce a low GPP and the tree will be out-competed for light due to faster-growing neighbors; however, at the other extreme, if the stomata behave very aggressively (riskily), this will produce a high GPP but also empty the subsurface storage quickly, consequently leading to a prolonged dry period of soil moisture. This would lead to substantial xylem damage (and/or root death), and then the carbon required to grow new sapwood (or roots) could exceed the benefits from the additional GPP. Thus, the optimal location along the gradient would probably be located slightly below the $K/K_{\text{max}}$ associated with that critical xylem water potential. Currently, the xylem refilling and the associated carbon cost are not incorporated into FATES-Hydro. These two processes should be implemented in the model to obtain an enhanced understanding of the water–carbon balance, which remains for future work.

In the deep-rooted scenario, the values of the pre-drought period and early drought stage were clustered in the upper-right corner, above a $K/K_{\text{max}}$ of 0.6 (Fig. 8a). In this region, the stress from the loss of xylem conductivity likely would not be high enough to cause severe consequences if we were to use a 50 % loss of xylem conductivity as the threshold for mortality and/or permanent xylem damage. The deep-rooted tree can thus benefit by trading less GPP for maintaining xylem conductivity with a risky or more productive stomatal strategy during normal years. However, during the late stages of the drought (2014 and 2015), the conductivity values became much lower. If this mega-drought had stopped earlier (e.g., if it had been a mild drought that only lasted for 2 years), then the additional GPP obtained from risky stomata might have outweighed the carbon cost of repairing xylem damage. This suggests that if the 2012–2015 drought was not common in California, then natural selection might favor the risky or more productive stomatal strategy for deep-rooted trees. However, the same strategy would also expose trees to a high mortality risk during severe droughts.

In the shallow-rooted case (Fig. 8b), the values were all clustered lower and to the left compared with the deep-rooted scenario, irrespective of the drought status. Thus, for shallow roots, risky or more productive stomatal behavior resulted in a similarly high mortality risk during both the pre-drought and drought periods. Accordingly, under the long-term climate conditions found at CZ2, regardless of whether severe droughts were frequent or not, the only shallow-rooted trees that could persist would have to follow the safe and less productive stomatal strategy. Therefore, safe and less productive stomata also protect shallow-rooted trees from mortality risk during drought.

The model outcome indicated that, under drier root zone soil conditions, if pines were to follow a shallow rooting strategy, they would benefit from a safer stomatal strategy with more conservative water use. By contrast, if they were to follow a deep rooting strategy, they would benefit from riskier stomata. This is consistent with the Anderegg et al. (2016) finding regarding the relative stomatal conductance (gs) across elevation. They found that, at a low-elevation (lower-precipitation) site, P. ponderosa had a lower relative stomatal conductance and less loss of xylem conductivity, equivalent to safer stomata in our study, while at a mid-elevation (higher-precipitation) site, it had a higher relative stomatal conductance and a greater loss of xylem conductivity, equivalent to risky stomata in our study. The simulation results are consistent with the idea that the CZ2 region is dominated by deep-rooted trees, which is supported by previous studies. In situ measurements of regolith structure (particularly porosity) indicate that, at CZ2, there is a layer of thick, semi-weathered bedrock that allows trees to grow deep roots (Holbrook et al., 2014). Growing deep roots to access rock moisture to support plant water use was also observed in the Eel River CZO catchment (Rempe et al., 2018), another Mediterranean-type ecosystem along the western coast. Observed net CO₂ exchange and ET during the pre-drought period suggested that, during a wet year, deep moisture supported summer transpiration and productivity when the upper-layer moisture was low (Goulden et al., 2015). Because the deep rooting strategy is sufficient in most cases for avoiding the main effects of dry seasons and short droughts and because, conditional on having deep roots, the risky stomatal strategy confers a productivity advantage at a little increased risk of vulnerability, we would expect plants with these traits to dominate. However, under extreme cases such as the 2012–2015 drought, which ranked as one of the most severe in California in the last 1200 years (Lu et al., 2019), we would expect plants with this deep-rooted, risky stomatal strategy to be highly vulnerable to drought, which is consistent with the approximately 90 % mortality of the pine observed at CZ2 during the drought (Fettig et al., 2019). The water balance of the catchment based on long-term observation from precipitation, streamflow, and ET (Bales et al., 2018a; Goulden and Bales, 2019) also supports this being the slow depletion of deep moisture that caused tree mortality in the late stage of the prolonged 2012–2015 drought.

The findings of our study indicate that future drought mortality will likely occur in ecosystems that are limited by water and other factors. In such ecosystems, trees can benefit from having more efficient but less safe hydraulic traits, as they allow them to be more competitive for water and to bring in a higher GPP. The extra carbon gain can be used to develop measures for dealing with other constraining factors, such as increasing stored carbon to lower the risk of carbon starvation, building thicker bark to resist fire, and growing more roots, which further enhances their capacity to compete for water.

5 Conclusions

Our analysis indicates that root distribution can affect the most competitive stomatal traits. In a Mediterranean-type climate where the supply of energy and water is desynchronized and the accessible subsurface water storage capacity is close to annual precipitation, deep roots combined with risky stomata represent a beneficial strategy for high productivity in normal years with a low mortality risk; however, this strategy exposes trees to a high mortality risk during multi-year droughts. While such a strategy enables trees to fully use subsurface storage and precipitation for productivity over regular years, the lack of deep water storage recharge during droughts exposes trees to high drought stress and makes this strategy unfavorable under severe and prolonged drought. By contrast, shallow roots combined with safe stomata represent a strategy for drought resistance, albeit at the cost of considerably reduced productivity, as such a combination only allows trees to use shallow subsurface storage while leaving deep moisture untouched. Thus, less precipitation is used for productivity. However, this strategy leaves trees less susceptible to drought-induced mortality should the deep reservoir be depleted. By contrast, shallow roots with risky stomata lead to high mortality even during non-drought years, making this an uncompetitive combination at the site. These results suggest that stomatal strategy is controlled by root zone soil moisture and regulated by root distribution in that region. Thus, our study underscores the importance of considering plant rooting and hydraulic strategies within the larger context of plant ecological strategies.