Uniform selection as a primary force reducing population genetic differentiation of cavitation resistance across a species range.

BACKGROUND: Cavitation resistance to water stress-induced embolism determines plant survival during drought. This adaptive trait has been described as highly variable in a wide range of tree species, but little is known about the extent of genetic and phenotypic variability within species. This information is essential to our understanding of the evolutionary forces that have shaped this trait, and for evaluation of its inclusion in breeding programs.
METHODOLOGY: We assessed cavitation resistance (P(50)), growth and carbon isotope composition in six Pinus pinaster populations in a provenance and progeny trial. We estimated the heritability of cavitation resistance and compared the distribution of neutral markers (F(ST)) and quantitative genetic differentiation (Q(ST)), for retrospective identification of the evolutionary forces acting on these traits. RESULTS/DISCUSSION: In contrast to growth and carbon isotope composition, no population differentiation was found for cavitation resistance. Heritability was higher than for the other traits, with a low additive genetic variance (h(2) (ns) = 0.43±0.18, CV(A) = 4.4%). Q(ST) was significantly lower than F(ST), indicating uniform selection for P(50), rather than genetic drift. Putative mechanisms underlying Q(ST)<F(ST) are discussed.

Introduction

The climatic niches of forest tree species are moving faster than the maximum rate of migration, measured by palynology or gene flow analysis, as a direct consequence of the current increase in temperatures due to global warming [1], [2]. Forest tree populations are thus facing new selection pressures and are unable to track their bioclimatic envelope [3] over the time scale at which these changes are occurring. The local adaptation and survival of tree populations in a rapidly changing environment with warmer temperatures and more frequent water shortage is a major concern in efforts to ensure the sustainability of forest ecosystem services. In addition to this trend, climate experts are predicting more extreme climatic events, such as periods of severe drought [4], which will increase mortality rates [5], [6], [7]. These effects on tree mortality highlight the way in which the impact of climate change may depend on the changes associated with extreme events rather than trends [8]. In this context, there is a need to investigate relevant drought tolerance-related traits, to quantify both genetic variation and phenotypic plasticity, which together define the capacity of tree populations to adapt.

From a physiological point of view cavitation resistance is an important trait to estimate drought tolerance. Indeed, dysfunctions of the vascular system of the tree, such as xylem embolism due to cavitation events, is likely to be a key factor governing the mortality of these long-lived organisms [9].When a cavity is formed in the xylem sap under tension (negative pressure), it may spread in the vascular system through intervessel or intertracheid pits, thus compromising the capacity of the plant to transport water [10].

A direct causal link between survival (fitness) and cavitation resistance during extreme drought has been highlighted, based on two lines of evidence suggesting that cavitation resistance is an important adaptive trait. Firstly, assessments of the correlation between cavitation resistance and lethal water potential [11], [12] demonstrated a highly significant linear relationship (r2 = 0.9) between these two traits. Secondly, global surveys of cavitation resistance in woody species showed that xeric species are more resistant to embolism than hydric species [13], [14], [15]. These interspecific studies, with adaptive inferences concerning cavitation resistance being rendered robust by the incorporation of phylogenetic information [15], [16], [17], concluded that cavitation resistance-related traits are under natural selection [15], [18]. To validate this evolutionary pattern, a population-level perspective is appropriate, because variation observed across species cannot be assumed to reflect patterns within species.

At the intraspecific level, cavitation resistance can be analyzed by provenance or progeny trials. The few studies carried out to date (reviewed in Table S1) have included only small numbers of individuals (<9) and populations (<5), and it has therefore not been possible to estimate environmental and genetic effects on phenotypic variation accurately. We therefore still know little about the genetic determinism and micro-evolutionary pattern of this hydraulic trait, but such information is absolutely necessary if this trait is to be incorporated into breeding programs and for a more fundamental understanding of the evolutionary basis of tolerance to severe drought.

The aim of this study was to provide the first estimates of heritability, additive variance and population differentiation for cavitation resistance-related traits. We carried out a case study of maritime pine (Pinus pinaster Ait.), a species with a fragmented distribution in the western part of the Mediterranean region. The scattered distribution of this species may have prevented or limited gene flow between different groups of populations, promoting high levels of genetic divergence between ecotypes due to genetic drift [19] and/or natural selection (Quezel and Barbero 1998 in [20]). Here, we took advantage of a new technology (high-throughput phenotyping platform for cavitation resistance) to screen for the first time a large number of genotypes from six ecotypes to test the hypothesis that Pinus pinaster populations have been subjected to diversifying selection for cavitation resistance. More specifically, this experiment aimed to address the following questions: what is the level intraspecific variation and heritability for cavitation resistance? Can we separate the relative roles of drift and selection in population differentiation for this trait?

