The legacy of C4 evolution in the hydraulics of C3 and C4 grasses

The anatomical reorganization required for optimal C4 photosynthesis should also impact plant hydraulics. Most C4 plants possess large bundle-sheath cells and high vein density, which should also lead to higher leaf hydraulic conductance (Kleaf) and capacitance. Paradoxically, the C4 pathway reduces water demand and increases water-use-efficiency, creating a potential mis-match between supply capacity and demand in C4 plant water relations. We use phylogenetic analyses, physiological measurements, and models to examine the reorganization of hydraulics in closely-related C4 and C3 grasses. Evolutionarily young C4 lineages have higher Kleaf, capacitance, turgor-loss-point, and lower stomatal conductance than their C3 relatives. In contrast, species from older C4 lineages show decreased Kleaf and capacitance, indicating that over time, C4 plants have evolved to optimize hydraulic investments while maintaining C4 anatomical requirements. The initial “over-plumbing” of C4 plants disrupts the positive correlation between maximal assimilation rate and Kleaf, decoupling a key relationship between hydraulics and photosynthesis generally observed in vascular plants.


Introduction 32
The evolution of C4 photosynthesis in the grasses-and the attendant fine-tuning of both 33 anatomical and biochemical components across changing selection landscapes [1,2,3] -likely 34 impacted leaf hydraulics and hydraulics-photosynthesis relationships, both within the grass 35 lineages in which C4 evolved independently > 20 times [4] , and as compared to closely-related 36 C3 [5,6] . C4 plants typically exhibit lower stomatal conductance (gs) and consequently greater 37 water-use efficiency than C3, because the concentration of CO2 inside bundle sheath cells permits 38 reduced intercellular CO2 concentrations and conservative stomatal behavior [7,8,9] . At the same 39 time, C4 plants require high bundle sheath to mesophyll ratios (BS:M), which are accomplished 40 with increased vein density and bundle sheath size as compared to C3 plants. In C3 species, leaf 41 hydraulic conductance (Kleaf) has a positive relationship with vein density [10,11,12,13] . The 42 decreased inter-veinal distance and consequently higher vein density in C4 species has been 43 predicted to lead to a higher Kleaf than closely-related C3 species [14,15] . Further, increased bundle 44 sheath size was proposed to lead to a higher leaf capacitance in C4 species [15,16], This would lead 45 to a potential physiological "mis-match", where the evolution of the C4 pathway simultaneously 46 increases a plant's hydraulic capacity while reducing its transpirational demand.

