Optimization of leaf morphology in relation to leaf water status: A theory

Abstract The leaf economic traits such as leaf area, maximum carbon assimilation rate, and venation are all correlated and related to water availability. Furthermore, leaves are often broad and large in humid areas and narrower in arid/semiarid and hot and cold areas. We use optimization theory to explain these patterns. We have created a constrained optimization leaf model linking leaf shape to vein structure that is integrated into coupled transpiration and carbon assimilation processes. The model maximizes net leaf carbon gain (NPPleaf) over the loss of xylem water potential. Modeled relations between leaf traits are consistent with empirically observed patterns. As the results of the leaf shape–venation relation, our model further predicts that a broadleaf has overall higher NPPleaf compared to a narrowleaf. In addition, a broadleaf has a lower stomatal resistance compared to a narrowleaf under the same level of constraint. With the same leaf area, a broadleaf will have, on average, larger conduits and lower total leaf xylem resistance and thus be more efficient in water transportation but less resistant to cavitation. By linking venation structure to leaf shape and using water potential as the constraint, our model provides a physical explanation for the general pattern of the covariance of leaf traits through the safety–efficiency trade‐off of leaf hydraulic design.

tion between several leaf traits that are believed to comprise the adaptive strategy that affects leaf economy, known as the leaf economic spectrum (Wright et al., 2004). Furthermore, leaf surface area has been found to be positively associated with soil moisture levels (Wright et al., 2017;Xu, Guo, Xu, Wei, & Wang, 2009), and net carbon assimilation rate and leaf area are found to decrease toward the top of a tree, with decreasing petiole xylem water potential (Ambrose, | 1511DING et al. Sillett, & Dawson, 2009Ishii, Jennings, Sillett, & Koch, 2008;Kenzo et al., 2015;Koch, Stillet, Jennings, & Davis, 2004;Ryan, Phillips, & Bond, 2006). This phenomenon suggests that water potential can also affect leaf shape. Most previous studies on leaf traits have focused on thermal regulation; less attention has been paid to the leaf water status. The purpose of this paper is to propose a mechanism for the leaf economic spectrum in relation to water status.
The leaf water status directly affects the leaf carbon budget through a series of biochemical processes. When the water potential of the leaf vein conduit reaches a critical low value, cavitation occurs. Cavitation in minor vein conduits blocks their ability to fill with water and causes the loss of cell turgor pressure, which in turn negatively affects the biochemical processes of photosynthesis by reducing the maximum carboxylation rate, damaging ATP, and increasing mesophyll resistance for CO 2 diffusion (Grassi & Magnani, 2005;Lawlor & Tezara, 2009;Niinemets, Cescatti, Rodeghiero, & Tosens, 2005). Leaf traits considered to control the leaf carbon budget related to water status include stomatal conductance, leaf area, leaf shape, and venation structure (Ball, Woodrow, & Berry, 1987;Berninger, Mäkelä, & Hari, 1996;Blonder et al., 2011;Cowan, 1982Cowan, , 1986Jarvis, 1976;Noblin et al., 2008;Wright et al., 2017;Xu et al., 2009).
Stomatal conductance and leaf area are the two fundamental leaf traits that are tightly associated with the leaf carbon budget, the leaf water status, and water loss. It has been shown that to avoid the negative impacts of low water potential, the leaf will increase its stomatal resistance to reduce the transpiration rate (Ball et al., 1987;Berninger et al., 1996;Cowan, 1982Cowan, , 1986Jarvis, 1976). Stomatal resistance can be regulated by instantaneously changing stomata size (opening and closure) and varying the stoma density of the leaf over a longer time scale (Sack & Buckley, 2016;Woodward & Kelly, 1995;Xu & Zhou, 2008). Carbon and water fluxes are coupled through stoma. This adjustment comes at the cost of a reduced carbon fixing rate, since increasing stomatal resistance slows the carbon intake rate. Leaf surface area has been found to be positively associated with soil moisture levels (Wright et al., 2017;Xu et al., 2009). Furthermore, the net carbon assimilation rate and leaf area are found to decrease with tree height, because of the decreasing xylem water potential associated with tree height (Ambrose et al., 2009;Ishii et al., 2008;Kenzo et al., 2015;Koch et al., 2004;Ryan et al., 2006). These studies suggest that the adjustment of leaf area to water availability over time can influence the optimization of stomatal resistance in response to fluctuations in temperature and vapor pressure (Maseda and Fernández, 2006).
Leaf shape and leaf venation are the two other important traits relating to the leaf carbon budget and water use strategy. The former determines external geometry, and the latter determines the structure of the internal transportation system. Both can affect the rate at which net leaf carbon gain (NPP leaf ) scales with leaf size. Leaf shape, a genetically determined trait, is found to be associated with habitat and is a trait that defines plant water use strategies. Leaf width is found to increase with increasing precipitation and/or soil moisture (Givnish, 1987;McDonald, Fonseca, Overton, & Westoby, 2003). Plants growing in arid or semiarid environments have smaller and narrower leaves, while plants growing in humid or semihumid environments have broader and larger leaves (Wright et al., 2017).
The leaf venation network has been connected to net carbon assimilation rate, life span, leaf mass per area, and nitrogen content, the traits that affect the leaf carbon budget (Blonder et al., 2011).
Evidence suggests that leaf shape and venation structure are related. To provide structural support for a leaf, the major vein must extend from the petiole to the end of the leaf, and second-order major veins must extend from the major vein to the boundary of the leaf (Cochard, Nardini, & Coll, 2004;Roth-Nebelsick, Uhl, Mosbrugger, & Kerp, 2001). Therefore, assuming the venation has a dendritic structure following the metabolic scaling theory West, Brown, & Enquist, 1997), the width-to-length ratio of a leaf regulates the vein network parameter γ (the ratio of the length of n + 1th order vein to nth order vein). Thus, a narrowleaf will have a smaller width-to-length ratio, and the consecutive veins will consequently shorten more quickly than those of a broadleaf (e.g., with the same leaf length, the second-order major vein of a narrowleaf will be shorter than that of a broad leaf). In addition, the optimal average distance between the two adjacent smallest minor veins approximates the distance from the minor vein to the leaf surface, which is approximately half of the leaf thickness (Noblin et al., 2008). These venation properties control the rate of total leaf xylem resistance scaling with leaf area (West, Brown, & Enquist, 1999).
The studies reviewed above suggest that venation structure, leaf shape, and stomatal resistance are all related and affect the tradeoffs between leaf size and the characteristics of the conductive system in leaves and transpiration. The total carbon uptake and the construction and maintenance costs of a leaf are proportional to leaf size. Therefore, these leaf traits affect leaf carbon budget and leaf water status through transpiration. The objective of this study is to develop a leaf physiology model to explore how the leaf deals with this trade-off to maximize its net carbon gain.
Optimality is suggested as an appropriate method to model the adaptation strategy of plants to their environments (Schymanski, Sivapalan, Roderick, Beringer, & Hutley, 2008). Recently, optimality has been used to predict the instantaneous response of stomatal conductance to environmental conditions such as variation of soil moisture and elevated atmospheric CO 2 concentration (Ball, Cowan, and Farquhar, 1988;Cowan, 1982;Franks and Casson, 2014;Manzoni et al., 2013;Medrano et al., 2015;Sperry et al., 2016Sperry et al., , 2017.
These models optimize instantaneous stomatal conductance either by water loss or the loss of plant hydraulic conductance  without explicitly considering the effect of other leaf properties, such as size, shape, and venation structure.
Here, we use an optimality approach to explain the coadjustment of leaf traits as an adaptive strategy in a given environment.
We present a constrained optimization leaf model that incorporates leaf venation and shape with coupled transpiration and carbon budget processes to address the coadjustment of these two leaf traits as an adaptive strategy. As mentioned above, stomatal conductance can be regulated through two means: the instantaneous opening and closing of stoma and the long-term adjustment of stoma density and aperture size. Our model focuses on the latter. Here, our model is not an optimization model of instantaneous stomatal conductance, in the sense that it does not accurately resolve the transpiration rate under temporally detailed realistic conditions; instead, the model provides an explanation of the adaptation strategy in an environment over the life span of the leaf rather than predicting/explaining the instantaneous response of stoma. The model also provides a physical link between stoma behavior and other leaf traits. Furthermore, our model uses water potential as the limiting factor, instead of water, as used in previous models (see Section 2.1 for justification).

| CON S TR AINED OP TIMIZ ATI ON LE AF MODEL
In general, a constrained optimization model has two components (functions): an objective function, describing the profit (net carbon gain) of a system (leaf) to be maximized, and a constraint function, defining the total cost of a factor allowed to make profit. Here, the constraint is the maintenance of the loss of water potential from the petiole to the terminal minor vein above a critical value, which depends on the petiole water potential. This is somewhat similar to imposing a hydraulic risk avoidance that works at the leaf level (i.e., without a description of the decrease in water potential from soil to petiole).
Below, we first justify the appropriateness of using water potential as the constraint factor and then describe the two components of the model: the objective function and the constraint function.

| Constraint: water or water potential
A critical issue with the optimality approach is the appropriate structuring of the model Wolf, Anderegg, & Pacala, 2016). Available water in soil has long been considered the limiting factor of plant growth. For example, the stomatal optimization theory assumes that the adjustment of stomatal resistance leads to the maximization of the cumulative carbon assimilation over a fixed loss of water through transpiration (Cowan & Farquhar, 1977, Katul, Oren, Manzoni, Higgins, & Parlange, 2012. However, in some regions (e.g., humid areas and topographic convergent areas), water is not limited. Furthermore, in places where plants can access groundwater, water itself is not limited, because no matter how much water the plant uptakes from the capillary fringe, it can be recharged from groundwater in the saturated zone. As such, water cannot be used as limiting factor in the optimization model for these regions. However, a cap of root zone water potential exists, because roots do not grow in the saturated zone, whether because it is biologically imperative to avoid anoxic conditions or because roots will physically rot if they are in water too long. Also, it is the lowering of water potential of the xylem that causes the cavitation and negative impact on plants. Therefore, the amount of water potential a plant can lose through transpiration is limited, and this limitation is valid not only in arid regions but also in humid ones. By using water potential as the limiting factor, optimization can be applied in these regions.
In addition, average root zone water potential is not solely controlled by precipitation-the input of water, but also the soil depth and the depth of groundwater. Thus, using water as the limiting factor cannot consider all of these factors, but using water potential does allow us to further integrate these aspects in the future. Furthermore, several leaf economic traits that vary with soil moisture are found globally. This indicates that water can regulate plant behavior not only in arid regions but also in humid areas. As such, using water potential, instead of water, as the limiting factor in the optimization model not only allows the optimization model to be valid for all regions but also holds the potential to integrate the processes of hillslope hydrology in the future. So, for all of these reasons, we use water potential as an alternative constraint.

| Constrained optimization leaf model
The constrained optimization leaf model has two components: an objective function, as the net leaf carbon gain (NPP leaf ), and a constraint function, the total loss of xylem water potential of a vein, starting at the petiole and terminating at the minor vein. The model maximizes the objective function, the NPP leaf , subject to the constraint function, the total loss of xylem water potential of a vein (petiole xylem water potential) over a leaf's life span.

| The objective function-net leaf carbon gain
NPP leaf is described here as the balance between the leaf carbon assimilation rate (the first term on the right side of Equation 1) and the rate of total carbon loss from respiration and construction over its life span (the second term on the right side of Equation 1). The leaf carbon assimilation rate is the gross primary productivity rate, or the net carbon fixation rate per unit of leaf area, GPP (μmol m −2 s −1 ), multiplied by the leaf area of a single leaf, LA (m 2 ): where NPP leaf (mol/s) is the net carbon gain of a mature leaf. R (mol m −1 s −1 ) is the coefficient of the total carbon cost of the leaf (the second term on the right side of Equation 1), the sum of the life span standardized construction cost (e.g., the total mass of the leaf divided by its life span) and the respiration loss, both are (1) proportional to total mass of a leaf. The total mass is the product of the average leaf mass density and the volume of the leaf.
Assuming the average density of the leaf does not change with leaf size, for a given leaf geometry, the total carbon cost of the leaf scales with leaf area, by an exponent of 3/2. LA (m 2 ) is the leaf area; C a (0.0004 mol/mol) is the CO 2 concentration of air; a 1 (mol m −2 s −1 ) is the maximum carboxylation capacity; and a 2 (510 mol/mol) is a given as K c (1 + Co a /K o ), where K c and K o are the Michaelis constants for CO 2 fixation and oxygen inhibition, and Co a is the oxygen concentration in the air. r es (s/m) is the mass stomatal resistance of water vapor; cons3 (0.025 m 3 /mol) is a constant for stomatal resistance conversion from mass to mol; S is the ratio of leaf CO 2 concentration to the CO 2 concentration of the air, which is considered a constant with a value around 0.6 (Katul, Manzoni, Palmroth, & Oren, 2010).
The formula of GPP used here is adapted from the linear model of Katul et al. (2010) at equilibrium. The linear model assumes that the ratio of inner leaf CO 2 concentration, C i , to the CO 2 concentration of the air in the denominator of the biochemical model is a constant, while allowing C i to vary in the numerator (Katul et al., 2010).
We use the linear model so that we can obtain a relatively simple and understandable expression.

| The constraint function-the total loss of xylem water potential (Ψ loss )
The constraint function describes the total loss of xylem water potential from the petiole to the terminal minor vein under the limitation that the total loss of xylem water potential is the product of the average mass flow rate per unit area of cross section of xylem, Jx (m 3 s −1 m −2 ), and the total xylem resistance per unit area of the cross section of a single flow path, r Xtotal (s/m) (see Appendix S1 for derivation): where Ψ loss (Pa) is the total loss of xylem water potential from the petiole to terminal minor vein; e sat is the saturated vapor pressure; e a is the atmospheric vapor pressure; cons1 (kg/mol) is the molar mass of H 2 O; cons2 (s 2 /m 2 ) is the inverse of the specific gas constant of water vapor; (2a) k Z is the empirical scaling coefficient providing the ratio of minor vein distance to leaf thickness; Z (m) is the thickness of the leaf; r Xtotal is the total xylem resistance of a single flow path from petiole to terminal minor vein; D m is the diffusivity of water in mesophyll; w 0 (mol/m 3 ) is the saturated water concentration of mesophyll; l T (m) is the length of the terminal minor vein; r T (m) is the radius of the terminal minor vein xylem conduit; µ (kg m −1 s −1 ) is the viscosity of water; ρ w (km/m 3 ) is the density of water; k SP is the shape parameter of the leaf (proportional to W/L ratio); γ and β are the parameters of the venation network (explained in the next paragraph); and ΔΨ Lmax (Pa) is the maximum total xylem water potential that is allowed to be lost along a single flow path from petiole to terminal minor vein without causing any negative impact on leaf health. It is the difference between the petiole xylem water potential Ψ P (Pa) and the turgor loss xylem water potential Ψ C (Pa) at which the percentage of the vein that forms cavitation will cause detrimental effects on biochemical processes. Ψ C is the hydraulic trait of plants; thus, the value depends on the plant species. Ψ 50 of the vulnerability curve of the sigmoid model (Markesteijn, Poorter, Paz, Sack, & Bongers, 2011;Pammenter & Willigen, 1998) will be a good estimation of Ψ C , because this is where the sharpest drop of conductance occurs (Venturas, Sperry, & Hacke, 2017).
The total xylem resistance is derived from a simplified hydraulic network structure following the metabolic scaling theory . The vein network of our leaf model has a dendritic structure, as described by Price and Enquist (2007). The veins contain xylem conduits, and water flows through the xylem conduits contained by veins ( Figure 1). Note that the vein network and hydraulic network composed by xylem conduits are not the same, but they do share some characteristics, described as follows. The veins are numbered in increasing order from major vein to terminal minor vein (e.g., the major vein is numbered as Order 0, the second major vein is numbered as Order 1, and so on). The xylem conduits branch where veins branch; therefore, the length ratio of two consecutive orders of veins and the length ratio of the two consecutive orders of xylem conduits are the same, which we denote as γ. However, the branching ratio of the vein and conduit differ. In our model, the branching ratio (the ratio of the number of conduits of two consecutive orders veins) of xylem conduit is denoted by m, and the vein branching ratio is denoted by n. The ratio of the radius of the two consecutive xylem conduits is denoted as β, assuming the conservation of total xylem conduit cross-sectional area West et al., 1997), xylem conduit branching ratio, m, equals 1/β 2 .
Based on the previous studies mentioned above, we make two connections between vein structure and leaf shape. First, the length ratio of two consecutive veins, γ, is proportional to half of the widthto-length ratio of the leaf. Second, the distance between two adjacent terminal minor veins is half of the leaf thickness ( Figure 1).
If the water potential of the xylem conduits-which we call xylem water potential in the rest of the paper-of the terminal minor veins drops below a threshold (Ψ C ), cavitation will occur, and thus leaf cell turgor pressure will no longer be maintained. The threshold xylem water potential (Ψ C ) is a functional hydraulic trait determined by the vulnerability curve (Markesteijn et al., 2011;Pammenter & Willigen, 1998) and cell physiology. The difference between the actual petiole xylem water potential and the critical xylem water potential determines how much total water potential can be lost through the flow pathway.
From Equation 2, we can see that both increasing leaf area and decreasing stomatal resistance can increase NPP leaf , but they also increase the total loss of xylem water potential. From an evolutionary perspective, plants evolved in a given environment will tend to maximize NPP leaf in that type of environment. Therefore, the model maximizes NPP leaf , the objective function (Equation 1), over the loss of total xylem water potential (Equation 2a), under the condition that the total loss of xylem water potential (Ψ loss ) must be less than the difference between the petiole xylem water potential and the critical xylem water potential (Equation 2b).
Below, we explain in detail how leaf shape (W/L ratio) affects the network structure of the leaf hydraulic system and hence the difference in the relation between total xylem resistance and leaf area.
The total number of xylem conduits of the major vein (0 Order vein) is the following: where M T is the number of xylem conduits of a single terminal minor vein, and m is the xylem conduits branching ratio (e.g., M k+1 /M k ).
Following the metabolic scaling theory, the vein network is assumed to have a space-filling, fractal structure (West et al., 1997). Thus, the total number of xylem conduits of the major vein, M 0 , scales to the power of the number of the total vein order, T, given by the following (see Appendix S1 for derivation): Expanding Equation 4, T can be expressed as the following: where k SP = W/L, and LA is the total leaf area (m 2 ). Both the numerator and denominator of Equation (5) decrease with increasing W/L, but the numerator decreases at a faster rate than the denominator. Therefore, T increases with increasing W/L ratio. M T is determined by leaf size, minor vein length, distance between minor veins, and number of xylem conduits contained by a single minor vein, assuming that the minor vein is invariant (e.g., similar in length and number of xylem conduits contained) and the minor vein distance is determined by leaf thickness. The number of minor veins is the total length of minor veins divided by the average length of a single minor vein, l T . Following the same logic of West et al. (1999) (see Appendix S1 for detail), the total length of minor veins is given as the leaf area divided by the average distance between two adjacent minor veins.
Thus, with the same LA, M T is the same for leaves differing in W/L ratios if they have the same thickness. With a larger T, the number of xylem conduits of the major vein of a broadleaf is less than that of a narrowleaf (Equation 3). With the conservation of the total area of the cross section of xylem conduit, the radius of xylem conduit of the major vein, r 0 , of a leaf having larger W/L is larger than of a leaf with a smaller W/L. Therefore, if the leaf size is the same, and assuming the conduit size of terminal minor vein is the same, leaves with different shapes (W/L ratio) will have the same total cross section area of xylem conduits at petiole, but leaves with small W/L ratios will have fewer orders of veins and a larger number of smaller size xylem conduits of the major vein, compared to a leaf with a large W/L (Table 1).
Here, using the link between leaf shape and venation, we examine the effect of leaf width-to-length ratio and leaf thickness on the coadjustment of the optimal stomatal resistance and leaf area with petiole xylem water potential and corresponding NPP leaf (as NPP leaf per unit of leaf area). For this purpose, we allow leaf width-  Table 2).
Equation (6)  The solution expresses the maximum NPP leaf , optimal leaf area and stomatal resistance, and corresponding GPP, by a set of biophysical and environment variables, as described in Equation 1.