Methods

Provenance trial and climatic data

We carried out a provenance-progeny trial, in which young trees (six-year-old plants) were planted in December 2003 at the INRA forestry station in the Aquitaine Region (44°44’N, 00°46’W). The mean annual temperature at this site is 13.2°C and mean annual rainfall is 836 mm (1984–2006). The soil is a sandy podzol with a water table rising to about 0.5 m below the surface in winter and descending to a depth of 2 m in late summer. Seedlings were grown in the nursery from open-pollinated seeds collected from 24 natural populations (or ecotypes) in France, Spain, Morocco and Tunisia, to cover the fragmented distribution of Pinus pinaster. Each population was represented by 20 to 30 half-sib families. The trial was arranged in a randomized block design (15 blocks) with single-tree plots. Each block contained at least one tree from each half-sib family. There were 600 seedlings per block, giving 9,000 seedlings for the entire experiment.

Choice of populations

The assessment of cavitation resistance is a intensively time-consuming process [21]. We therefore designed a procedure for the selection of a subset of populations representing all the variability of climatic envelope of maritime pine. For a total of 754 grid points covering the entire natural range of the species [22], we first extracted climatic data from the CRU CL 2.0 10′ global dataset for the period 1961–1990 [23], [24], [25]. These data included monthly average precipitation, mean, minimum and maximum temperature, diurnal temperature range, water vapor pressure, cloud cover, wet day frequency, ground frost frequency, radiation, wind, Martonne index, Turc's potential evapotranspiration and soil water deficit. We also derived the vapor pressure deficit from these parameters. Principal component analysis (PCA) was then used to reduce the number of dimensions for the whole set of climate variables (Figure 1). The data were centered and scaled before PCA. The 14 populations were finally placed on the main plane of the PCA (accounting for 76% of the variation, Table S2) and six of these populations (Table 1) were selected as a representative subset of the climatic envelope covered by Pinus pinaster species. In each population, eight families (5 half-sibs/family/block) were randomly selected for further analysis.

Figure 1

Principal component analysis (PCA) on the [763 population locations x 14 climatic variables] data matrix.

Upper panel: The contour plot represents the presence's probability (kernel density estimate) of Pinus pinaster population (small black dot) within the bioclimatic envelope representing by PC1 and PC2 axes, accounted for 54% and 21% of the variance, respectively. The studied populations and provenance test are indicated by black circles. PCA was performed with the variables indicated in the methods section. See Table S2, for additional information about the relative contribution of climatic variables to the axes. Lower panel: projection of 14 climatic variables on the subspace spanned by the first two eigenvectors (correlation circle).

Table 1

Climatic data, location and elevation of the studied maritime pine populations.

Sampling location	Longitude(°)	Latitude(°)	n	Elevation(m)	P i(mm)	T m(°C)	VPDmax(hPa)	ETP( mm)
Bayubas de Abajo (Central Spain)	−2.87	41.52	39	955	561	10.5	11.42	882.9
Coca (Central Spain)	−4.08	41.37	40	788	452	11.9	14.23	718.7
Mimizan (South-western France)	−1.30	44.13	40	37	1176	13.2	7.26	751.59
Oria (South-eastern Spain)	−2.62	37.87	40	1232	451	13.4	14.29	922.59
San Cipriano (Northern Spain)	−8.70	42.13	40	310	1625	13.8	8.54	721.91
Tamrabta (Southern Morocco)	−5.02	33.66	40	1760	550	15.1	18.56	976.54

n, number of sampled individuals for hydraulic measurements; P i , mean annual precipitation; T m, mean annual air temperature; VPD, maximal of water vapor pressure deficit (in July for all the provenance); ETP, annual sum of potential evapotranspiration.

Sample preparation for the assessment of cavitation resistance