48
The significance of such a potential physiological mismatch depends on the potential costs and 49 tradeoffs associated with the building of an 'over-plumbed' leaf. If the costs are high [12,17] , then 50 one would expect to see a reduction of Kleaf over evolutionary time, as continued selection works 51 to optimize the C4 metabolism [5,18] . Alternatively, a maintenance of high Kleaf over time could 52 result from either a lack of strong selection to reduce Kleaf, or a strong evolutionary constraint 53 imposed by the anatomical requirements of C4 photosynthesis. In other words, the high BS:M 54 ratio required for an efficient C4 system may directly limit the ability of C4 plants to optimize 55 their hydraulic architecture. 56 57 The evolution of a new photosynthetic pathway that results in multiple potential changes to the 58 plant hydraulic system represents the ideal platform to expand our understanding of the 59 relationship between photosynthesis and water transport. It is generally thought that maximum 60 photosynthetic rate (Amax) and hydraulic capacity (Kleaf) are tightly linked, because the ability to 61 transport water through leaves to the sites of evaporation at a high rate allows for the 62 maximization of carbon gain. Studies have documented a positive correlation between Amax and 63 Kleaf across many scales, from a broad phylogenetic spectrum of species spanning vascular 64 plants [11] , to smaller clades of closely related species [13] . Grasses are largely absent from previous 65 efforts to examine this relationship, which is unfortunate because of the parallel venation found 66 in grasses and other monocots. With over 20 origins of C4 photosynthesis with ages that span ~ 67 30 million years, grasses also present a unique opportunity to examine the influence of C4 68 evolution on Amax-Kleaf relationships. Using a broad sampling of grasses ( Fig. 1 leaf turgor loss point, Amax and lower or equivalent gs than their closest C3 relatives (Fig. 2). The 80 one C3-C4 intermediate species, Steinchisma decipiens, in our analysis had Kleaf similar or 81 equivalent to C4, but leaf capacitance, leaf turgor loss point, gs and Amax equivalent to C3 (Fig. 2). 82 By analyzing our data in the context of the evolutionary models (Supplementary Table S1 We also looked for evolutionary trends in hydraulic traits after the evolution of a C4 system to 93 probe for an extended 'optimization' phase of C4 evolution [3,20] . Identifying directional trends in 94 continuous character evolution is difficult without fossil taxa, and it is impossible to directly 95 measure hydraulic traits for fossils; however, we can test for trends indirectly using extant 96 species. For example, if reduction in Kleaf is selected for subsequent to C4 evolution we expect 97 older C4 lineages to have lower Kleaf values than younger C4 lineages. We extracted the 98 evolutionary age of C4 origin for each of our lineages from the dated phylogeny [19] . Regressions 99 of evolutionary age versus hydraulic traits provide strong evidence for a long-term directional 100 trend in hydraulic evolution following the origin of C4 photosynthesis (Fig. 3). Kleaf, leaf turgor 101 loss point and capacitance showed significant negative correlations with evolutionary age, while 102 Amax had a significant positive correlation. In contrast, there was no significant relationship 103 between gs and evolutionary age. No evolutionary relationships were detected in C3 species, 104 which indicated the correlations between evolutionary age and hydraulic traits were unique to C4 105 species. We also tested for an evolutionary trend by modelling hydraulic trait evolution using a 106 phylogeny with branch lengths scaled to molecular substitutions/site, which provides an estimate 107 of differences in evolutionary rates between lineages [4] . While the second approach requires 108 many assumptions that are likely violated, the results also provide additional support to a 109 directional trend in Kleaf and capacitance in C4 lineages: comparing 12 different types of models 110 with or without evolutionary trends (supplementary Table S7), we found Kleaf and leaf 111 capacitance were best fitted by the Brownian motion model with a significant negative trend for 112 C4 (Supplementary Table S8, Table S9    stomatal conductance (gs) and maximal assimilation rate (Amax) vs. the evolutionary age for the 153 nine origins of C4 to show the evolutionary trend within C4 and within their closely-related C3 154 species. The evolutionary age for each sampled origin is derived from the dated phylogeny [19] . 155 156 157   We used our mechanism-based physiological model [32] to consider how the evolution of higher 164 Kleaf would affect the optimal gs and photosynthesis in C3 and C4 plants. An increase in Kleaf in 165 the C3 ancestor selects for higher gs and increases the steady-state leaf water potential to a 166 limited extent (Fig. 5, S1). Changing Kleaf has a smaller effect on the photosynthesis rate of C4 167 than that of C3 (Fig. 6, Table S25), Decreasing Kleaf by half or doubling it changes the 168 photosynthesis rate of a C4 plant by an average of -4.27% and 3.48%, respectively. In contrast, 169 the same shifts in Kleaf has average effects of -10.07% and 9.14% on the assimilation rate of a C3 170 plant. The sensitivity of the assimilation rate to changes in Kleaf decreases with increasing CO2 171 concentration and increasing water-limitation for both C3 and C4 plants (Table S25). These 172 differences in sensitivity to Kleaf were robust to differences in physiological properties between 173 C3 and C4 (specifically, the temperature response properties and Jmax/Vcmax ratio; Table S25). The 174 assimilation rate of C4 plants was still less sensitive to Kleaf than that of C3 species under 175 different CO2 concentration and water-limited conditions (Table S25). The physiological 176 modeling results indicates that C4 species maintain lower gs and higher leaf water potential 177 compared to closely related C3 species because the CCM reduces transpirational demand. The 178 modeling effects of varying Kleaf on photosynthesis confirmed the diminished returns for high-179 efficiency water transport in C4 species mentioned above.  The evolution of the C4 pathway in the grasses caused a series of shifts in hydraulic properties as 218 compared to closely-related C3 grasses. The anatomical requirements of C4 initially increased 219 Kleaf and leaf capacitance, as predicted by previous studies [14,15,16] ; however, Kleaf and leaf 220 within a phylogenetic framework when comparing multiple species [21,22] , and phylogenetic 224 studies have assumed trait evolution as simple Brownian motion [23,24] . Hydraulic traits, however, 225 may have evolved along different trajectories before and after the evolution of the C4 pathway 226 and associated anatomical reorganization, resulting in more complicated evolutionary dynamics. 227 Our evolutionary models indicated C4 grasses initially had higher Kleaf, leaf capacitance, turgor 228 loss point than corresponding C3, and a lower stomatal conductance (gs) than grasses consistent 229 with previous studies [25,26] . Decreased vein distance and increased bundle sheath size are thought 230 to be anatomical precursors to the evolution of C4 [27,28] , and both are thought to increase Kleaf 231 and/or leaf capacitance [14,15] . Therefore, the shifts of Kleaf and leaf capacitance likely occurred 232 before, or at the initial formation of, the C4 CCM. After the full formation of C4, Kleaf and/or leaf 233 capacitance started to decrease, which led to higher or equivalent Kleaf and leaf capacitance in the 234 current C3 and C4 species (Fig. 2) Amax across and within plant taxa [11,13,31] . We found that Amax and Kleaf are positively correlated 250 in our C3 species but not in C4 (Fig 4). Ocheltree et al. (2016) [22] similarly found no relationship 251 between Kleaf and Amax in a set of nine C4 species. We see possible explanations that are not 252 necessarily mutually exclusive. First, the positive relationship of Amax and Kleaf is weakened 253 under high Kleaf, possibly due to diminished returns of further increasing the efficiency of water 254 transport [11 ,31] , a conclusion supported by our physiological modeling results below. As Kleaf 255 tends to be lower in grasses than in other species, it is possible that the diminishing returns from 256 increasing Kleaf manifest at lower values in grasses, and the initial high Kleaf resulting from C4 257 anatomy could be in the Amax "saturation" zone. Lastly, we see evidence here that the time-since-258 C4-evolution affects several hydraulic traits across and within lineages, and it could be that a 259 walk towards Amax-Kleaf optimality is slowly occurring within C4 grass lineages in relatively 260 newfound ecological niches. However, the similar correlations of gs vs. Amax in C3 and C4 and 261 lack of evolutionary trend in gs indicated the evolutionary processes of gs might be already near 262 the optimal condition or stabilized quickly. Other hydraulic traits of leaf capacitance and leaf 263 turgor loss point do not seem to contribute to the Amax directly because of weak correlations. 264

265
We identified the mode and direction of evolution for hydraulic traits in C3 and C4 lineages and 266 found evidence that different traits followed different evolutionary processes. Hydraulic 267 conductance and leaf capacitance could therefore evolve with directions in a step-wise fashion 268 due to anatomical constraints, but gs and leaf turgor loss point might have a more quick process 269 of readjustments, which allows them to stabilize soon. This suggests that there could be greater 270 diversification of Kleaf and leaf capacitance in the existing C4 species and maybe in the future. 271 Also, these rearrangements of hydraulic properties interacted with each other throughout the 272 evolutionary trajectory. For example, increased Kleaf and leaf capacitance would lead to an 273 increased water transport efficiency, which enabled greater gs of the C4 ancestor (either a C3 274 grass or a C3-C4 intermediate), but the formation of the full C4 CCM enables a decrease of gs. 275 Therefore, observed gs in C4 grasses reflects a balance of these two contrasting physiologies 276 playing out in a given ecological and phenological background, which may explain why although 277 C4 gs was lower than the C3, the difference was not large. This line of reasoning might also 278 explain the inconsistent observations of gs comparisons between C3 and C4. Most previous 279 studies found that C4 grasses had lower gs than C3 grasses in both closely related and unrelated 280 species [25,33] , yet Taylor et al. (2014) found that C4 grasses maintained a higher or equivalent gs 281 to closely-related C3 grasses [34] . Likewise, artificial selection or genetic engineering might have 282 more success in adjusting these hydraulic traits in advance. Consciously selecting or 283 manipulating narrower xylem, decreasing the expression of aquaporins, or other mechanisms of 284 decreasing leaf conductance while maintain high bundle sheath to mesophyll ratio, together with 285 CCM may increase the water use efficiency of C4 species further. Our phylogenetic analyses can 286 thus inform both the evolutionary history of C4 plants and future efforts to modify C4 crops. high Kleaf are still significant in C4 plants [12,35,36,37,38] . The most efficient way to reduce Kleaf costs 298 would be to reduce venation density, as veins come with high construction costs [12,17] , and also 299 reduce the leaf area that is available for carbon fixation. Yet the anatomical requirements of the 300 C4 system preclude this option: reducing vein density would result in a highly inefficient C4 301 system [15] , which would negatively impact the plant's carbon budget, presumably to a much 302 greater extent than the cost of an overbuilt venation system. As vein construction is a primary 303 contribution to the cost of a high Kleaf, and high vein densities are now linked to a new function 304 (C4 carbon fixation), the cost-benefit calculations in optimizing Kleaf have shifted, and the 305 tradeoff is in favor of overplumbing in order to maintain a highly efficient new carbon fixation 306 system. In evolutionary vocabulary, what emerges is a new constraint -and in this example, it is 307 clear that the emergence of a new constraint to organismal evolution is simply due to a shift in 308 the tradeoffs associated with characters that influence multiple aspects of organismal function. In 309 other words, we assume a low vein density is a phenotype that is still developmentally 310 achievable for C4 grasses; what has prevented its emergence is the shift in functional costs 311 associated with reduced vein densities. 312 313 And yet, we documented a gradual reduction in Kleaf over time, which we presume was 314 accomplished via changes in other factors that influence leaf hydraulic capacity-perhaps by 315 changing xylem conduit diameters, shifts in extra-xylary mesophyll conductance, decreased 316 expression of aquaporins, and reorganization of internal air spaces [6,12,37,39,40] . It is possible that 317 these changes resulted from a continued and direct selection pressure to reduce investment in an 318 underutilized hydraulic system. An alternative explanation is that all of the traits that influence 319 Kleaf also play important roles in other aspects of leaf function -and the emergent of a new 320 constraint (a high vein density to maintain C4 function) has released still other constraints on 321 other traits so that they may be optimized for their other functions. A striking pattern in our data 322 is that older C4 lineages have achieved both lower Kleaf and higher Amax -suggesting that they 323 are continuing to optimize their photosynthetic capacity, long after the initial origin of C4. We 324 suspect that the slow evolutionary decline in Kleaf is due in large part to the optimization of traits 325 to increase Amax at the expense of Kleaf, which is possible only because hydraulic capacity was 326 already "buffered" by the vein density requirements of C4 -allowing for continued reductions of 327 Kleaf at no functional cost. Increased suberization of bundle sheath cells is one example of a 328 potential release of constraint [22] : it allows C4 plants to gain higher Amax through reducing bundle 329 sheath leakiness, but it likely simultaneously reduces water flow from veins out into the 330 mesophyll. Since C4 plants are already operating in hydraulic excess, bundle sheath suberization 331 may be optimized for C4 function without any negative repercussions for plant water relations. 332 This hypothesis could also explain the opposing trends in Amax and Kleaf when viewed as a 333 function of evolutionary age. The examination of C4 evolution in grasses provides an exciting 334 system to study the evolutionary dynamics of constraints highlighted by the interplay between 335 photosynthesis and plant hydraulics.

Hydraulic traits 357
Leaf hydraulic conductance (Kleaf) was measured using the evaporative flux method [41] , with 358 some adjustments to maintain stability of the evaporative environment to which the leaf was 359 exposed (Supplementary Methods). The evening before measurements, potted plants were 360 brought to the laboratory, watered, and then covered by black plastic bags filled with wet paper 361 towels to rehydrate overnight. For the leaf gasket, a 1 cm diameter, ~ 1 cm long solid silicone 362 rubber cylinder was cut nearly in two, leaving a hinge on one end. The cylinder was placed 363 around the leaf blade near the ligule and glued shut with superglue [42] . The leaf was cut from the 364 plant with a razor blade while submerged in a 15 mmol L -1 KCl solution; the rubber gasket was 365 then attached to tubing filled with the same KCl solution. The other end of the tubing was inside 366 a graduated cylinder that sat on a digital balance (Mettler-Toledo). The leaf was then placed 367 inside a custom, environmentally controlled cuvette that allowed for the measurement of entire 368 grass blades. Throughout measurements, cuvette temperature was controlled at 25 o C and the 369 humidity was 55-65% (VPD range of 1.1-1.4 kPa) across measurements, but remained constant 370 during a particular measurement. Photosynthetically active radiation in the system is 1000 µmol 371 m -2 s -1 . Flow from the balance was monitored for 45 m to 1h until the flow rates reach steady 372 state. After the measurements, the leaf was detached and was put into a plastic bag to equilibrate 373 for 20 minutes to measure the leaf water potential (Model 1000, PMS Instrument, USA). Kleaf 374 values were further standardized to 25 o C and leaf area to make the Kleaf comparable among 375 studies and across species. Data indicating a sudden change of flow and whose leaf water 376 potential was an obvious outlier were deleted. 377

378
We measured pressure-volume (PV) curves for six leaves per species using the bench-drying 379 method [43,44] . A leaf was cut directly from the same plants rehydrated in the lab (as described 380 above) using a razor blade and leaf water potential was measured immediately. Then, the leaf 381 weight was recorded. The leaf was initially allowed to dry on the bench for 2-minute intervals 382 and put into a ziplock bag and under darkness for 10-minute equilibration before measuring the 383 leaf water potential and leaf weight again. Then, the waiting intervals could be adjusted based on 384 the decrease of the leaf water potential (from 2 minutes-1h). Ideally, a decreasing gradient of -385 0.2MPa for leaf water potential was obtained for the curves, until the leaf weight reached a 386 steady state. At the end of the experiment, leaves were dried in the oven at 70 o C for 48h to obtain 387 the dry weight. The PV curves were used in curve fitting to obtain leaf capacitance, and leaf 388 turgor loss point using an excel program from Sack and Pasquet-Kok (2010) [

Phylogenetic analysis 401
Phylogenetic analysis for C3 and C4. We pruned the dated phylogeny from a published grass 402 phylogeny to include only the species in our physiological experiments [19] (Fig. 1). Using the 403 dated phylogeny, for each of the hydraulic traits, we fitted evolutionary models to test which 404 evolutionary model best explains observed distribution of traits along the phylogeny and how 405 these models differ between C3 and C4 (Table S1). We fitted evolutionary models belonging 406 Brownian Motion model and Ornstein-Uhlenbeck Model using the package "mvMORPH" in 407 R [47] . To determine the best fitted evolutionary model, we compared two criteria, the small-408 sample-size corrected version of Akaike information criterion (AICc, the lower AICc, the better 409 fit) and Akaike weights (AICw, the higher AICw, the better fit) [ C4 origin from the dated phylogenetic trees. Then, we regressed the hydraulic traits with 419 evolutionary age. A significant negative correlation between evolutionary age and hydraulic trait 420 will indicate a further decreasing evolutionary direction after C4 evolved. We also performed an 421 additional analysis to test the original states and further direction together. We extracted 422 molecular phylogeny for all the species from Edwards, GPWG II (2012) [4] . Except for the six 423 evolutionary models mentioned above, the molecular phylogeny allows us to fit for additional six 424 Brownian motion models with trend (Supplementary Table S7). Likewise, if Brownian motion 425 model with trend fits the phylogenetic patterns better than Brownian motion model without trend 426 it means there is an evolutionary trend, and a significant LRT test for a two-trend model suggests 427 that C3 and C4 lineages differ in the speed or direction of hydraulic evolution. We also mapped 428 the traits on the phylogeny for potential further references (Fig. S2-S5). 429 To further test whether there are significant differences among C4 subtypes, evolutionary models 430 with subtypes (Table S1) were used to fit the data. We again used AICc, AICw and LRT 431 methods to find the best model variants: whether there are significant differences for hydraulic 432 shifts and evolutionary trends among three different subtypes. For the leaf capacitance analysis, 433 Dichanthelium clandestinum is deleted as it is an obvious outlier. 434 Phylogenetic analysis for correlations among traits. Multivariate analysis in "mvMORPH" 435 was used to estimate the correlations between Amax and each of the hydraulic traits and to test the 436 hypotheses that whether such correlations are different between C3 and C4. The process of 437 brownian motion with different root for C3 and C4 was used for Kleaf, gs and leaf turgor loss and 438 brownian motion with the same root was used for leaf capacitance. Since the Ornstein-439 Uhlenbeck process is difficult to take the root state difference into consideration, here we used 440 Brownian motion assumptions as approximation for leaf turgor loss. Seven different correlation 441 models are fitted (Table S19). We used LRT for the seven correlation models to test whether the 442 correlation of the two traits is significantly different from 0 and whether the correlation of two 443 traits is significantly different between C3 and C4. Such correlation analysis is similar to PGLS 444 considering C3 and C4, but with more varieties on the setting of variance and covariance matrix. Furthermore, we used physiological models that couples the photosynthesis systems and 449 hydraulic systems to predict the effect of changing Kleaf on assimilation rate [32] . The change of 450 Kleaf was assumed to change the plant hydraulic conductance (Kplant) proportionally in the 451 modeling process. We double or reduce by half Kleaf relative to the original value to predict the 452 effects on assimilation rates for C3 and C4 pathways. We assumed C4 had the same 453 photosynthetic properties with C3 species (e.g., Rubisco affinity and specificity, Supplementary 454 Table S24) other than the carbon concentration mechanism, which mimics the initial evolution of 455 C4 and the closely-related C3-C4 system. We also model the additional scenarios in which C4 had 456 different photosynthetic properties to support the above condition further (Supplementary