Quantifying in situ phenotypic variability in the hydraulic properties of four tree species across their distribution range in Europe

RESEARCH ARTICLE Quantifying in situ phenotypic variability in the hydraulic properties of four tree species across their distribution range in Europe N. Gonza ´lez-Muñoz 1 *, F. Sterck 2 , J. M. Torres-Ruiz 1 , G. Petit 3 , H. Cochard 4 , G. von Arx 5,6 , A. Lintunen 7 , M. C. Caldeira 8 , G. Capdeville 1 , P. Copini 2,9 , R. Gebauer 10 , L. Gro ¨nlund 7 , T. Ho ¨ltta ¨ 7 , R. Lobo-do-Vale 8 , M. Peltoniemi 11 , A. Stritih 12 , J. Urban 10 , S. Delzon 1 1BIOGECO, INRA, Universite ´de Bordeaux, Pessac, France, 2Forest Ecology and Forest Management Group, Wageningen University & Research, Wageningen, The Netherlands, 3Universitàdegli Studi di Padova, Dep. TeSAF, Legnaro (PD), Italy, 4PIAF, INRA, Universite ´Clermont-Auvergne, Clermont-Ferrand, France, 5Swiss Federal Institute for Forest, Snow and Landscape Research WSL, Birmensdorf, Switzerland, 6Climatic Change and Climate Impacts, Institute for Environmental Sciences, Geneva, Switzerland, 7Department of Forest Sciences, University of Helsinki, Helsinki, Finland, 8Forest Research Centre, School of Agriculture, University of Lisbon, Tapada da Ajuda, Lisboa, Portugal, 9Wageningen Environmental Research (Alterra), Wageningen, The Netherlands, 10 Department of Forest Botany, Dendrology and Geobiocoenology, Mendel University, Zemědělska ´3, Brno, Czech Republic, 11 Natural Resources Institute Finland (Luke), Latokartanonkaari 9, Helsinki, Finland, 12 Swiss Federal Institute of Technology ETH, Planning of Landscape and Urban Systems, Zurich, Switzerland *[email protected] Abstract Many studies have reported that hydraulic properties vary considerably between tree species, but little is known about their intraspecific variation and, therefore, their capacity to adapt to a warmer and drier climate. Here, we quantify phenotypic divergence and clinal variation for embolism resistance, hydraulic conductivity and branch growth, in four tree species, two angiosperms (Betula pendula,Populus tremula) and two conifers (Picea abies, Pinus sylvestris), across their latitudinal distribution in Europe. Growth and hydraulic efficiency varied widely within species and between populations. The variability of embolism resistance was in general weaker than that of growth and hydraulic efficiency, and very low for all species but Populus tremula. In addition, no and weak support for a safety vs. efficiency trade-off was observed for the angiosperm and conifer species, respectively. The limited variability of embolism resistance observed here for all species except Populus tremula, suggests that forest populations will unlikely be able to adapt hydraulically to drier conditions through the evolution of embolism resistance. Introduction Massive forest mortality events due to drought stress and rising temperatures have been observed at the global and regional scales [1–5]. Considering that climate change models predict further increases in mean temperature and in the frequency and severity of extreme drought events [6], more negative impacts on tree survival are expected [7]. In this context, PLOS ONE | https://doi.org/10.1371/journal.pone.0196075 May 1, 2018 1 / 17 a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 OPEN ACCESS Citation: Gonza ´lez-Muñoz N, Sterck F, Torres-Ruiz JM, Petit G, Cochard H, von Arx G, et al. (2018) Quantifying in situ phenotypic variability in the hydraulic properties of four tree species across their distribution range in Europe. PLoS ONE 13(5): e0196075. https://doi.org/10.1371/journal. pone.0196075 Editor: Berthold Heinze, Austrian Federal Research Centre for Forests BFW, AUSTRIA Received: July 5, 2017 Accepted: April 5, 2018 Published: May 1, 2018 Copyright: ©2018 Gonza ´lez-Muñoz et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited. Data Availability Statement: All data underlying the study are available at the Dryad Digital Repository (https://doi.org/10.5061/dryad. b2pg468). Funding: This article is based on work from COST Action FP1106 STReESS, supported by COST (European Cooperation in Science and Technology). This study was funded in part by the “Investments for the Future” programme (grant no. ANR-10-EQPX-16, XYLOFOREST) from the French assessments of the properties associated with drought resistance in trees, and of the capacity of species to deal with environmental changes, may help us to anticipate the impact of climate change on forest tree species. Hydraulic failure due to xylem embolism is now considered one of the main causes of drought-induced tree mortality [8–10]. When soil water potential drops due to water shortage, the tension of the xylem water column increases, promoting the formation of embolisms that reduce the hydraulic functioning of the plant [11,12]. In cases of prolonged drought, soil water potential continues to fall, triggering the spread of embolisms throughout the xylem conduit network, leading to the hydraulic dysfunction of the plant vascular system and, finally, to lethal damage to the plant [13,14,15]. Therefore, determining the resistance to embolism of the species is crucial for evaluating the consequences that the expected increase in drought event frequency can have on a given population, forest or biome. P 50 is the xylem pressure at which 50% of conductivity is lost due to embolism formation, and it is widely used to assess plant hydraulic safety to embolism. Xylem-specific hydraulic conductivity (K S ), i.e. the rate of water transport through a given area of sapwood per unit pressure difference and per unit length, is commonly used to assess hydraulic efficiency [12]. Across species, literature shows a weak correlation between hydraulic safety and hydraulic efficiency, but the absence of species displaying both high hydraulic efficiency and safety suggests a possible safety-efficiency trade-off [16]. In conifers, P 50 and K s are only weakly correlated, as embolism resistance is driven mostly by the torus-aperture overlap in pit pairs [17–19], whereas xylem hydraulic efficiency is not influenced by this pit trait. By contrast, in angiosperms, both P 50 and K s are associated with pit membrane structure [20–22] and thickness [23], as well as with the perforation structure [24,25]. Differences in resistance to embolism, i.e. in P 50 , across species have been widely reported [17,19,26,27,28]. However, less attention has been paid to within-species phenotypic variation in this hydraulic property. Phenotypic variability results from a combination of genetic variation (differences in genotype among different individuals within the population and between populations) and phenotypic plasticity (genotype property to render different phenotypes in different environments [29]), and defines the capacity of populations to succeed under changing environmental conditions [30,31]. Low levels of phenotypic variability across large spatial scales may indicate a low potential of species to adapt to ongoing climate change. Contrary to other plant functional properties (see for instance [32] for leaf phenology, [33] for leaf functional traits), previous works on hydraulic properties show that phenotypic differences within species are by far lower than those found across species [34–38], and these differences are even smaller in gymnosperms than in angiosperms [39]. For instance, Lamy et al. 2011, 2014 found neither phenotypic variability in situ and nor genetic differentiation between maritime pine populations and suggested uniform selection rather than genetic drift, for P 50 . However, whether the phenotypic variation of hydraulic properties varies across species distribution ranges remains largely unexplored. Furthermore, studies assessing the extent to which phenotypic variability in hydraulic properties is lower than that of other key species traits over large scales are also lacking. The main aim of this study was to evaluate phenotypic variability in the functional hydraulic safety and efficiency (P 50 and K S ) of four European tree species (two conifers and two angiosperms) along a latitudinal gradient covering most of their distribution range. We assessed the capacity of the species to adapt to changing environmental conditions by exploring the links between hydraulic properties and latitude and climate. We also evaluated the phenotypic variability of branch growth in trees from the same populations to assess the extent to which hydraulic properties were conserved relative to other key traits. Finally, we assessed the safetyefficiency trade-off at intraspecific level. We hypothesized 1) a weaker phenotypic variation for Phenotypic variability in hydraulic traits PLOS ONE | https://doi.org/10.1371/journal.pone.0196075 May 1, 2018 2 / 17 National Agency for Research and the Cluster of Excellence COTE (ANR-10-LABX-45, within the DEFI project) to SD. NGM was supported by the Agreenskills+ Fellowship Programme, which has received funding from the EU’s Seventh Framework Programme under grant agreement No. FP7-26719 (Agreenskills+ contract). GvA was supported by a grant from the Swiss State Secretariat for Education, Research and Innovation SERI (SBFI C14.0104). RLdV was funded by a postdoctoral fellowship from the Portuguese Fundac¸ão para a Ciência e a Tecnologia (FCT; SFRH/BPD/86938/2012). RG and JU were supported by the Ministry of Education, Sports and Youth of the Czech Republic (COST-LD13017). MP was supported by EU Life Programme (LIFE12 ENV/FI/000409). Forest Research Centre (CEF, School of Agriculture, University of Lisbon) is a research unit funded by FCT (UID/AGR/00239/ 2013). Competing interests: The authors have declared that no competing interests exist. hydraulic safety than for hydraulic efficiency and branch growth, given the highly conserved evolutionary nature of P 50 [37,38]; 2) a phenotypic cline -a gradual change of a phenotypic character in a species over a geographical areain both hydraulic safety and efficiency and 3) a weak safety-efficiency trade-off within species. This study provides for the first time a multispecies assessment of inter and intra-specific phenotypic variability in functional hydraulic properties along a large latitudinal gradient. Our results will help to characterize the adaptive capacities of European forests, which will have to face drier and warmer climatic conditions in the future. Materials and methods Study species and populations We focused on four widely distributed European species, with different water-transport structures, from diffuse porous with scalariform perforation plates (Betula pendula Roth) or with simple perforation plates (Populus tremula L.) to softwood (two tracheid-bearing species, Picea abies (L.) Karst and Pinus sylvestris L.). For each species, four to six populations were selected across their distribution range (see the distribution range of the species and the location of the populations in S1 Fig). The mean annual temperature and total annual rainfall across selected populations ranged from -1.8 to 9.5˚C and 538 to 1739 mm, respectively (S1 Table). We also selected two different sites a few kilometres apart, for each population. Climatic data Data for mean annual temperature (MAT) and total annual precipitation (MAP) were obtained from WorldClim original 30-s data (http://www.worldclim.org/bioclim) [40] downscaled to 100-m resolution based on a high-resolution digital elevation model (DEM) and moving window regression technique [41] for all but the Italian (IT) and Swiss populations (SW-LOE and SW-PFY). The MAP data for the IT (Italy), SW-LOE (Switzerland-Loetschental) and SW-PFY (SwitzerlandPfynwald) populations were obtained from nearby weather stations at San Vito di Cadore (Centre for Alpine Environment Studies) and Sierre (www. meteoswiss.ch), respectively, due to considerable variations in topography. The aridity index (AI) was calculated as MAP/PET (total annual precipitation/annual potential evapotranspiration). PET was extracted from the Global Aridity and PET Database (http://www.cgiar-csi. org). We averaged the mean temperatures (T_Sum) or aridity indices (AI_Sum) of June, July and August to obtain mean values for the summer (see S1 Table for the climatic conditions of the populations studied). Xylem vulnerability to embolism We collected branches from five to 11 healthy mature trees per population in the early morning during the wet season (spring 2015). One or two branches with three to five functional rings were sampled at mid-crown and south oriented. Samples had a standard length of 45 cm. Transpiration losses were prevented by removing the leaves or needles immediately after sampling and wrapping the branches in moist paper to keep them humid and cool (3˚C) until the measurement of embolism resistance (within three weeks of sampling). The bark was removed from conifer branches to prevent resin to fill the cavitron reservoirs (see below, [17]), and all branches were recut with a razor blade, under water, to a standard length of 0.27 m. For each angiosperm species, 10 samples per species were used to test the open vessel artefact [42] by injecting air at 2 bars at one end and no open vessels were detected for any of them. Phenotypic variability in hydraulic traits PLOS ONE | https://doi.org/10.1371/journal.pone.0196075 May 1, 2018 3 / 17 Vulnerability to drought-induced embolism was determined at the Caviplace (University of Bordeaux, Talence, France; http://sylvain-delzon.com/caviplace) and INRA-Clermont-Ferrand facilities, with the Cavitron technique [43,44]. Samples were infiltrated with a reference ionic solution of 10 mm 25 KCl and 1 mm CaCl 2 in deionized ultrapure water. Centrifugal force was used to generate negative pressure into the xylem and induce cavitation. This method allows to measure xylem conductance under negative pressure using the custom software Cavisoft 4.0 (Univ. Bordeaux, Pessac, France). Initially, the maximum conductance of stem (K max , in m 2 MPa -1 s -1 ) was calculated under low xylem pressures. The percentage loss of conductance (PLC) of the stems was calculated at different xylem pressures (P i ) from -0.8 to -5 MPa with the following equation: PLC ¼100 1 K Kmax   We obtained one vulnerability curve per tree by measuring one or two of the collected branches. These vulnerability curves show the percentage loss of xylem conductance as a function of xylem pressure [17]. For each branch, the relationship between PLC and xylem water pressure was fitted with the following sigmoidal equation [45]: PLC ¼100 1þexp S 25 ðPiP50Þ    where P 50 (MPa) is the xylem pressure inducing a 50% loss of conductivity and S(% MPa -1 ) is the slope of the vulnerability curve at the inflection point. All sigmoidal functions were significant and fitted with the NLIN procedure in SAS (version 9.4 SAS Institute, Cary, NC, USA). The xylem-specific hydraulic conductivity (K s , kg m 1 MPa -1 s -1 ) was calculated by dividing the hydraulic conductivity measured at low speed by the sapwood area of the sample. Branch growth measurements We also collected one branch per tree from three to five trees per population and site. We selected straight branches and did not keep any sample with reaction wood for our measurements. The allometric relationship between branch radius (mm) and xylem age (number of years) was used as a surrogate for tree growth, as radial branch growth and tree growth patterns are highly correlated [46]. The branch surface area and the number of tree rings were systematically measured at 70 cm from the branch apex. Statistical analyses We assessed the phenotypic variability of functional hydraulic properties (P 50 and K s ) and branch growth (branch radius/xylem age) in each species, by testing the effect of population and site with nested ANOVAs, in which the population and the site nested in population were considered factors. If statistically significant differences were observed, post-hoc Tukey tests were conducted for multiple comparisons between populations. Before running the ANOVAs, we checked that the data satisfied the assumptions of normality and homoscedasticity. As vessel size can rapidly increase with branch size during early years of tree growth, and then may have a potential effect on hydraulic conductivity, we tested any potential correlation between K s and branch diameter. We also calculated the inter-population and intraspecific coefficients of variation (% CV inter and CV sp , respectively). For each species, Spearman´s or Pearson´s correlation coefficients (depending on the linearity condition) were calculated between the averaged by site P 50 ,K s and branch growth and the latitude and the five climatic variables Phenotypic variability in hydraulic traits PLOS ONE | https://doi.org/10.1371/journal.pone.0196075 May 1, 2018 4 / 17 mentioned above. Finally, we also checked for intraspecific safety-efficiency trade-offs, with Spearman´s or Pearson´s correlation tests. Statistical analyses were performed with the R project for statistical computing (R Development Core Team, 2016) [47]. Results Phenotypic variability across species distribution ranges Xylem vulnerability curves followed a sigmoid function in all species (Fig 1,S2 Fig), showing the lack of an open vessel artefact and the accuracy of the results obtained here. Betula pendula and Picea abies showed, respectively, the lowest and highest resistance to embolism of the four species evaluated. The mean P 50 ±SE (MPa) was -1.78 ±0.02 for Betula pendula, -2.45 ±0.08 for Populus tremula, -3.16 ±0.03 for Pinus sylvestris, and -3.58 ±0.02 for Picea abies (Fig 1). Differences in P 50 between populations were observed for all species, whereas the site (nested in population) had an effect on P 50 in all species but Picea abies (Table 1,Fig 2). The CV in P 50 was low for all species other than Populus tremula. The variability in P 50 of Betula pendula,Picea abies and Pinus sylvestris ranged from 4.15 (CV inter of Picea abies) to 10.23% (CV sp of Pinus sylvestris), whereas that of Populus tremula ranged from 24.82 (CV inter ) to 25.07% (CV sp ) (Table 2). The high variability observed for Populus tremula was mostly due to the population of Finland-Va¨rrio¨(FI-VA), which had the least negative P 50 values of any of the populations studied (Fig 2). Fig 1. Xylem vulnerability curves for each population of the four species studied (Betulapendula,Populus tremula,Picea abies and Pinus sylvestris). The shaded band represents the standard deviation. CR: Czech Republic; Fi-RU: Finland-Ruotsinkyla¨; NE: The Netherlands; PO: Portugal; Fi-HYY: Finland-Hyytia¨la¨; Fi-VA: Finland-Va¨rrio¨; IT: Italy; SW-LOE: Switzerland-Loetschental; SW-PFY: SwitzerlandPfynwald. https://doi.org/10.1371/journal.pone.0196075.g001 Phenotypic variability in hydraulic traits PLOS ONE | https://doi.org/10.1371/journal.pone.0196075 May 1, 2018 5 / 17 Variability levels were much higher for K s than for P 50 (Table 2,Fig 2). K s differed significantly between populations, for all species other than Populus tremula (Table 2,Fig 2). For this Table 1. Effects of population and site on P 50 (MPa), xylem-specific hydraulic conductivity (K s , kg m -1 MPa -1 s -1 ) and branch growth (BG, estimated as branch radius/xylem age (mm/year)) of study species, according to nested ANOVAs. The F, p-values and degrees of freedom are shown. Pop: population. Angiosperms Conifers Betula pendula Populus tremula Picea abies Pinus syvestris df F p df F p df F p df F p Population 3 9.978 <0.001 3124.885 <0.001 57.971 <0.001 515.284 <0.001 P 50 Site (Pop) 4 7.894 <0.001 49.905 <0.001 6 1.042 0.405 6 6.505 <0.001 Population 3 34.949 <0.001 2 0.249 0.781 5 31.381 <0.001 54.810 <0.001 K s Site (Pop) 4 1.022 0.405 3 2.809 0.051 6 1.743 0.123 6 2.019 0.070 Population 3 10.957 <0.001 312.999 <0.001 537.755 <0.001 550.702 <0.001 BG Site (Pop) 4 2.140 0.107 4 1.021 0.416 6 5.539 <0.001 614.674 <0.001 https://doi.org/10.1371/journal.pone.0196075.t001 Fig 2. Mean P 50 (MPa), xylem-specific hydraulic conductivity (K s , kg m -1 MPa -1 s -1 ) and branch growth (mm/year) per species, population and site. The two sites are represented in different colours (white and grey). The bars represent the nominal range of data variation, with the upper and lower ends showing the upper quartile plus 1.5 times the interquartile range and the lower quartile minus 1.5 times the interquartile range, respectively. Values beyond these limits are plotted as circles. CR: Czech Republic; PO: Portugal; NE: The Netherlands; Fi-RU: Finland-Ruotsinkyla¨; Fi-VA: Finland-Va¨rrio¨; Fi-HYY: Finland-Hyytia¨la¨; IT: Italy; SW-LOW: SwitzerlandLoetschental; SW-PFY: SwitzerlandPfywald. Different letters indicate statistically significant differences between populations. https://doi.org/10.1371/journal.pone.0196075.g002 Phenotypic variability in hydraulic traits PLOS ONE | https://doi.org/10.1371/journal.pone.0196075 May 1, 2018 6 / 17 species, we could not obtain absolute values of K s for the FI-HYY population, due to software recording issues. No differences between sites were observed for K s (Table 1). Populus tremula had the smallest CV inter of the four species studied (14.16%), but the largest CV sp (80.04%, Table 2). Betula pendula had the largest CV inter , with a mean difference of up to 2.8 kg m -1 MPa -1 s -1 between the populations located at the extreme ends of its latitudinal distribution range (Table 2,Fig 2). When significant, the correlations between K s and branch diameter were weak (Betula pendula rho = -0.348, p = 0.008; Pinus sylvestris rho = 0.267, p = 0.005; Populus tremula rho = 0.090, p = 0.542; Picea abies rho = -0.210, p = 0.055). Finally, branch growth differed between populations for all species, whereas site (nested in population) had a significant effect on branch growth only for conifers (Table 1,Fig 2). The phenotypic variability of branch growth was greater than that of P 50 (Table 2,Fig 2). Furthermore, the phenotypic variability of branch growth was greater than that of K s in most cases (Table 2,Fig 2). Pinus sylvestris had the largest CV sp and CV inter in branch growth (85.81 and 76.67%, respectively), whereas these two coefficients were the lowest in Populus tremula (43.29 and 37.66%, respectively) (Table 2). Phenotypic clines with climate and latitudinal gradients Populus tremula presented strong significant clines in P 50 , as five out of the six climatic variables studied here showed significant correlations with P 50 (Table 3,Fig 3). P 50 values for this species were positively correlated with latitude and aridity index (AI and AI_Sum), but negatively correlated with MAT and T_sum (Table 3,Fig 3). There was also a statistically significant negative correlation between P 50 and T_Sum in Betula pendula (Table 3). By contrast, no significant clines in P 50 were observed for conifers (Table 3,Fig 3). K s was less strongly related to climate than P 50 . The K s /climate correlation was statistically significant only for Betula pendula, with lower K s values at sites with higher MAT values (Table 3,Fig 4). We found steeper clines for branch growth than for hydraulic properties, with all species showing at least one statistically significant correlation between branch growth and latitude/ climate variables (Table 3,Fig 5). Latitude was correlated with branch growth only in Populus tremula, for which the lowest branch growth values were obtained for the northernmost population (Table 3,Fig 5). In general, when statistically significant, branch growth was positively correlated with MAT and T_sum, and negatively correlated with AI, AI_sum and MAP (Table 3,Fig 5). Safety-efficiency trade-off At the intraspecific level, we found statistically significant but weak positive correlations between P 50 and K s for the conifers studied, with the most vulnerable individuals having the Table 2. Intraspecific (CV sp ) and inter-population (CV inter ) coefficient of variability (%) for the xylem pressure inducing a 50% loss of conductance (P 50 , MPa), xylem-specific hydraulic conductivity (K s , kg m -1 MPa -1 s -1 ) and branch growth (BG, estimated as branch radius/xylem age (mm/year)) for each study species. CV sp CV inter Species P 50 K s BG P 50 K s BG Betula pendula 9.67 53.47 58.94 5.56 50.28 47.05 Populus tremula 25.07 80.04 43.29 24.82 14.16 37.66 Picea abies 6.57 49.07 57.49 4.15 42.68 53.48 Pinus sylvestris 10.23 48.41 85.81 6.45 23.80 76.67 https://doi.org/10.1371/journal.pone.0196075.t002 Phenotypic variability in hydraulic traits PLOS ONE | https://doi.org/10.1371/journal.pone.0196075 May 1, 2018 7 / 17 largest hydraulic conductivities (Table 4,S3 Fig). No significant correlation was found between P 50 and K s for either of the angiosperms studied (Table 4). Discussion We assessed the phenotypic variability of hydraulic safety and efficiency traits (P 50 and K s , respectively) and branch growth in four tree species across a long latitudinal gradient covering most of their distribution range in Europe. P 50 displayed lower phenotypic variability than K s and branch growth, consistent with our initial hypothesis. The low variability of P 50 across populations has been related to uniform evolutionary selection or canalization [37,38]. Indeed, these studies provided evidence of natural selection acting on this trait. This uniform selection reproduces trait conservatism and eventually leads to stasis [37]. In contrast, it has been suggested that K s variability is related to the interaction between genotype and environment [35]. K s may also vary significantly with sampling position along the branch axis [48,49], although we tried to overcome this limitation through the use of systematic sample preparation procedures. We also expected branch growth to be more variable than P 50 , because branch growth is strongly influenced by multiple interacting factors, such as the availability of nutrients, light, water and temperature [46,50,51,52,53] and biotic interactions [54]. The limited embolism resistance variability observed here, in all species other than Populus tremula, suggests that forest populations of the studied species will potentially find difficulties to cope with a warmer and drier conditions by increasing their embolism resistance. However, considering the differences in climatic ranges between the studied species, we have to be cautious when interpreting these patterns. Further studies investigating larger precipitation gradients and/or marginal Table 3. Correlation coefficients (Pearson or Spearman) and p-values for the relationships between the mean xylem pressure inducing a 50% loss of conductance (P 50 , MPa), xylem-specific hydraulic conductivity (K s , kg m -1 MPa -1 s -1 ) and branch growth (BG, estimated as branch radius/xylem age (mm/year)) and the climatic variables for each sampling site. Betula pendula Populus tremula Picea abies Pinus sylvestris Cor. p Cor. p Cor. p Cor. p P 50 Latitude 0.167 0.692 0.750 0.032 -0.568 0.054 0.193 0.547 MAT -0.547 0.161 -0.770 0.025 0.420 0.174 0.056 0.862 MAP 0.539 0.168 -0.214 0.610 0.288 0.364 -0.466 0.127 AI 0.460 0.251 0.886 0.003 -0.112 0.728 -0.462 0.130 T_Sum -0.793 0.019 -0.909 0.002 0.341 0.278 -0.120 0.711 AI_Sum -0.289 0.487 0.934 0.001 0.098 0.761 -0.578 0.049 K s Latitude 0.228 0.586 -0.166 0.753 0.001 1.000 -0.122 0.704 MAT -0.886 0.003 0.308 0.553 -0.147 0.649 0.408 0.187 MAP -0.119 0.779 0.086 0.872 -0.414 0.181 0.276 0.384 AI 0.231 0.582 -0.206 0.695 -0.239 0.455 -0.279 0.379 T_Sum -0.428 0.290 0.292 0.575 0.082 0.799 0.225 0.481 AI_Sum -0.180 0.670 -0.439 0.383 -0.028 0.931 0.019 0.952 BG Latitude 0.497 0.210 -0.894 0.003 -0.386 0.215 -0.224 0.484 MAT -0.423 0.296 0.876 0.004 0.755 0.004 0.627 0.029 MAP -0.786 0.021 0.452 0.260 0.133 0.679 0.027 0.934 AI -0.669 0.069 -0.904 0.002 -0.779 0.003 -0.593 0.042 T_Sum 0.364 0.375 0.909 0.002 0.870 0.002 0.394 0.205 AI_Sum 0.527 0.179 -0.848 0.008 -0.394 0.205 -0.387 0.214 Statistically significant correlations are highlighted in bold. MAT: mean annual temperature (˚C); MAP: total annual precipitation (mm); AI: aridity index (MAP/PET or potential evapotranspiration). T_Sum (˚C) and AI_Sum: averaged mean temperature and aridity indices, respectively, for June, July and August. https://doi.org/10.1371/journal.pone.0196075.t003 Phenotypic variability in hydraulic traits PLOS ONE | https://doi.org/10.1371/journal.pone.0196075 May 1, 2018 8 / 17 populations are needed. Indeed, a recent study showed that marginal populations of beech significantly differed in embolism resistance [55] while core populations exhibited similar P 50 values [36]. Our results also show that within species phenotypic variability in K s and growth are large, and in general larger than that of P 50 , suggesting that intra-population variability should not be neglected in further studies at local scales. Fig 3. Mean P 50 (MPa) per population plotted against latitude (3.a, decimal degrees) and the climatic variables for each sampled population and site: 3.b. mean annual temperature (MAT, ˚C); 3.c. total annual precipitation; (MAP, mm); 3.d. AI: aridity index (MAP/PET or potential evapotranspiration). https://doi.org/10.1371/journal.pone.0196075.g003 Phenotypic variability in hydraulic traits PLOS ONE | https://doi.org/10.1371/journal.pone.0196075 May 1, 2018 9 / 17 39. Anderegg WRL (2015) Spatial and temporal variation in plant hydraulic traits and their relevance for climate change impacts on vegetation. New Phytol 205(3):1008–1014. PMID: 25729797 40. Hijmans RJ, Cameron SE, Parra JL, Jones PG, Jarvis A (2005) Very high resolution interpolated climate surfaces for global land areas. Int J Clim 25:1965–1978. 41. Zimmermann NE, Roberts DW (2001) Final Report of the MLP climate and biophysical mapping project, Birmensdorf. 42. Torres-Ruiz JM, Cochard H, Choat B, Jansen S, Lo ´pez R, Toma ´s ˇkova ´I, et al. (2017) Xylem resistance to embolism: presenting a simple diagnostic test for the open vessel artefact. New Phytol 215:489– 499. https://doi.org/10.1111/nph.14589 PMID: 28467616 43. Cochard H (2002) A technique for measuring xylem hydraulic conductance under high negative pressures. Plant, Cell and Environment 25:815–819. 44. Cochard H, Damour G, Bodet C, Tharwat I, Poirier M, Ameglio T (2005) Evaluation of a new centrifuge technique for rapid generation of xylem vulnerability curves. Physiol Plantarum 124:410–418. 45. Pammenter NW, Vander Willigen C (1998) A mathematical and statistical analysis of the curves illustrating vulnerability of xylem to cavitation. Tree Physiol 18:589–59. PMID: 12651346 46. Weiskittel AR, Maguire DA, Monserud RA (2007) Response of branch growth and mortality to silvicultural treatments in coastal Douglas-fir plantations: Implications for predicting tree growth. For Ecol Manage 251:182–194. 47. R Development Core Team (2016) R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna. 48. Yang S, Tyree MT (1993) Hydraulic resistance in the shoots of Acer saccharum and its influence on leaf water potential and transpiration. Tree Physiol 12:231–242. PMID: 14969914 49. Petit G, Anfodillo T, Mencuccini M (2008) Tapering of xylem conduits and hydraulic limitations in sycamore (Acer pseudoplatanus) trees. New Phytol 177(3):653–64. https://doi.org/10.1111/j.1469-8137. 2007.02291.x PMID: 18069964 50. Makinen H, Saranpaa P, Linder S (2001) Effect of nutrient optimization on branch characteristics in Picea abies. Scand J Forest Res 16:354–362. 51. Makinen H (2002) Effect of stand density on the branch development of silver birch (Betula pendula Roth.) in central Finland. Trees Struct Funct 16:346–353. 52. Hatfield JL, Prueger JH (2015) Temperature extremes: Effect on plant growth and development. Weather Climate Extremes 10(A):4–10. 53. Lipiec J, Doussan C, Nosalewicz A, Kondracka K (2013) Effect of drought and heat stresses on plant growth and yield: a review. Int Agrophys 27(4):463–477. 54. Makinen H (1999) Effect of stand density on radial growth of branches of Scots pine in southern and central Finland. Can J Forest Res 29:1216–1224. 55. StojnićS, Suchocka M, Benito-Garzo ´n M, Torres-Ruiz JM, Cochard H, Bolte A, et al. (2017) Variation in xylem vulnerability to embolism in European beech from geographically marginal populations. Tree Physiol in press. 56. Mencuccini M, Comstock J (1997) Vulnerability to cavitation in populations of two desert species, Hymenoclea salsola and Ambrosia dumosa, from different climatic regions. J Exp Botany 48:1323– 1334. 57. Vander Willigen C, Pammenter NW (1998) Relationship between growth and xylem hydraulic characteristics of clones of Eucalyptus spp. at contrasting sites. Tree Physiol 18: 595–600. PMID: 12651347 58. Schuldt B, Knutzen F, Delzon S, Jansen S, Mu¨ller-Haubold H, Burlett R, Clough Y, Leuschner C (2016) How adaptable is the hydraulic system of European beech in the face of climate change-related precipitation reduction? New Phytol 210(2):443–458. https://doi.org/10.1111/nph.13798 PMID: 26720626 59. Maherali H, De Lucia EH (2000) Xylem conductivity and vulnerability to cavitation of ponderosa pine growing in contrasting climates. Tree Physiol 20:859–867. PMID: 11303576 60. Martı ´nez-Vilalta J, Piñol J (2002) Drought-induced mortality and hydraulic architecture in pine populations of the NE Iberian Peninsula. For Ecol Manage 161: 247–256. 61. Cornwell WK, Bhaskar R, Sack L, Cordell D, Lunch CK (2007) Adjustment of structure and function of Hawaiian Metrosideros polymorpha at high vs. low precipitation. Funct Ecol 21:1063–1071. 62. Hajek P, Kurjak D, von Wu¨hlisch G, Delzon S, Schuldt B (2016) Intraspecific variation in wood anatomical, hydraulic, and foliar traits in ten European beech provenances differing in growth yield. Front Plant Sci 7:791. https://doi.org/10.3389/fpls.2016.00791 PMID: 27379112 63. Ahmad HB, Lens F, Capdeville G, Burlett R, Lamarque LJ, Delzon S (2017) Intraspecific variation in embolism resistance and stem anatomy across four sunflower (Helianthus annuus L.) accessions. Physiol Plant in press. Phenotypic variability in hydraulic traits PLOS ONE | https://doi.org/10.1371/journal.pone.0196075 May 1, 2018 16 / 17 64. Lexer C, Fay M, Joseph J, Nica MS, Heinze B (2005) Barrier to gene flow between two ecologically divergent Populus species, P. alba (white poplar) and P. tremula (European aspen): the role of ecology and life history in gene introgression. Mol Ecol 14: 1045–1057. https://doi.org/10.1111/j.1365-294X. 2005.02469.x PMID: 15773935 65. Lexer C, Joseph J, van Loo M, Prenner G, Heinze B, Chase MW, et al. (2009) The use of digital imagebased morphometrics to study the phenotypic mosaic in taxa with porous genomes. Taxon 58:349– 364. 66. Delzon S, Sartore M, Burlett R, Dewar R, Loustau D (2004) Hydraulic responses to height growth in maritime pine trees. Plant Cell Environ 27(9):1077–1087. Phenotypic variability in hydraulic traits PLOS ONE | https://doi.org/10.1371/journal.pone.0196075 May 1, 2018 17 / 17