We collected branches, according to the sampling procedure described below, during winter 2009, before 10 am, to avoid native embolism. Needle water potential was lower than -1 MPa, far from the minimum needle water potential in summer (−2 MPa) of Pinus pinaster [26]. The branch sample corresponded to the 2007 and 2008 growth units on the 2007 whorl when possible, in order to measure on the same number of rings in each sample. Sampled branches were fully exposed to the sun, longer than 40 cm and with a diameter of 0.3 to 1 cm (<4 years of age). The current needles were removed and the branches were wrapped in wet paper towels and bagged upon collection to prevent dehydration. In the lab, samples were cut under water just before measurement to obtain 0.28 m segments (i.e., much longer than the longest tracheid). Bark was removed from all segments before measurements.

Assessment of cavitation resistance

Cavitation resistance was measured on 240 genotypes (6 populations * 8 families * 5 offsprings), with the Cavitron technique [21], [27], [28]. Centrifugal force was used to establish negative pressure in the xylem and to provoke water stress-induced cavitation, using a custom-built honeycomb rotor (Precis 2000, Bordeaux, France) mounted on a high-speed centrifuge (HS18, MSE Scientific, London, UK). This technique enables measurement of the hydraulic conductance of a branch under negative pressure. Xylem pressure (P x) was first set to a reference pressure (−0.5 MPa) and maximal conductance (k max) was determined by measuring the flux of a reference ionic solution (10 mmol dm−3 KCl and dm−3 mmol dm−3 CaCl2 in deionized water) through the sample. The centrifugation speed was then set to a higher value for 3 min to expose the sample at a more negative pressure. Conductance (k) were measured 4 times for each step, and the average was used to compute percent of loss of xylem of conductance (PLC in %) following PLC = 100 (1−k/k max). The procedure was repeated for at least eight pressure step with a −0.5 MPa step increment until PLC reached at least 90%. The percent loss of xylem conductance as a function of xylem pressure (MPa) represents the sample's vulnerability curve (Figure 2). Rotor velocity was monitored with a 10 rpm resolution electronic tachymeter (A2108-LSR232, Compact Inst, Bolton, UK) and xylem pressure was adjusted to about ±0.02 MPa. We used Cavi_soft software (version1.5, University of Bordeaux) for conductance measurements and computation of all vulnerability curves (VCs). The 10,800 measurements of conductance were performed at the new high-throughput phenotyping platform for hydraulic traits (CavitPlace, University of Bordeaux, Talence, France).

Figure 2

Vulnerability curves of one genotype for each studied population.

Black dot are the raw measure of percent of loss of conductance (PLC in %) along the negative pressure gradient (in MPa). The grey line is the Weibull reparameterized model and the black line is the sigmoid reparameterized Model.