| Model interpretation-optimization
To help in understanding why the model predicts certain patterns, we first interpret the optimization model graphically (Figure 2). Both is the isoline of the production function that touches the isoline of the constraint function (red dot). A leaf can only produce on and below the isoline of that ΔΨ Lmax , because above that line, the xylem water potential of terminal minor will be lower than the turgor loss value, and the leaf will be malfunction and no longer be able to properly perform photosynthesis. Setting a lower constraint (e.g., moving from the red to blue line) increases both stomatal conductance and leaf area, as well as NPP leaf (blue dot) (Figure 2c). This explains the general trends. The parameters of the model (Equation 2) affect the shapes of the isolines of the two functions. With different values, the exact location of the optima varies. Next, we discuss the effect of leaf shape as parameters on the leaf optima.

| Impact of leaf shape-venation relation on leaf optimum
First, we examine the effect of leaf width-to-length ratio on the change of leaf optimum, GPP, and NPP leaf with increasing petiole xylem potential (Figures 3 and 4  increase with ΔΨ Lmax at a much slower rate than those of a broadleaf ( Figure 3). This means that a broadleaf is more sensitive to the change of ΔΨ Lmax and will have higher marginal gain of productivity than a narrowleaf. However, on the other hand, if ΔΨ Lmax decreases, a broadleaf will also suffer greater productivity loss than a narrowleaf.
The optimal leaf sizes of narrow and broad leaves will increase with increasing ΔΨ Lmax , but that of a narrowleaf will increase with ΔΨ Lmax more slowly than that of a broad leaf; also, the optimal leaf size of a leaf with a small W/L will, overall, be smaller than that of a leaf with a large W/L (Figure 4a). The optimal stomatal conductance of leaves with different W/L ratios will increase with increasing ΔΨ Lmax ; the optimal stomatal conductance of narrow leaves will be lower than that of broad leaves at any petiole xylem water potential (Figure 4b). Assuming Ѱ C , a plant hydraulic trait, does not change for a given plant, the change of total allowable xylem water potential loss (ΔѰ Lmax ) represents the change of petiole xylem water potential (Ѱ P ); the result of the optimization shows that the GPP, NPP leaf , and optimal leaf area will all increase with petiole xylem water potential, while stomatal resistance will decrease, each at different rates. The marginal increases of the GPP and NPP leaf of a narrowleaf will be smaller than those of a broad leaf.

| Effect of leaf thickness
The model predicts that with different leaf thicknesses, NPP leaf , net carbon gain per leaf area (NPP), GPP, and optimal leaf area will all increase with petiole xylem water potential, while optimal stomatal resistance will decrease (Figures 3 and 4, green lines).
However, overall, along the gradient of maximum allowable water potential loss (ΔѰ Lmax ), a thicker leaf will have lower productivities, which also change with ΔѰ Lmax at a slower rate (Figure 3, green lines). This will result in a more significant effect of leaf thickness on NPP leaf and leaf optimums when petiole xylem water potential is high. In other words, NPP leaf , GPP, optimal stoma resistance, and leaf area are more sensitive to changes of petiole xylem water potential if the leaf is thin. A thicker leaf will also have higher optimal stomatal resistance and lower optimal leaf area than those of a thin leaf, at any given ΔѰ Lmax (Figure 4, green lines). The effects of leaf thickness on the above leaf properties are not as strong as the effects of leaf W/L ratio.

| Implications of the model
Our leaf model predicts two patterns: (a) stomatal conductance (1/res, inverse of stomatal resistance), leaf area, and GPP increase with increasing petiole xylem water potential; (b) with different leaf shapes, the aforementioned variables change at different rates with petiole xylem water potential. This is due to different venation structures.
Our first prediction agrees with other previous studies on the changing of leaf traits with tree height. These studies, conducted in different places, have found that net carbon assimilation rate and leaf area decrease with tree height due to decreased xylem water potential along the stem (Ambrose et al., 2009;Ishii et al., 2008;Kenzo et al., 2015;Koch et al., 2004;Ryan et al., 2006). Ishii et al. (2008) have found that bulk leaf water potential decreases with the height of the leaf for Sequoia sempervirens. They have also found that leaf area decreases from ~9-1 cm 2 from the bottom to top of the tree (Figure 3d in Ishii et al., 2008) and that net carbon assimilation rate first increases with height as a result of increasing light intensity, reaching the maximum at 80 m, then decreases with height as the xylem water potential decreases with height along the tree. Stomatal conductance also decreases with height from 250 mm mol m −2 s −1 at the base to 800 mm mol m −2 s −1 at the top of the tree (Figure 5c in Niinemets, 2002). Koch et al. (2004) have found that xylem water potential decreases with height linearly.
Leaf size and shape also change with height: leaves become smaller and more expanded. Using carbon isotopes, Koch et al. (2004) have confirmed that the reduction in photosynthesis rates with height is caused by the reduced stomatal conductance. These studies support the argument that leaf morphology and stomatal conductance are equally limited by leaf water status (Koch et al., 2004;Niinemets, 2002;Thomas & Winner, 2002;Woodruff, Bond, & Meinzer, 2004;Woodruff, Meinzer, Lachenbrunch, & Johnson, 2009). Empirical study has found the general pattern that species having higher maximum stomatal conductance quickly close their stoma at low vapor pressure deficit (VPD; Oren et al., 1999). Recently, Henry et al.
(2019) have developed a stoma response model that explains this pattern through four mechanisms. The optimal stomatal resistance in our model is equivalent to the maximum stomatal conductance, or gsmax.
Our model can be connected to their stoma model through optimal  Our second prediction is based on the assumption of the relation of leaf shape and hydraulic systems. Next, we discuss in detail the connection between leaf shape and hydraulic conductance through the regulation of vein structure.

| Effect of leaf shape on leaf hydraulic conductance
Leaf shape affects NPP leaf , optimal stomatal resistance, and leaf area through the vein network structure parameters of the constraint function-γ and β within our leaf model-and leaf width-to-length ratio determines the leaf hydraulic network parameter γ (i.e., the length ratio of daughter-to-mother conduits [and veins]). This parameter determines the rate at which total xylem resistance r Xtotal scales with leaf area. Increasing W/L ratio causes at the same level Ψ loss isoline to move northwest ( Figure 5). This means that given the same amount of water potential that can be lost, a leaf with a larger W/L ratio can support higher carbon production factors, stomatal conductance, and leaf area than a leaf with a smaller W/L ratio ( Figure 5). This occurs because the rate at which the total xylem resistance increases with leaf area is much faster when γ is smaller (Equation 2, Figure 6). This causes the narrowleaf to have smaller optimal leaf area but higher optimal stomatal resistance than those of a broadleaf and thus a lower maximum F I G U R E 2 Graphic interpretation of the optimization model: (a) isolines of objective function (Equation 1); (b) isolines of constraint function, and (Equation 2); (c) coadjustment of stomatal resistance and leaf area to maximized net leaf carbon gain when the level of the constraint function increases from a low level (red dash line) to a higher level (blue dash line)

F I G U R E 3
Change of (a) net leaf carbon gain (NPP leaf ); (b) net leaf carbon gain rate per unit leaf area (NPP); and (c) net carbon assimilation rate per unit leaf area (GPP) with total allowable loss of xylem water potential NPP leaf . The stomatal resistance of a narrowleaf is consequently higher than that of a broadleaf with the same leaf area.
It should be noted that the minor vein density controls the relation between the covariance of stomatal conductance, leaf area, and total xylem water potential loss. Leaf thickness affects optimal stomatal conductance, leaf area, and NPP leaf if the minor vein distance is proportional to leaf thickness, as observed in angiosperms (Noblin et al., 2008). In our model, minor vein density is linked to leaf thickness through the term k Z ⋅ Z∕2, which represents the distance between the two adjacent minor veins (i.e., the inverse of minor vein density). For a given level constraint, increasing minor vein distance not only shifts the isoline of the constraint function inward but also increases its curvature (Figure 7). This means that given the same amount of total xylem water potential loss, increases in the distance between minor veins (decreases in minor vein density) will result in lower stomatal conductance and smaller leaf area, which consequently lowers the maximum NPP leaf that can be reached. For plants whose minor vein density is independent of leaf thickness, leaf thickness will have no effect on the leaf optima or NPP leaf from a hydraulic perspective. However, increasing leaf thickness may result in higher respiration costs and thus have a negative impact on NPP leaf .

F I G U R E 4
Change of (a) optimal leaf area and (b) stomatal resistance with increasing allowable loss of total xylem water potential F I G U R E 5 Effect of leaf width-to-length ratio on the isoline line of the constraint function. All four lines represent the same level of total loss of xylem water potential. With larger W/L ratio, the isoline line (dot line) is located in the northeast; it migrates southwest as the W/L ratio decreases F I G U R E 6 Change of total xylem resistance with leaf area of different leaf shape (W/L ratio), lines represent different W/L F I G U R E 7 Effect of leaf thickness on the location of the isoline line of the constraint function. The isoline lines have the same level of total loss of xylem water potential, but with different leaf thickness, Z, they differ in location and curvature The physiological explanation of our model is that with the same amount of transpiration per unit leaf area (e.g., with the same stomatal conductance per unit leaf area), a thicker leaf with lower minor vein density will require increasing the flow rate in the conduit, thus causing more loss of total xylem water potential. This effect is greater when petiole xylem water potential is high.
Increasing leaf thickness by increasing volume of mesophyll can decrease leaf mass per unit area (LMA; Enrique, Olmo, Poorter, Ubera, & Villar, 2016;Griffith, Quigley, & Anderson, 2016). Our model predicts that a thick leaf will also have higher stomatal resistance, which will result in lower maximum carbon assimilation rate per unit leaf area (and per unit leaf mass); therefore, a negative relation between LMA and maximum carbon assimilation rate results from changing leaf thickness. This negative relation between LMA and maximum carbon assimilation rate has been observed worldwide (Wright et al., 2017).

| Linking the leaf to the whole plant and landscape
The individual leaf is the fundamental unit of a whole tree.
Understanding the carbon balance of an individual leaf is the first step and is critical for understanding the whole plant's carbon balance. In this study, we have connected vein structure and leaf shape, which is further incorporated into the coupling of the transpiration and photosynthesis processes. This has provided us a to upscale biophysical processes occurring at the stoma scale to whole plants and, further, to the landscape.
Leaf shape and area also affect the light penetration and distribution of the canopy (Bonan, 1996;Dickinson, 1983;Sellers, 1985).
Smith, Sperry, and Adler (2016)  The hydraulic safety-efficiency trade-off-namely, plants whose vascular systems are more resistant to cavitation but inefficient at water transport-is a widely recognized phenomenon F I G U R E 8 Theoretical change of (a) leaf area, (b) net leaf carbon gain, and (c) net leaf carbon gain per unit leaf area with petiole xylem water potential due to safety-efficiency trade-off of plant hydraulic system and potential impact on plant distribution along hillslope from topographically induced soil moisture gradient (d) (Hacke, Sperry, Wheeler, & Castro, 2006;Manzoni et al., 2013;Markesteijn et al., 2011;McElrone, Pockman, Martínez-Vilalta, & Jackson, 2004). It is found that, in general, species with lower xylem conductance have safer xylem . Both hydraulic traits (pit and xylem network properties), and nonhydraulic traits (e.g., wood density) can contribute to the hydraulic safety-efficiency trade-off (Gleason et al., 2016). Among the hydraulic traits, pit properties have strong effects on the safety of xylem, while conduit diameter is a xylem network property found to be related to large P50 in xylem (safer xylem) in some species (Blackman, Brodribb, & Jordan, 2010;Brodribb, Bienaimé, & Marmottant, 2016;Scoffoni et al., 2017). For example, Blackman et al. (2010) have found the P50 in the leaves of 20 phylogenetically disparate woody angiosperm species from montane rainforest and dry sclerophyll forest in cool, temperate Australia to be significantly correlated with the average lumen width of the vein.
In our model, a narrowleaf has a small γ (the length ratio of two consecutive veins) and a small xylem conduit diameter and thus a higher total xylem resistance compared to a broadleaf with the same leaf area. Following the scaling relationship in WBE theory (West et al., 1997), small xylem conduits located at the end of the petiole will result in smaller xylem conduits throughout the en- Inefficiency in water transport also implies inefficiency in carbon fixation; thus, a narrowleaf plant will have an advantage when soil moisture is low but gradually lose that advantage with increasing soil moisture. So, the net primary productivities of narrowleaf plants increase with soil moisture more slowly than those of broadleaf plants. As a consequence, with their size increasing more slowly with increasing soil moisture, narrowleaf plants will gradually lose their advantage to broadleaf plants in wetter soils (Figure 8b, c).
Soil moisture is regulated by topography (Barling, Moore, & Grayson, 1994;Beven & Kirkby, 1979;Iorgulescu & Musy, 1997;O'Loughlin, 1986); it is lowest at ridge top and increases along the flow path, moving down the hillslope toward the channel. Thus, in general, along the hillslope, narrowleaf plants dominate the ridge area where soil moisture is low. Moving down along the hillslope, the soil becomes wetter, and broadleaf plants become more dominant ( Figure 8d).
Recently, the recognition of the interaction between plants and earth surface processes has led to interdisciplinary studies of the coevolution of plants and their physical environments (Corenblit, Steiger, Gurnell, & Naiman, 2009;Fisher, Heffernan, Sponseller, & Welter, 2007;Murray, Knaapen, Tal, & Kirwan, 2008;Reinhardt, Jerolmack, Cardinale, Vanacker, & Wright, 2010). However, an appropriate means to explicitly integrate the feedback mechanism between plant and earth surface processes remains unclear. The leaf model here can further facilitate these studies. Ding, Johnson, & Martin (2018) have developed a formula that can be used to directly connect groundwater regulated soil moisture to advection and diffusion erosion processes (i.e., landscape development processes) in low order watersheds. Once the leaf model is scaled to the whole plant and integrated into the advection and diffusion landscape processes, we can predict the abundance and distribution of a given plant throughout the landscape using the two erosion coefficients and the leaf shape.

| CON CLUS ION
We have connected leaf vein structure to leaf shape, which has provided the scaling relation of the total xylem resistance of veins with leaf area. This relation has been incorporated into the coupling of the transpiration and leaf carbon budget of a leaf model to investigate the impact of leaf shape on how plants coadjust stomatal resistance and leaf area to maximize NPP leaf . The optimality model predicts that both NPP leaf and leaf area increase with petiole xylem water potential, while stomatal resistance decreases, each at different rates. A narrowleaf has an overall lower NPP leaf and leaf area but higher stomatal resistance compared to a broad leaf. This occurs because a broadleaf has a larger vein length ratio, γ as suggested by our model.
With the same leaf area, a broadleaf will have, on average, larger xylem conduits and consequently a more efficient transportation of water and less resistance to cavitation than a narrowleaf. Our study indicates that when incorporating the impact of vein structure on the vulnerability curve, a trade-off exists between a higher marginal leaf carbon fixing rate with respect to xylem water potential and the ability of the leaf to resist cavitation. This explains how plants adapted to dry and/or hot conditions will gradually lose their advantages with increasing soil moisture. Using the link between vein structure and leaf shape, we can further connect plant hydraulic traits to other traits affecting strategies of leaf thermal regulation and canopy light utility, thereby providing a method to upscale from stoma to whole plants and the ecosystem.

ACK N OWLED G M ENT
This study was funded by a Discovery Grant from the Natural Sciences and Engineering Research Council of Canada (NSERC) to E.A. Johnson.

CO N FLI C T O F I NTE R E S T
The authors declare no conflicts of interest.