Based on a sensitivity analysis (data not shown) and graphical checking (Figure 2), we retained the reparameterized sigmoid function fitted to the conductance data (k) (see [29] and [30] for an exhaustive review) rather than the Weibull model, for determination of the pressure at which the sample lost 50% of its conductance (P) and the slope of the curve at P (S).where k max is the highest conductance measured for each sample (equivalent to k sat in the original Ogle's model), k is the mean conductance at a given pressure concerned, X is the percentage loss of conductance of interest (in %), P (in MPa) is the pressure inducing X% loss of conductance and S (in MPa.% −1) is the slope of the tangent at the P abscissa point on the curve. Analysis has been performed for pressures and slopes corresponding to X = 12, 50 and 88 % loss of conductivity (P 12, P 50, P 88, S 12, S 50 and S 88 respectively).

Carbon isotope ratio and growth measurement

Carbon isotope ratio (δ 13C in %o) was obtained as previously described [32], [33]. Needles of the growth unit used for cavitation analysis were harvested and 20 needles were sampled at random. The needles were dried and ground to a powder and 3 mg sample was analyzed with an isotope ratio mass spectrometer (FISONS Isochrom, Manchester, UK) at INRA facility in Reims (France). Total height was measured at the ages of two (2004) and three (2005) years, on the same six populations and eight families, for all 15 blocks. The annual increase in height (Δh) was calculated as the difference between these two measurements (in 2004 and 2005).

Quantitative genetic analysis

Genetic analysis was conducted with the following mixed model:where is the vector of observation for a trait, is the vector (number of block) of fixed block effects, is the vector (number of populations) of random population effects, is the vector (number of mother trees) of the random genetic effects of mother tree within the population, ε is the vector (number of individuals) of residuals, X is called the design matrix, Z and Z are the incidence matrices linking the observations to the effects. A variance was fitted for each random effect: is the genetic variance between populations, is the genetic variance between mother trees nested within a population and is the residual variance.

Variance or covariance components were estimated by the restricted maximum likelihood (REML) method, assuming a normal distribution of the random effects. The significance of variance components were tested using log-likelihood ratio tests. We included population as a random effect to draw inference at species levels [34] and to obtain an unbiased estimate of heritability and genetic population differentiation [35]. The normality, identity and independency of residuals of each trait were graphically checked by plotting studentized marginal and conditional residuals (available on request), which confirmed that the data match with the assumption of mixed model. We estimated narrow-sense heritability as follows:where is the within-population additive variance. In our study, was estimated by  =  as trees from the same family were presumed to be half-sibs (open-pollinated seeds). We did not include the population effect in the heritability calculation, because natural selection appeared to occur within each population [36]. The standard deviation of heritability was calculated with the equations of delta method (see Appendix in [37]).

Variance components were standardized by the trait mean [38] as follows, CVX = 100√(Variance)/MeanX where X is the trait considered, and CV is the coefficient of variation. Each variance component is expressed with a CV (CVA: additive coefficient of variation; CVBP (VBP = ): coefficient of variation between populations; CVP: phenotypic coefficient of variation; CVR: residual coefficient of variation). The variance of each component was extracted from the asymptotic covariance matrix. The significance of mean population difference was estimated using the same model (Eq. 2) with a proc GLM with a Student-Neuman-Keuls post hoc test.

Correlation between traits

To facilitate interpretation of correlation, negative value of P 50 and δ 13C were converted from negatives to positives. Genetic correlations between traits were evaluated by calculating Pearson's coefficient on the family Best Linear Unbiased Predictor estimation (BLUP, for additional information see [37] p745). BLUP estimation ensures that data are corrected for block effect. We will refer to these correlations as genetic correlations. For phenotypic correlation, all Pearson correlations were computed over the BLUP family plus BLUP population and the grand-mean.

Estimation of population differentiation

The estimate of phenotypic differentiation between populations, Q [39], was calculated asPutatively neutral nuclear microsatellites (nuSSRs) were used to account for genetic differentiation (F) caused by demographic and other processes not related to selection (e.g., genetic drift resulting from geographic isolation or population expansion). Eight polymorphic nuSSRs were selected from those previously developed by [40] (NZPR413, NZPR1078, ctg64), [41], [42] (ctg275, FRPP91, FRPP94, ITPH4516) and [43] (A6F03). The markers were selected to be evenly distributed over the various linkage groups of the maritime pine genetic linkage map, with, at least, 4 alleles and multiplexing capacity. Genotyping was performed on genomic DNA isolated from the needles of 20 to 30 individuals from each of the six selected populations, as previously described [44], [45]. We used F ST (which is estimated from the allelic frequency) rather than R ST (which also takes into account allele size) because genetic drift affects allele frequency but not mutation rate. F values for each locus were estimated with Genepop [46] using the framework developed by [47] adapted for SSR data [48].

F ST and Q ST comparison

For the comparison of Q and F, to disentangle the effects of genetic drift from those of selection, we develop a new test to avoid previously reported limitations [49], [50], [51], [52] and allow to test Q >F and Q 51], [53]. We will refer to this distribution of F ST * as the “drift distribution”.

We estimated the sampling variance of Q ST (Eq 4), by simulating the distribution of each variance component (,) with a parametric bootstrap [50], using the Satterthwaite approximation [54]. This distribution is highly conservative and takes into account the deviation from homogeneity of variance [54]. is a parametric bootstrap replicate of variance component of i factor. It is obtained by multiplying (observed variance component) by a chi-square distribution scaled with an “effective” degree of freedom (dfe). dfe i is the effective degree of freedom for the variance component of i factor. df i is the observed degree of freedom. is observed variance component due to i factor. is observed variance component due to i-1 factor (nested or residuals factor). n is the total size of the sample. j is the number of level of factor i. We calculate a Q ST * for each replicate from and (Eq. 6). We will refer to this distribution of Q ST * as the “phenotypic distribution”, although Q ST is a standardized measurement of additive genetic variance between populations.

Finally, we compared the F ST * and Q distributions, using nonparametric and free distribution two-sample test for equality of the 2.5 and the 97.5 quantile (see [55] for theoretical proof) with a Bonferroni correction for multiple comparisons. We also applied a studentized bootstrap, which gave similar results but required more computation time [56], [57], [58]. All the analyses were performed with SAS version 9.2. Codes are available on request.

Results

1. Between-population variation

For each population, vulnerability curves showed similar sigmoid shape with the air-entry (P 12) around −3.25±0.006 MPa (see Figure 2). Linear curves were discarded from the analysis [59]. The between-population effect (VBP) was significant for δ13C and Δh but not for P 50 (Table 2). Similarly, no difference was found for the other cavitation resistance-related traits (S 12, S 88, S 50, P 12, P 88, data not shown). Cavitation resistance-related traits had much lower coefficients of variation than Δh (Table 2). This was particularly true for the between-population coefficient of variation (CVBP = 1% and 18% for P 50 and Δh, respectively). It should be noted that CVs for δ 13C are not comparable with those of other traits, because they are estimated relative to a standard [60] and are therefore independent of scale change but not of origin. The fixed block effect was significant for all the traits studied, indicating that some of the environmental variation was taken into account by the experimental design.

Table 2

Variance components (VP, VBP, VA, VR), narrow-sense heritability (h2 ns), coefficient of variation (CVP, CVA, CVBP, CVR) and population differentiation (Q ST) for all studied maritime pine populations.

Traits	VP	VBP	VA	VR	h2 ns±SE	CVP	CVA	CVBP	CVR	QST
P 50	0.067	0.002 ns	0.028*	0.058	0.438±0.18	6.6	4.4	1	6.2	0.027
δ13C	0.284	0.030 **	0.059*	0.269	0.213±0.10	1.7a	0.8a	0.6a	1.7a	0.197
Δh	112.7	55.0 ***	40.96***	102.5	0.363±0.06	26.9	16.2	18.8	25.7	0.188

h is the narrow-sense heritability and SE is the standard error of heritability, V is the phenotypic genetic variance, V is the additive genetic variance, V is the between-population variance, V is the residual variance. CV is the variation coefficient of additive variance after adjustment for the block effect. CV is the variation coefficient of phenotypic variance after adjustment for the block effect. CV is the residual coefficient of variation. CV is between-population coefficients of variation. Q is the genetic quantitative variation between populations (Spitze, 1993). The significance of random effects is indicated after each variance estimator: ns P > 0.05, * P<0.05, ** P<0.01, *** P<0.001. a CVs for δ13C are not comparable with other traits as they are estimated relative to a standard.

The populations from the wettest areas (Mimizan and San Cipriano) had the highest Δh values (Figure 3a), whereas Tamrabta population (from Morocco) presented the lowest value. Iberian populations from very different climatic areas (Coca, Bayubas, Oria) had intermediate values, with no detectable trend as a function of environmental aridity.

Figure 3

Mean values of height increment (Δh, (a)) (n = 297 per population).

Mean values of cavitation resistance (P 50, (b)) and carbon isotope composition (δ13C, (c)) for each studied population (n = 40 per population). The error bars represents the standard errors. Different letters indicate significant differences between populations at α = 0.05.

No significant difference between populations was detected for P 50 (Figure 3b), although Tamrabta surpassed the other populations and was the most cavitation-resistant population. Tamrabta also presented the lowest δ 13C value (Figure 3c), demonstrating a significantly lower water-use efficiency than the other populations, all other populations presented similar δ 13C values.

2. Within-population variation

Heritabilities and normalized measurements of trait dispersion (i.e. CVs) were estimated to evaluate the within-population additive variance, evolvability (through the analysis of CVA) and micro-environmental sensitivity (through the analysis of CVR) (Table 2). Narrow-sense heritability (h2 ns) for P 50 was higher (0.44±0.18) than those estimated for δ 13C (0.21±0.10) and Δh (0.35±0.06), showing that cavitation resistance was genetically controlled, although the standard error was high, probably due to the small number of progenies per mother tree analyzed. The CVs of P 50 and Δh presented contrasting patterns, with a lower coefficient of additive variation for P 50 (CVA = 4.4%) than for Δh (CVA = 16.2%), suggesting limited evolvability of P 50.

3. Evolutionary forces driving population differentiation

Q ST and F ST comparisons have three possible outcomes [39]: (i) if Q ST > F ST, the degree of differentiation for quantitative traits exceeds that attainable by genetic drift alone (ii) if Q ST and F ST are not significantly different, the observed degree of differentiation for quantitative traits could have been reached by genetic drift alone, and (iii) if Q ST50], [61] , we found that F ST * and Q ST * presented skewed distributions (Figure 4). Only Δh and P 50 had a Q ST * distribution different from the F ST * distribution (P = 0.003 and P <0.0001 respectively). For, δ 13C, the difference between Q ST * and F ST * values was centered on 0 (see Figure 4c right panel), and it was therefore not possible to distinguish between drift and selection (P = 0.88). Conversely, differences between the Q ST * and F ST * distributions for P 50 were centered on -0.18, suggesting that the studied populations were less differentiated than would be expected in the presence of drift alone (Figure 4b), which means that natural selection favored the same mean phenotype in different populations (consequence of uniform selection) . For Δh, the difference between Q ST * and F ST * distributions was centered around 0.27, suggesting that the studied populations displayed more differentiation than would be expected with drift alone (Figure 4c) which is interpreted as a consequence of diversifying selection.

Figure 4

Comparison between F ST (histogram in gray) and Q ST (histogram in black) distributions for growth rate (Δh, (a)), cavitation resistance (P 50, (b)) and carbon isotope composition (δ 13C, (c)) in the left panel.

The observed distribution (gray histogram) and the kernel density (black curves) of the Q ST-F ST difference are represented in the right panel for each trait. On the right panel, we also show the integral probability of the distribution (using the kernel density estimator) above (see “Integral p(x>0)” on the right panel) and below (see “Integral p(x<0)” on the right panel) zero (marked with the tick and dotted line).

4. Correlation between traits

We found a significant positive phenotypic correlation between absolute value of P 50 and δ 13C at the phenotypic level (r = 0.30, P = 0.035 see Figure S1), indicating that the more cavitation-resistant genotypes tended to be less water-use efficient. However, this relationship was not significant at the genetic level (r = 0.14, P = 0.320). No relationship between P 50 and Δh was found at either phenotypic or genetic level. A significant negative phenotypic correlation between Δh and δ 13C was detected (r = −0.68, P<0.0001). This correlation was barely significant at the genetic level (r = −0.29, P = 0.053).

Discussion

We reliably estimated for the first time the genetic variability of cavitation resistance, a functional trait that allows plants to survive under severe drought. We also provided evidence of natural selection acting on this trait. These results were based on the greatest number of genotypes ever measured to date in an experimental design (Table S1 and supplementary references). Despite the high level of variation of cavitation resistance between species [15], we detected no significant differences between maritime pine populations from a wide range of environments. Moreover, the between-population variability of cavitation resistance was significantly lower than would be expected under a hypothesis of genetic drift alone (Q ST- F ST comparison). We can therefore reject the hypothesis of diversifying selection. We suggest instead that uniform selection has shaped the phenotypic variability of this trait. Uniform selection could be seen as a stabilizing selection acting within each population with the same selection optimum in each population despite the steep climatic gradient [62], [63]. Conversely, growth and water-use efficiency displayed different patterns and were found to be subject to strong diversifying selection and genetic drift, respectively. Quantitative genetics analysis also showed that cavitation resistance presented a significant heritability, higher than that estimated for growth and water-use efficiency. This is the first evidence of uniform selection in woody plants and the underlying mechanisms are discussed below from a micro-evolutionary point of view.

Intra- vs. interspecific variability of cavitation resistance

Despite the steepest climatic gradient (precipitation ranging from 400 to 1,600 mm in the sampled populations) and strong phylogeographic structure between the six studied populations [22], [64]. Very low between-population variance for cavitation resistance were found (CVBP = 1%). The few studies published to date (reviewed in Table S1) tended to skim over the issue of intraspecific variation of cavitation resistance in provenance or progeny trials and reported little or no difference between populations [65], [66]. As these studies were not designed to assess the genetic component of phenotypic variation, further investigations are required, to generalize our finding to other species. In addition, phenotypic variation for cavitation resistance was low (CVP = 6.6%), but consistent with the range reported for wood properties of maritime pine, such as mean ring density, lignin content and fiber morphology ([67], [68] Lamy, unpunblished data). However, we are lacking information about the intraspecific variation of pit pair anatomical traits that are known to be implicated in cavitation resistance [69], [16]. The low within-species variability for cavitation resistance is remarkable (P 50 = −3.93±0.04 MPa, estimated over the whole dataset), given that substantial variability has been described between species. For instance, Delzon et al (2010) showed that cavitation resistance ranged from −3 to −12 MPa in a sample of 40 coniferous species. This variability was interpreted as the effect of natural selection rather than phylogenetic legacy [15].

In contrast, the population differentiation observed for growth and water-use efficiency (WUE) was significant and consistent with previous results for this species [70]. The Moroccan population had the lowest WUE, consistent with previous findings based on both gas exchange measurements and carbon isotope discrimination [71], [72], [73], [74]. In a provenance trial carried out in south-western France, this Moroccan population displayed lower stomatal sensitivity to water stress (delayed stomatal closure), leading to greater water loss throughout the summer period and a lower WUE (as reflected by carbon isotope composition). Genotype × environment interaction could potentially alter differences between populations [75], [76]. Our results therefore require confirmation in provenance trials carried out in drier climates.

Relationships between traits

The weak but significant positive correlation found between absolute value of cavitation resistance and water-use efficiency (carbon isotope composition) suggested that drought-tolerant genotypes had lower water-use efficiency. In dry environments, genotypes that allocate more carbon to the construction of cavitation-resistant wood in order to avoid runaway embolism might be able to maintain higher stomatal conductance and hydraulic conductance at low leaf water potential, resulting in a decrease in water-use efficiency. Our results are consistent with previous findings [77] of a strong and positive relationship between these two traits in two cedar species. However, little or no correlation has generally been reported [78], [79], [80]. The negative correlation between carbon isotope composition and growth has been reported in previous studies [77], [81], [82], assuming that growth is a function of carbon assimilation and carbon isotope composition is an index of retrospective gas exchanges.

Evidence of uniform selection for cavitation resistance

The phenotypic distribution of cavitation resistance was significantly lower than the expected distribution under the drift hypothesis. This may be interpreted as a consequence of uniform selection (also called homogenous, spatially homogenizing, convergent selection, uniform stabilizing selection or stabilizing selection across population). This inference (Q ST52], leading to a false positive result. However, the F ST estimate was more than five times greater that the Q ST estimate. We thus believe that this difference is biologically meaningful and not due to a statistical artifact [52]. The robustness of this result is, also, supported by fourth lines of evidence: (i) the different patterns obtained for growth and WUE in the same experimental design. (ii) Selection procedure of population (see methods) increase the probability to find diversifying selection because we selected extreme populations in term of climatic origin and evolutionary history (different mitotypes and chlorotypes, see [64]), consequently uniform selection could not be interpreted as a sampling bias. (iii) Wood density (measured by X ray on the same data set) showed exactly the same pattern (Lamy, unpunblished). (iv) Willson et al (2008) showed from interspecific data with narrow taxon sampling (limited to the Juniperus genus), that cavitation resistance gave strong phylogenetic conservatism, suggestive of uniform selection for the maintenance of ancestral traits.

For growth, diversifying selection was highlighted by Q ST being greater than F ST [83]. For WUE, we detected no signature of selection, as the phenotypic and drift distributions did not differ significantly, as previous reported by [84] for Quercus suber (P ST/F ST comparison). However, the distribution (Q ST-F ST) of this trait was shifted to the right (integral probability above 0 > integral probability below 0, Figure 4c right panel), suggesting that diversifying selection was slightly more pronounced than drift.

What are the mechanisms behind “uniform selection” for cavitation resistance?

The causal mechanisms underlying uniform selection or leading to evolutionary stasis are not well understood [63], [85], [86], [87], [88]. We discuss here only the processes most likely to account for the observed pattern (Q ST

Weak molecular variation

Lethal or sublethal mutations may limit the variability of the genes they affect, thereby controlling trait variation. In conifers, cavitation resistance is known to be determined by xylem anatomy, including, in particular, the characteristics of intertracheid pits [16], [69]. As knowledge about the nucleotide diversity of different functional categories of genes accumulates, it may become possible to test the hypothesis that genes involved in intertracheid pit formation (once these genes have been identified) display lower levels of diversity.

Genetic constraints

If selection acts on a trait that is negatively correlated with another trait (or traits) also under selection, than the decrease of rate of evolution for the first trait is proportional to the strength of the correlation [94]. A multi-trait approach could be used to explore this hypothesis indirectly [95], but could fail if the trait is canalized.

Canalized trait

This hypothesis suggest that cavitation resistance is canalized to buffer the variation of this key hydraulic trait against all kinds of disturbance, being of genetic (mutation, hybridization, recombination) and/or environmental nature [89], [90], [91]. Emergent properties of molecular networks could buffer molecular variability [92], to maintain phenotypic function, in accordance with the robustness theory [93]. In zoology, dipterian wings shape or centroid size (or mammalian body temperature) are the best known cases of canalized traits [90]. Indeed they reported a similar wing shape between species despite a great variability of climatic niche and a low additive genetic variance between populations for this trait. Except for leaf shape in Arabidopsis thaliana, there is no evidence of canalized trait in plants nowadays. For cavitation resistance, two arguments lead us to consider canalization as the most likely mechanism: (i) low additive genetic variance between populations (V BP) and (ii) the similarity of cavitation resistance values among all the Pinus species [69].

Variance component analysis

Variance component analysis for cavitation resistance resulted in the first estimates of the heritability of this trait (h2 ns = 0.4) and its additive coefficient of variation (CVA = 4.4%). These values suggest that this trait may respond to truncation selection frequently practiced in breeding. However, for a given selection intensity, genetic gain for cavitation resistance would be limited by the low additive variance, although long-term artificial selection experiments have shown that quantitative traits have a non negligible mutational variance [91], [96], [97], which could supply further additive variance at each generation. However, due to the small number of half-sib families and progenies within each family, which could inflate the value for heritability [37], this estimate should be interpreted with caution. Additional studies with a larger sample size are required for further exploration of the genetic determinism of this hydraulically important trait.

The much higher CVA value (16.2%) for height increment is consistent with previous reports and accounts for the genetic gain achieved for this trait over successive generations in breeding programs [98], [99], [100]. For δ13C, a previous study [31] reported a slightly lower heritability (h2 ns = 0.17 vs. 0.21 here), but with estimation based on a diallel cross of limited size, with a narrow genetic background restricted to 12 elite trees from south-western France.

Future directions

Phenotypic variation is a fundamental prerequisite for evolution because natural selection acts on phenotype. Adaptation and evolution via natural selection requires the presence of genetic variation among individuals in a population upon which natural selection can act. Intra-population genetic variability can thus be seen as the fuel for future adaptation. However, an environmentally induced shift in phenotype is also a major component of the variation we see in nature. Recent studies [80], [101] showed a weak but significant phenotypic variability for cavitation resistance. Our results suggest that this between populations variability might be under environmental control rather than genetic determinism. These considerations call for more research (ongoing) aiming at quantifying the in situ phenotypic variability of cavitation resistance and the extent of phenotypic plasticity using provenance trials installed under different edapho-climatic environments. Further studies are also being pursued to dissect the genetic architecture of cavitation resistance to determine the number, map location and effects of Quantitative Trait Loci controlling part of the variation of this trait.

Genetic (right panel) and phenotypic (left panel) correlation between traits. For ease of interpretation, we have converted all the negative values to positive values (P50, δC13). For the genetic correlation, all Pearson correlations (r) were computed over the best linear unbiased prediction (BLUP) (n = 48 for P50, δ13 and n = 151 for Δh). For phenotypic correlation, all Pearson correlations were computed over the BLUP family plus BLUP population and the grand-mean, to ensure that the order of degree of freedom remained the same and the block effects are removed. P50, pressure at 50 % loss of conductivity in MPa, Δh the annual increment between 2004 and 2005, in mm, δ13C is the isotope discrimination for carbon 13 in ‰.

(TIF)

Click here for additional data file.

Review of intraspecific studies for cavitation resistance estimated using Npop: number of populations used, Nind: number of individuals per population used to assess cavitation resistance. The table is divided in two parts, the first part corresponds to provenance or progeny trials, and the second to “in situ” studies.

(DOC)

Click here for additional data file.

Result of the principal component (PC) analysis (PCA) for climatic data of populations (listed in the Contributions to the first, second, third and fourth axes are indicated for each variable (PC1, PC2, PC3, PC4). The eigenvalues of PC1 = 7.65, PC2 = 3.059, PC3 = 0.97, PC4 = 0.86. W is mean wet ground days (days). I is mean Martonne's index (Pi/(Ta+10)). P is the mean precipitation (mm.days−1). C is percent of cloud cover (%). S is the mean of wind speed (m.s−1). V is the water vapor pressure in air (hPa). VPD is the water vapor pressure deficit of air (hPa). T min is the minimum temperature (°C). ΔDT is the mean diurnal temperature range (°C). T m is mean temperature (°C). R is mean global radiation (W.m2). T max is the maximum temperature (°C). H is mean soil water deficit (Pi-ETP, in mm). ETP is mean Truc's potential evapotranspiration (mm).

(DOC)

Click here for additional data file.