Molecular Phylogenetics and Evolution 68 (2013) 555–566

Contents lists available at SciVerse ScienceDirect

Molecular Phylogenetics and Evolution journal homepage: www.elsevier.com/locate/ympev

Adaptive evolution of Mediterranean pines Delphine Grivet a,⇑, José Climent a, Mario Zabal-Aguirre a, David B. Neale b, Giovanni G. Vendramin c, Santiago C. González-Martínez a a

Department of Forest Ecology and Genetics, Forest Research Centre (CIFOR), INIA, Madrid, Spain Department of Plant Sciences, University of California at Davis, Davis, CA, USA c Plant Genetics Institute, Division of Florence, National Research Council, Sesto Fiorentino (FI), Italy b

a r t i c l e

i n f o

Article history: Received 21 August 2012 Revised 10 December 2012 Accepted 31 March 2013 Available online 13 April 2013 Keywords: Phylogeny Adaptation Candidate genes Coevolution Phenotypic traits Pinus

a b s t r a c t Mediterranean pines represent an extremely heterogeneous assembly. Although they have evolved under similar environmental conditions, they diversified long ago, ca. 10 Mya, and present distinct biogeographic and demographic histories. Therefore, it is of special interest to understand whether and to what extent they have developed specific strategies of adaptive evolution through time and space. To explore evolutionary patterns, the Mediterranean pines’ phylogeny was first reconstructed analyzing a new set of 21 low-copy nuclear genes with multilocus Bayesian tree reconstruction methods. Secondly, a phylogenetic approach was used to search for footprints of natural selection and to examine the evolution of multiple phenotypic traits. We identified two genes (involved in pines’ defense and stress responses) that have likely played a role in the adaptation of Mediterranean pines to their environment. Moreover, few life-history traits showed historical or evolutionary adaptive convergence in Mediterranean lineages, while patterns of character evolution revealed various evolutionary trade-offs linking growth-development, reproduction and fire-related traits. Assessing the evolutionary path of important life-history traits, as well as the genomic basis of adaptive variation is central to understanding the past evolutionary success of Mediterranean pines and their future response to environmental changes. Ó 2013 Elsevier Inc. All rights reserved.

1. Introduction Mediterranean pines (section Pinus, subsection Pinaster) represent an extremely heterogeneous and interesting assembly, which ranges from shore and island pines to mountain pines in the areas surrounding the Mediterranean (Klaus, 1989; Richardson, 1998; Fig. 1). These pines are adapted to extremely variable, and frequently extreme, environments. Very often, in Mediterranean forests, different pine species co-occur and thus, they may have adapted to similar environmental conditions. As these taxa are also phylogenetically close, we could expect them to share common adaptive strategies. However, since Mediterranean pines diversified c. 10 Mya (divergence time between the subsections Pinaster and Pinus; Gernandt et al., 2008), each taxon has distinct biogeographic and demographic histories (e.g. Grivet et al., 2011; Fady, 2012), and may thus have responded quite differently to selective pressures. Therefore, we could also expect each individual taxon to

⇑ Corresponding author. Address: Department of Forest Ecology and Genetics, Forest Research Centre (CIFOR), National Institute for Agriculture and Food Research and Technology (INIA), 28040 Madrid, Spain. Fax: +34 91 347 6767. E-mail addresses: [email protected] (D. Grivet), [email protected] (J. Climent), [email protected] (M. Zabal-Aguirre), [email protected] (D.B. Neale), giovanni. [email protected] (G.G. Vendramin), [email protected] (S.C. González-Martínez). 1055-7903/$ - see front matter Ó 2013 Elsevier Inc. All rights reserved. http://dx.doi.org/10.1016/j.ympev.2013.03.032

have evolved following its own trajectory. It is thus of special interest to understand whether and how Mediterranean pines have developed specific strategies of adaptive evolution through time and space, and which genes and phenotypic traits are involved in this process. To explore evolution patterns in the Mediterranean pines, two (independent) sets of common characters were examined, one at the molecular level (candidate genes) and the other at the phenotypic level (life-history traits). Although some of these genes may be directly involved in a given phenotype under study (for example, 4-coumarate-CoA ligase (4cl) with wood density; see Fig. 2 in González-Martínez et al., 2007) currently there is not enough evidence to establish a clear association between genes and phenotypes. Moreover, allelic effects on phenotypes would most probably be polygenic, and strong association between genes and phenotypes is not expected (Rochman, 2012). Genes involved in species adaptation are still largely unknown, but examining genes that are differentially expressed under biotic (e.g. pathogen attack) and abiotic stress (such as high/low temperatures, and deficit/surplus of water, minerals or micronutrients), may help determine those that play a role in adaptive evolution. Examining the same set of genes across the Mediterranean pine taxa should allow discovering not only nucleotide similarities or the genetic variation that underlies differences among pine species, but also identifying

556

D. Grivet et al. / Molecular Phylogenetics and Evolution 68 (2013) 555–566

Fig. 1. Phylogenic trees inferred using 21 low-copy nuclear genes for the seven Mediterranean pine species and the American outgroup. Numbers in bold indicate the probability of the partition for each branch; decimals indicate the sample-wide mean concordance factor (with 95% credibility); numbers in parenthesis correspond to the standard deviation of sample-wide mean concordance factor across runs; circled digits indicate nodes’ identity. The biogeographic region of origin (after Richardson, 1998) is also given.

selection events that have affected one particular branch of the gene tree (Yang, 1998), or/and particular sites under selection (Nielsen and Yang, 1998). The phylogenetic approach employed to detect footprints of selection complements previous population-genetic studies performed at the within species level for different Mediterranean pines (Eveno et al., 2008 and Pot et al., 2005 for P. pinaster; Grivet et al., 2009 and Grivet et al., 2011 for P. pinaster and P. halepensis). This approach examines adaptive evolution of candidate genes at a wider evolutionary scale across several taxa within the subsection Pinaster, and sheds some light on the features that make species similar to each other or unique. Life-history traits have a direct bearing on fitness and thus direct consequences on species adaptive evolution (Stearns, 1977). Because pines of the Mediterranean Basin grow in rather extreme environments, they have developed specific strategies to face them, especially against drought and fire, the main disturbance factors in this region (Tapias et al., 2004, and references therein). Adaptation of Mediterranean conifers to their environment involves a suite of phenotypic traits that have been postulated as central for understanding their life-history strategies. These traits are related to vegetative growth (Grotkopp et al., 2002), seedling allometry and development (e.g. Strauss and Ledig, 1985), sexual

reproduction (e.g. Tapias et al., 2001; Ne’eman et al., 2011), seed dispersal ability (e.g. Benkman, 1995) and, particularly, tolerance to fire (Keeley and Zedler, 1998; Schwilk and Ackerly, 2001; Pausas et al., 2004; He et al., 2012) and drought (e.g. Martínez-Vilalta et al., 2004). Wood physical properties have also received much attention recently, since a close link with contrasted life histories has been shown in different tree species (e.g. Poorter et al., 2010). Finally, genome size has been suggested to play a pivotal role in determining life-history traits in many organisms (e.g. Grotkopp et al., 2004). Our knowledge about the correlated evolution of these traits and their underlying molecular bases are still fairly limited. Furthermore, the lack of a consistent phylogeny for Mediterranean pines has prevented unveiling the evolutionary relevance of phenotypic traits. Reconstructing underlying evolutionary processes requires robust inference of phylogenetic relationships. In the last decades, numerous studies have attempted the construction of pine phylogenies (including all or part of the Mediterranean pines), using chloroplast markers (e.g. Wang et al., 1999; Eckert and Hall, 2006; Parks et al., 2009), nuclear markers (e.g. Liston et al., 1999; Syring et al., 2005), a combination of molecular and non-molecular data (e.g. Gernandt et al., 2008), or combining phylogenetic trees

Fig. 2. Non-synonymous (dN) and synonymous (dS) values along with their ratio, x = dN/dS for locus lp31, plotted for each branch of the consensus phylogenetic tree of the Mediterranean pines’ clade. In bold: branches for which x > 1; div 0: the ratio could not be computed because the denominator equals zero.

D. Grivet et al. / Molecular Phylogenetics and Evolution 68 (2013) 555–566

(supertrees) based on molecular and morphological data (e.g. Schwilk and Ackerly, 2001; Grotkopp et al., 2004). Although Mediterranean pines have been the object of study with regard to the number of taxa and their taxonomic organization, the relationships among Mediterranean pines are still being debated (Willyard et al., 2009). Mediterranean pine species are placed in the subgenus Pinus (diploxylon or hard pines), section Pinus, subsection Pinaster (Gernandt et al., 2005). Previous attempts to reconstruct a well-supported phylogeny for the seven species belonging to the subsection Pinaster – namely, P. brutia, P. canariensis, P. halepensis, P. heldreichii, P. pinaster, P. pinea and P. roxburghii– have been unsuccessful. Because multiple independent markers have the power to disentangle the phylogenetic relationships at low taxonomic levels (see Syring et al., 2005), a new set of 21 unlinked low-copy nuclear genes was used that, when analyzed with multilocus Bayesian tree reconstruction (Ané et al., 2007; Liu and Pearl, 2007), allowed robust phylogenetic inference in Mediterranean pines. Having a well-supported phylogeny in hand opens the possibility of examining adaptive evolution at the interspecific level and over a large evolutionary scale, at both the molecular (see De Mita et al., 2007 and Palmé et al., 2009 for examples in Medicago and conifers, respectively) and the phenotypic (see Guzmán et al., 2009 and He et al., 2012 for examples in Erica and conifers, respectively) levels. The objectives of the present study are: (1) to first propose a robust phylogeny of the subsection Pinaster, and (2) to search for footprints of natural selection in candidate genes for adaptive traits (as shown by expression profiles), as well as at the evolution of phenotypic traits across the phylogeny to gain insight into Mediterranean pines adaptive evolution.

557

of gene expression. Loci were re-sequenced in all seven species of Mediterranean pines using the same primers and produced clear, single bands. Together with the high gene and amino acid sequence homology, this suggests low-copy genes and orthology across species. All DNA sequences were obtained by direct sequencing from haploid seed megagametophytes. Sequencing of haploid tissue minimizes the co-amplification of paralogs, a common problem in plant species with large genomes.

2.2. Phenotypic traits For each species, 12 phenotypic traits related to growth and development (ontogenetic heteroblasty – i.e. juvenile’s vegetative development–age at maturity, maximum size, maximum lifespan, and wood density), reproduction (average seed mass and cone size, seed dispersal ability by wind, and size at maturity) and fire tolerance (bark thickness, cone serotiny, and sprouting ability) (see details in Table 1) were collected. Phenotypic data for all traits, as well as genome size, came from a broad bibliographic survey combined with our own field data, in cases of ambiguous citations or missing data (see details in Table S2). One exception is the seedling ontogenetic heteroblasty index that was evaluated in our nursery common garden, where plants were cultivated in 7-L pots under standard cold greenhouse conditions for 18 months and then harvested for separate dry mass measures (for details about the methods see Climent et al., 2011). Because the present study specifically focuses on Mediterranean pines, traits’ values and evolution patterns are relative to this group of taxa. To gain insight on the adaptive evolution of these traits at a broader evolutionary scale we compared our results with published works that examined the same traits across an extended (but not well-supported for Mediterranean pines) Pinus (subgenus) phylogeny (see Section 4.3.).

2. Materials and methods 2.1. Species and sequence production

2.3. Low-copy genes and sequence editing

We obtained DNA sequence and phenotypic data from the seven pine species forming section Pinus, subsection Pinaster: the Turkish or Calabrian pine (P. brutia), the Canary Island pine (P. canariensis), the Aleppo pine (P. halepensis), the Bosnian pine (P. heldreichii), the maritime or cluster pine (P. pinaster), the stone pine (P. pinea) and the Chir pine (P. roxburghii). Data were also collected for one outgroup from section Trifoliae, subsection Australes, the loblolly pine (P. taeda) (Fig. 1). Mediterranean pines’ biogeographic distribution encompasses the Mediterranean Basin, the Canary Islands and the Himalayas (the Himalayan pine P. roxburghii is included in subsection Pinaster as previous research showed that this taxon is closely related to the Mediterranean pines; see Section 4.1.), while the native distribution of P. taeda is restricted to the southern United States, and this taxon was selected as outgroup because it is the reference species from which most candidate genes were originally isolated. DNA sequences used to infer the pine phylogeny were obtained from various sources and represent random samples covering species’ ranges (see details in Table S1). Most sequences were produced either in our laboratory following protocols from Grivet et al., 2011 (all P. roxburghii sequences and 4 out of 21 loci for all species) or using services at Agencourt Biosciences, MA, USA (the other 17 loci except for P. roxburghii). DNA sequences from P. taeda were retrieved from TreeGenes (http://dendrome.ucdavis.edu/ treegenes/; see Wegrzyn et al., 2008) at University of California, Davis. Loci have been selected based on their expression in wood forming, drought and salt stressed, and/or pathogen-challenged tissues in pine, as well as on knowledge of the genes’ role in xylem formation, disease resistance, abiotic stress response or regulation

Twenty-one low-copy nuclear loci were examined in the seven Mediterranean and the outgroup pine taxa. Four loci studied in our previous work in P. pinaster and P. halepensis (Grivet et al., 2009 and Grivet et al., 2011) were also sequenced in the other species included in this study: the 4-coumarate: CoA ligase gene (4cl-Pt), the dehydrin 2 (dhn2-Pp), and two genes from the ABA stress and ripening family (ASR), lp31-Pt and lp33-Pp. The 17 other loci were specifically selected for this study following a 2-step process: primers were originally designed in P. taeda, and tested in P. pinaster. Then, the successful loci were sequenced in all Mediterranean pine taxa (Table 2). Finally, genes were integrated in the phylogenetic analysis progressively until reaching a well-supported phylogeny (see Section 2.4.). Raw sequences were processed with the Sequal version 2.2 bioinformatic pipeline (Lang T. and Garnier-Géré P., personal communication; program available upon request to PGG), and the final sequence was obtained manually using CodonCode Aligner version 3.7.1 (CodonCode, Dedham, MA, USA). For phylogeny reconstruction, when sequences from multiple individuals were available, each species was represented by the consensus sequence (see details in Table S1 and File S1). The structural gene annotation was determined based on homology against the P. taeda EST database and the reference GenBank protein database using Geneious Pro software (Drummond et al., 2011). The biological function of genes was determined based on their homology with the Arabidopsis thaliana and Pinus taeda protein database (E-value are reported in Table S5) and the Protein Knowledgebase – UniProtKB (http://www.uniprot.org/uniprot/).

558

D. Grivet et al. / Molecular Phylogenetics and Evolution 68 (2013) 555–566

Table 1 Description of the 13 phenotypic traits estimated for the seven Mediterranean pine species. Last column lists relevant references providing evidence of the trait adaptive value in pines or other woody plants. Actual values, measurement units and data sources are given in Table S2. Code

Full name

Description

References

ONT

HEIGHT

Maximum size

Ratio between secondary needles dry mass and total needle dry mass in seedlings at a nursery common garden Mean age for the first effective seed cone bearing, equivalent to minimum generation time Maximum tree height

Jones, 1999; Jordan et al., 2000

AAM

Ontogenetic heteroblasty index Age at maturity

LIFESP WOOD SEED CONE

SAM

Maximum life span Wood density Average seed mass Average volume of female cones Dispersal ability by wind Size at maturity

BARK

Age of oldest trees Wood specific gravity: ratio of air-dried weight to fresh volume Average mass of a single seed without its wing Calculated from mean cone length and width, assumed an ellipsoid shape

Grotkopp et al., 2002, 2004; Charnov and Berrigan, 1990; Verdú, 2002 Keeley and Zedler, 1998; Tapias et al., 2004 Charnov and Berrigan, 1990 Van Gelder et al., 2006; Poorter et al., 2010 Cornelissen, 1999; Leishman et al., 2000 Benkman and Miller, 1996

Ratio between average wing length and average seed mass

Benkman, 1995

Bark thickness

Average height at the onset of reproduction, either from direct records (own data) or calculated by crossing AAM with age-height growth curves Average breast-height bark thickness at age of 40 yr

SEROT

Cone serotiny

Average proportion of serotinous cones in wild stands

SPROUT

Sprouting ability

GEN

Mean genome size

Categorical index: 0 no sprouting; 1 juvenile sprouting from dwarf shoots; 2 adult epicormic sprouting Genome size (in pg) estimated using laser flow cytometry

Berrigan and Charnov, 1994; Bonser and Aarssen, 2009 Keeley and Zedler, 1998; Tapias et al., 2004 Keeley and Zedler,1998; Schwilk and Ackerly, 2001 Keeley and Zedler, 1998; Tapias et al., 2004 Grotkopp et al., 2004

DISPERS

Table 2 Polymorphism levels of 21 low-copy nuclear genes used to infer the phylogeny of the eight pine species, along with the putative protein they code for. Biological function of the putative protein is provided in Table S5. Locus 0_259_01 0_771_01 2_2702_01 0_2885_01 0_3261_01 0_3600_01 0_4394_01 0_6293_01 0_7916_01 0_8479_01 0_8531_01 0_9383_01 2_9930_01 0_10054_01 0_11919_01 0_13929_02 0_16860_01 4cl dhn2 lp31 lp33 Total

Length (bp)a Total

Coding regionb

SNP Total

Coding regionb

555 383 379 395 421 370 391 363 398 545 466 437 399 396 423 336 509 868 823 811 365 10 033

117 na 330 102 420 369 na 222 348 378 375 435 399 258 0 219 na 684 513 375 195 5739

25 40 27 22 20 16 30 29 35 24 22 25 11 20 34 19 45 52 65 66 38 665

22 na 22 4 20 16 na 16 25 11 15 25 11 14 na 8 na 28 32 25 22 316

Accession number Putative protein Metacaspase 1 na Glyoxal oxidase-related protein Inorganic phosphate and probable anion transporter Pentatricopeptide repeat-containing protein Phenylalanine ammonia-lyase na Uncharacterized protein Xyloglucanendotransglucosylase/hydrolase protein Probable inositol transporter Monocopper oxidase-like protein Ubiquitin carboxyl-terminal hydrolase-like protein Transcription factor Transcription factor na Protein basic pentacysteine na 4-coumarate-CoA ligase Dehydrin ABA and WDS induced proteinc ABA and WDS induced proteinc

na: not available. a Estimates based on sequences aligned across species, including indels. b Size given as to correspond to full amino acids (triplet of nucleotides). c Family of plant proteins induced by water deficit stress (WDS), or abscisic acid (ABA) stress and ripening.

2.4. Phylogenetic analysis To infer the pine genealogy, a 2-step method was used. First, a Bayesian hierarchical analysis was performed to estimate species trees from multilocus sequences using BEST version 3.1.2 (Liu and Pearl, 2007). The evolutionary model of each locus was estimated independently using PAUP version 4b10 (Swofford, 1991) and MrModeltest2 version 2.3 (Nylander, 2004). Each partition (i.e. each locus) was set to have its own parameters of evolution using the function unlink for tree topology, branch lengths, state frequencies (nucleotide frequencies) and gene mutation rates. Pos-

terior probability of phylogenetic trees was obtained using Markov Chain Monte Carlo (MCMC) procedures that approximate the posterior probability of trees by drawing samples from the posterior distribution. Parameters for the MCMC were ngen = 30  106 (number of generations), samplefreq = 1000 (sampling frequency of the chain), nchains = 4 (1 cold chain and 3 heated chains) and nruns = 2 (number of runs in which the partition was sampled). Two runs were considered to have converged when they met two criteria: (i) when the average standard deviation of split frequencies for the tree corresponding to each gene was below 0.10; and (ii) when the log probability of the generation versus

D. Grivet et al. / Molecular Phylogenetics and Evolution 68 (2013) 555–566

the log probability of observing the data began to plateau, which suggests that the two runs mixed. Burn-in (i.e. number of samples discarded from the beginning of the chain) was determined by counting the samples for which the two runs did not mix. Finally, the analysis was performed twice (analysis 1 and analysis 2) starting from different seeds to ensure that the true phylogenetic tree was inferred. Second, a Bayesian concordance analysis (BCA) was performed using Bucky 1.2 (Larget et al., 2010). This software takes as input the complete tree files of each individual locus generated by BEST and performs a Bayesian analysis that produces an output consisting of a sample of gene trees from their joint distribution, from which concordance factors (the proportion of genes that have the same clade) are estimated with credibility intervals. Finally, these concordance factors (CFs) are used to produce the concordance tree. Burn-in was identical to that for BEST. With BCA, each locus is assumed to have a unique genealogy, with different loci having different genealogies. The a priori level of discordance among loci is controlled by the parameter a that was set to one (high discordance), based on the reticulate evolution generally found in pines (Willyard et al., 2009). 2.5. Detecting selection in gene trees As selection tests are based on comparing coding and noncoding regions, adaptive patterns could be examined only in loci for which the structural annotation was possible (see Section 3.2.). Moreover, because we used consensus sequences as input, we chose to ignore ambiguous nucleotides only when included in a pairwise comparison (i.e. they were not ignored in the entire data set). 2.5.1. Neutral substitution test This test was performed to assess whether there are more or less replacement or silent substitutions than expected under neutrality. A significant test may be evidence for: (i) positive selection reflected by an increased number of observed replacements, but it is not known whether it is directional or non-directional selection; or (ii) negative selection whereby an increased number of silent substitutions is observed. The neutral substitution test was performed using CRANN 1.04 (Creevey and McInerney, 2003) for each amplicon individually (coding region) and using the best species tree inferred from Bayesian analyses (see Section 2.4.). To assess whether the observed number of replacement or silent substitutions was significantly greater than expected from neutrality, two tests, the G-test and the Fisher’s exact test, were carried out. Following the recommendations of the program manual, if the Fisher’s exact test returned any results then it was taken as the true P-value; otherwise the G-test was used. 2.5.2. Relative rate ratio test This test allows assessing whether there is directional selection (high number of replacement-invariable substitutions; i.e. those replacement substitutions where the new character-state is preserved in all subsequent lineages) or non-directional selection (high number of replacement-variable substitutions; i.e. the replacement substitution is not preserved in all lineages and has changed at least once more in a subsequent lineage). The relative rate ratio test (Creevey and McInerney, 2002) was conducted with CRANN 1.04 (Creevey and McInerney, 2003) for each amplicon individually (coding region) and using the best inferred phylogenetic tree, as reported above. 2.5.3. PAML analysis The program codeml implemented in the package PAML 4.4 (Yang, 2007) was used to compare different codon-based models

559

that are based on the x ratio (the ratio of non-synonymous dN to synonymous dS substitution rates). The x ratio, if heterogeneous among lineages, indicates a violation of the neutral model of sequence evolution. Likelihood ratio tests were used to compare the models. The best species’ tree inferred from the Bayesian analyses (see Section 2.4.) was used for the input tree topology, while coding regions of each gene were analyzed for synonymous and non-synonymous sites. Two analyses were performed: 2.5.3.1. Branch model This model allows one to test whether the x ratios are different among lineages by comparing the likelihood under the one-ratio model (same x for all branches) and the free-ratio model (one x for each branch). The degrees of freedom (df) used equals the total number of branches within the phylogeny (df of the free-ratio model) minus one (df of the one-ratio model). 2.5.3.2. Site specific model This model allows the x ratio to vary among sites, i.e. among codons of a gene. Two pairs of models were examined: the M1a (nearly neutral) – M2a (positive selection) model and the M7 (beta) – M8 (beta and x) models (Yang et al., 2000). In the M7 and M8 models, beta represents the shape of the x distribution; M7 does not allow for positively selected sites, while M8, by adding one extra class of sites to the beta model, accounts for the possible occurrence of positively selected sites. Comparing the two models of each model pairs allows testing for the presence of positive sites using a likelihood-ratio test. In both tests df = 2 was used. When a likelihood test was significant, the Bayes Empirical Bayes (BEB; Yang et al., 2005) implemented in PAML was used to compute the posterior probabilities for site classes and identify sites under positive selection. 2.6. Phenotypic traits evolution To reconstruct the evolutionary history of life-history traits, the proportion of the likelihood associated with each of the alternative states at the nodes was estimated using the BayesMultiState program implemented in the BayesTraits 1.0 package (available at www.evolution.rdg.ac.uk; Pagel and Meade, 2006). The tree file contains the tree with branch lengths inferred from the Bayesian analyses; while the trait files contain individual traits coded as in Table 1. To estimate the model of evolution and ancestral states of the traits, Markov-chain Monte Carlo (MCMC) methods were used. The model using the reversible-jump MCMC method was used as suggested in the program documentation with a gamma distribution of interval 0–10 for the mean and the variance, 107 reiterations, sampled every 300 iterations and a burn-in of 105 (to produce 30,000 sampled points). The model was run several times to reach a level of acceptance rate (ratedev) between 20% and 40% (deviation of the normal distribution) as recommended in the program documentation. The reconstruction of the ancestral state at each node was performed using the Addnode command. The correlated evolution among pairs of traits was inferred using Felsenstein’s (1985) method of phylogenetic independent contrasts (PIC) implemented in Mesquite v. 2.75 (Maddison and Maddison, 2011), using the PDAP module. This method allows detecting whether species traits co-vary once correlation due to ancestry is removed, by integrating the phylogenetic relationships of the species. Before performing the correlation analysis we ensured that each individual trait met the assumptions required to calculate the PIC (many of which are related to the Brownian motion model of evolution; see Table S3). Two traits (HEIGHT and WOOD) violated the most critical assumption and were not included in the analysis; two other traits (DISPERS and SPROUT) violated another assumption that is not too critical to performing the

560

D. Grivet et al. / Molecular Phylogenetics and Evolution 68 (2013) 555–566

PIC, so these traits were included in the analysis (Table S3). For all eleven traits we computed pairwise phylogenetic regressions and reported the 2-tailed P-value for the least squares regression line, as well as the direction (slope) and the strength (R-squared) of the significant correlations. In addition to PIC, we performed Pearson’s correlations uncorrected by the phylogeny between all pairs of variables (based on the original values in Table S2, not on the codes used for PICs) in order to compare both the sign and strength of the correlations, both considering the phylogeny, or not.

3. Results 3.1. Species phylogeny The two phylogenetic trees inferred with the Bayesian method BEST using 21 low-copy genes (10,033 bp and 665 SNP; Table 2) showed posterior probability values of the partitions ranging from 96% to 100% for the first analysis (Fig. 1) and from 91% to 100% for the second analysis (Fig. S1). The two analyses gave the same species trees, and the two runs converged for each analysis (see details in Table S4). The concordance trees inferred from the Bayesian concordance analysis (BCA) grouped the eight species identically to the BEST output and each node was supported by sample-wide posterior mean concordance factors ranging from 0.75 to 1 (first analysis; Fig. 1) and 0.81 to 1 (second analysis; Fig. S1). These results led to a well-supported phylogeny and indicated that the dataset contains enough information for inferring the Mediterranean pine species tree. The inferred phylogenetic tree for the seven Mediterranean pines was retained for further analyses investigating adaptive evolution of candidate genes and life-history traits.

3.2. Detecting selection in gene trees Three amplicons were not included in the analysis because their coding and noncoding regions could not be determined (0_771_01, 0_4394_01 and 0_16860_01), and one amplicon was discarded because it did not display any coding region (0_11919_01) (see Table 2).

The Neutral Substitution test detected one gene (0_259_01), out of the 17 genes that could be tested, that showed negative selection (Table 3). This event of selection was detected at the outer branches of the trees. More specifically, the amplicon 0_259_01 displayed a significant departure from neutrality (i.e. fewer observed non-synonymous sites than expected) at internal branches connecting node 2 (Fisher’s P-value = 0.0336) and node 3 (Fisher’s P-value = 0.0043), and was almost significant at the internal branch connecting node 4 (G-test with P-value <0.100). This amplicon that gave significant results with the Neutral Substitution test did not show any significant results with the Relative Rate Ratio test of directionality. The Branch Model test gave, for individual genes, dN/dS for the entire gene tree under the null model (x0) ranging from 0.0001 (amplicon 2_9930_01) to 0.9705 (amplicon 4cl) (Table 3). Based on global dN/dS, all these genes would have undergone negative selection. When looking at dN/dS across lineages, the Branch Model test revealed a significant variation of x for lp31 (P-value = 0.0453), and almost significant for 0_259_01 and dhn2 (P-value = 0.0956 and 0.0563, respectively), but no specific pattern emerged in the phylogeny. However, locus lp31 is interesting as x for P. brutia and P. roxburghii showed evidence of positive selection (x > 1) (Fig. 2). Finally, the Site Specific Model revealed an almost significant variation of dN/dS for one codon (at position 40 of the alignment, coding for amino acid Glu/Gln) within the amplicon 0_259_01 that could be under positive selection (P-value = 0.0654 when examining the M7 and M8 models; Table 3). 3.3. Evolution of life-history traits The reconstruction of evolutionary history of the 12 pine lifehistory traits and the genome size on the Bayesian consensus tree is shown in Fig. 3, and the most relevant results can be summarized as follows: 3.3.1. Growth and development Mediterranean pines showed a trend towards relatively delayed ancestral shoot ontogeny, as reflected by a low proportion of secondary needles (ONT; state 1), with P. heldreichii showing the high-

Table 3 Detection of selective events in Mediterranean pines. 3 out of 21 gene loci could not be annotated and were not included in these analyses.

a b

Locus

x0a

dN/dS

Type of selectionb

P-value

0_259_01

0.0615

0.0077/0.1256

Negative (internal branch 2) Negative (internal branch 3) Negative (internal branch 4) Branch site (pos. 40 E/Q)c

2_2702_01 0_2885_01 0_3261_01 0_3600_01 0_6293_01 0_7916_01 0_8479_01 0_8531_01 0_9383_01 2_9930_01 0_10054_01 0_13929_02 4cl dhn2 lp31 lp33

0.2874 0.2481 0.7132 0.4186 0.1347 0.1330 0.4989 0.4671 0.0001 0.0743 0.5079 0.9705 0.1922 0.1387 0.5692

0.0120/0.0416 0.0000/0.0000 0.0036/0.0147 0.0226/0.0317 0.0260/0.0621 0.0149/0.1108 0.0038/0.0285 0.0192/0.0385 0.0185/0.0397 0.0000/0.0348 0.0048/0.0648 0.0187/0.0369 0.0100/0.0103 0.0186/0.0967 0.0315/0.2268 0.0349/0.0613

0.0336 0.0043 0.100 > P > 0.050 0.0956 0.0654 ns ns ns ns ns ns ns ns ns ns ns ns ns 0.0563 0.0453 ns

Branch Branch

x0: Estimation of global dN/dS under the null model (estimated using PAML).

Type of selection: negative (negative selection detected with the Neutral Substitution test using CRANN); branch (dN/dS variation across branches detected using PAML) and site (dN/dS variation across sites detected using PAML). c Position of the triplet coding for the amino acid E (Glu, glutamic acid)/Q (Gln, glutamine). ns: non-significant at P-value = 0.10.

D. Grivet et al. / Molecular Phylogenetics and Evolution 68 (2013) 555–566

561

Fig. 3. Evolution of life-history traits on the consensus tree. Posterior probabilities of Bayesian inference character state evolution are represented by pie charts. See Table S2 for details on trait state coding (zero reflects the absence of the trait for a given taxon, and 1 the smallest value).

est value (state 5). Age at Maturity’s ancestral state (AAM) is intermediate (state 2) with P. halepensis and P. brutia showing the

smallest values (state 1) and P. heldreichii the highest (state 5). Maximum size (HEIGHT) ancestral state was moderate for the

562

D. Grivet et al. / Molecular Phylogenetics and Evolution 68 (2013) 555–566

Mediterranean pines’s cluster (state 3), but showed higher values in the path leading to P. canariensis and P. roxburghii. Maximum life span (LIFESP) ancestral state was moderate (state 3), P. halepensis showing the smallest value (state 1) and P. heldreichii the highest (state 5). Wood density (WOOD) ancestral state showed low values (state 1), and no clear evolutionary pattern emerged across Mediterranean pines. 3.3.2. Reproduction A contrasting pattern of evolution was found within the Mediterranean pines with the cluster P. canariensis-roxburghii-pinea showing higher seed mass (SEED) and cone size (CONE), and lower dispersal ability by wind (DISPERS) than P.halepensis/brutia or P. heldreichii. For size at maturity (SAM) the ancestral state was towards high values (state 4), P. brutia and P. halepensis displaying the lowest values (state 1). 3.3.3. Fire tolerance Mediterranean pines showed medium to high levels of fire tolerance for 2 fire traits, serotiny (SEROT) and sprouting ability (SPROUT), with the exceptions of P. pinea and P. heldreichii (that do not produce serotinous cones) for the former, and P. heldreichii (that does not resprout, even at a juvenile stage) for the latter. Bark thickness (BARK) displayed a high value as the ancestral trait (state 3), and this state was conserved in most Mediterranean pine taxa, albeit to a lower extent in P. heldreichii (state 2) and particularly in P. halepensis (state 1). Overall, there is no specific pattern of life-history trait evolution within the Mediterranean pines’ clade, i.e. they have followed distinct evolutionary pathways in the different species. Furthermore, we did not find a common pattern across Mediterranean pine species, although there was an overall trend towards higher fire tolerance (thick bark, as well as high serotiny and sprouting ability) as compared with the extended pine phylogeny in He et al. (2012) (see Section 4.3). Moreover, with respect to genome size, no clear pattern of evolution was found, with most species exhibiting intermediate values (average of 32.68 pg), but still higher than the average in subgenus Pinus (27.02 pg; Grotkopp et al., 2004).

Remarkably, P. halepensis displayed distinct phenotypic values compared to its sister species P. brutia for 9 traits, with P. brutia generally displaying values more in line with other Mediterranean pines than to P. halepensis. Six (Pearson, uncorrected) to 15 (PIC, corrected) correlations (2– 5 per trait) were detected, revealing that all traits except SPROUT and genome size co-varied with other traits (Fig. 4). When integrating the phylogenetic relationships (i.e. PIC analysis), 4 out of 6 pairwise correlations persisted, while two correlations did not show up due to coancestry (AAM-LIFESP and SAM-SEROT); moreover, 11 new correlations were detected. Particularly, several significant correlations occurred between traits linked to dispersal (DISPERS covaried with 5 traits) and fire tolerance (SEROT and BARK covaried with 5 traits each). 4. Discussion 4.1. Species phylogeny The set of 21 nuclear low-copy genes simultaneously analyzed using Bayesian methods provided a well-supported species phylogeny of the closely related Mediterranean pines. It has been shown previously that different genes sampled from the same set of taxa can produce quite distinct phylogenies (see references in Ané et al., 2007). Indeed, individual genes can lead to discordant phylogenic trees due to processes such as horizontal transfer, incomplete lineage sorting, introgression, hybrid speciation, or/ and selection. In Mediterranean pines, gene flow among taxa has been documented but is relatively rare (e.g. between P. brutia and P. halepensis; Barbéro et al., 1998), which does not preclude, however, that gene flow could have blurred the phylogenetic signal when this group radiated 10 Mya. Moreover, incomplete lineage sorting of ancestral polymorphisms could have produced discordances between gene and species trees in the subsection Pinaster. To take into account this latter process when estimating the underlying species phylogeny from the set of 21 low-copy genes, we used a method that explicitly incorporates incomplete lineage sorting (Liu and Pearl, 2007). Summarizing the phylogenetic signal

Fig. 4. Significant pairwise correlations between 10 life-history traits and genome size (P < 0.05). Grey line: standard Pearson correlations (not corrected by phylogenic relationships); black line: Felsenstein’s phylogenetic correlations; double black line: correlations that were significant with both methods. Thick lines represent significant correlations with R-squared >0.50 and thin lines significant correlations with R-squared <0.50. Trait abbreviations are given in Table 2. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

D. Grivet et al. / Molecular Phylogenetics and Evolution 68 (2013) 555–566

from multilocus sequence data – by jointly estimating gene trees and incorporating uncertainty in individual gene tree estimates– was successful in producing a well-supported phylogenetic tree providing the basis to investigate the complex evolutionary history of the Mediterranean pine assembly. Our multilocus low-copy gene data confirmed previously known evolutionary relationships established with molecular markers, and in addition, resolved uncertainties found in previous studies. Accordingly, we confirmed the grouping of P. brutia and P. halepensis as previously defined with chloroplast markers (e.g. López et al., 2002). These two mostly parapatric circum-Mediterranean species form a uniform group with high taxonomic affinity (Quézel, 2000), some ecological and genetic similarities, and when they (rarely) co-occur they can form natural hybrids (Barbéro et al., 1998; Quézel, 2000; Schiller, 2000). The present study also confirmed the position of P. halepensis and P. brutia as sister group to the clade formed by P. canariensis, P. roxburghii, P. pinea, and P. pinaster. Finally, it corroborated the early-diverging position of P. heldreichii in subsection Pinaster, as previously inferred with chloroplast markers (e.g. Gernandt et al., 2005). The Bosnian pine (P. heldreichii) has been considered by some authors to be more closely related to other Asian hard pines than to Mediterranean pines, while other studies supported its affinity with the Mediterranean pines (see discussion in Wang et al., 1999). The incorporation of two species belonging to the subsection Pinus (P. sylvestris and P. nigra) in our nuclear-based phylogeny confirms the grouping of P. heldreichii with subsection Pinaster (see Fig. S2). Our analyses could establish with strong support the position of four taxa whose position with respect to other Mediterranean pines was ambiguous according to previous molecular and morphological studies. More specifically, P. canariensis and P. roxburghii formed a tight group, despite these species being geographically distant by more than 7000 km. Pinus canariensis and P. roxburghii have been previously considered close relatives by some authors (e.g. Page, 1974; Eckert and Hall, 2006), while for others this relation was not so clear (e.g. Liston et al., 1999). Two scenarios could explain the actual geographic location of these two taxa: (i) a very early split of P. roxburghii lineage from the Mediterranean pines, implying a common ancestor with a wide Eurasian distribution that would have disappeared during the climatic changes at the end of the Miocene (Page, 1974; Klaus, 1989), and (ii) P. roxburghii as an ancestral lineage to the Eurasian hard pines (Mirov, 1967; Page, 1974; Klaus, 1989; Wang et al., 1999), a scenario that violates tree-thinking as it would designate a currently extant species as the ancestor. Most studies based on chloroplast data (e.g. Eckert and Hall, 2006; Gernandt et al., 2008) and our own results based on nuclear loci support the former hypothesis, because they place Eurasian pines (e.g. P. sylvestris and P. nigra) in an early-diverging position compared to P. roxburghii. Pinus pinea’s position in our well-supported phylogenetic tree was sister to P. canariensis and P. roxburghii. This species is considered by many authors to be enigmatic and isolated (Mirov, 1967; Klaus, 1989), and its taxonomic position was unclear (compare, for instance, Liston et al., 1999 with Wang et al., 1999). Finally, P. pinaster’s position was found sister to that of P. pinea + P. canariensis + P. roxburghii. This species presents as many specific morphological features that separate it from other Mediterranean pines as those that link it to them (Klaus, 1989). The position of P. pinaster and P. pinea in the phylogenetic tree was the most difficult to assess (i.e. it required increasing the number of low-copy genes to obtain a well-supported phylogeny). Using multiple genes allows minimizing the misleading signal that some individual genes may generate when reconstructing a phylogenetic tree and allows gathering a set of loci able to resolve both deep and shallow nodes (Edwards, 2009a). In the present study, a new set of low-copy nuclear genes permitted examining

563

multiple evolutionarily independent events and resolving the relationships among the recently diverged taxa that define subsection Pinaster, providing a solid phylogeny to study the evolution of molecular variants and phenotypic traits relevant to the adaptation of Mediterranean pines. 4.2. Detecting selection in gene trees Estimating selective events in a phylogeny gives access to long evolutionary scales that may have affected one or more lineages. In the present study, 2 genes out of 17 showed different significant selective patterns (at P = 0.05). (i) The sequence of locus 0_259_01 is highly homologous to that of a metacaspase (Table S5) which is involved in plant defense response by acting as a positive regulator of cell death (Coll et al., 2010). This gene displayed purifying selection in the lineages leading to the group formed by P. pinea-roxburghii-canariensis. Purifying selection in this group can be explained by the action of negative selection in an ancestral branch (e.g. ancestral lineage connecting node 4 that is almost significant) and the persistence of its footprint in all subsequent descendants (i.e. ancestral lineages connecting nodes 3 and 2) (Fig. 1); another scenario would involve the action of negative selection in successive ancestral lineages (i.e. ancestral lineages leading to nodes 4, 3, and 2). It is not surprising to find that this gene is experiencing purifying selection, since many defense response genes are functionally important and have been conserved (see, for instance, Bakker et al., 2008 where 27 defense response genes in Arabidopsis thaliana were analyzed). (ii) Locus lp31 belongs to the ASR family (abscisic acid/water stress/ripening induced; Table S5) and is involved in stress response (Frankel et al., 2006). Locus lp31 showed a highly variable dN/dS ratio among evolutionary lineages that is incompatible with a neutral model of sequence evolution. Results for individual lineages (Fig. 2) indicated two cases of positive selection (x > 1 for P. brutia and P. roxburghii). Moreover, differences in the amount of selective constraint among tree branches cannot be ruled out to explain the variation in dN/dS among lineages. Interestingly, at the intra-specific level lp31 shows signatures of homogenizing selection in P. pinaster at range-wide distribution scale (Eveno et al., 2008), of positive selection in P. halepensis eastern populations (Grivet et al., 2009) and of diversifying selection in P. sylvestris western Scotland populations (Wachowiak et al., 2011). This gene seems thus to have been the target of distinct selective forces at both the intra- and interspecific levels depending on species and environment. The power of detecting signals of selection at the branch and codon levels is based on nucleotide diversity at synonymous (dS) and non-synonymous (dN) sites. In our study, the moderate level of synonymous and non-synonymous substitutions (77% of the 175 branches of the entire gene tree across all genes displayed dS = 0 and dN = 0; Fig. S3), is compensated by the relatively high number of lineages analyzed simultaneously, which increases the power of sites-based approaches (Yang and Bielawski, 2000). Our analyses highlight the challenge in detecting selection by looking at dN/dS ratios averaged over all sites and all lineages, as selection rarely affects all sites of a gene over a prolonged time (Nielsen and Yang, 1998). Accordingly, the assumption that all sites in a sequence and all branches in a phylogeny are under the same underlying dN/dS ratio is unrealistic, and analyses at branch and site levels provide a more effective approach to detect selective events in a phylogeny (as exemplified in this study). In the present study, 12% of the genes displayed signatures of selection in the Mediterranean pine phylogeny, of which half was due to purifying selection. This proportion is rather low compared to a previous study based on 21 nuclear genes in the phylogeny of 13 pine taxa from the subgenus Pinus and Strobus (Palmé et al., 2009). These authors found little evidence of positive selection

564

D. Grivet et al. / Molecular Phylogenetics and Evolution 68 (2013) 555–566

but observed extensive action of negative selection, i.e. in nearly all the genes examined. The low proportion of genes under selection found in our study can be related to the relatively low taxonomic level at which selection tests were performed, which reduces their power due to low divergence among lineages. This lack of divergence was reflected by the difficulties in resolving the phylogenetic relationships within the Mediterranean pine section. A lack of positive selection was also shown within the genus Pinus, in a multilocus analyses of major clades of soft pines (subgenus Strobus) (Eckert A.J., personal communication). The same amplicons examined for signatures of selection were also used to infer the Mediterranean pine phylogeny. This approach may be questioned as genes under selection may produce incorrect species phylogenies. However, in our study a high proportion of loci produced phylogenetic trees with the same clades (i.e. a high value of sample-wide mean concordance factor), which supports a unique multilocus tree topology. Moreover, a recent review suggests that natural selection should not constitute a problem for phylogenetic analysis (Edwards, 2009b). Indeed, while some selective forces would have minor effect in producing misleading phylogenies (e.g. stabilizing selection or directional selection), others could hinder reconstructing a phylogeny correctly (e.g. balancing selection or selection-driven convergence) but would need to act at the genome scale to impact the phylogenetic reconstruction and can therefore be uncovered by looking at many genes simultaneously (Edwards, 2009b). Accordingly, mulilocus data sets allow recovering the signals of many sites across many genes and minimize the erroneous phylogenetic trees that could result from selective events. These views are supported by the well-supported phylogeny obtained in this study for the Mediterranean pines, even though two genes (out of the 21 used in the phylogeny) were found potentially under selection. 4.3. Evolution of life-history traits Several insights emerged from our analysis of adaptive trait evolution. Overall, there was not concerted evolution of traits across Mediterranean pine taxa, a result that reflects the distinct life histories of this heterogeneous assembly as already stated by other authors (e.g. Klaus, 1989; Keeley and Zedler, 1998; Tapias et al., 2004). However, there were a few remarkable exceptions, as exemplified by the evolution towards decreasing size at maturity in the lineage leading to P. halepensis/brutia, which suggests that heterochrony (i.e. the heritable alteration of the rates of development; e.g. Smith, 2003) has played some role in the divergence of these taxa, probably as an adaptation to severe perturbations under limited resources. Moreover, two adaptive fire-related traits showed a specific pattern of evolution: high bark thickness was maintained in most species and decreased in only two taxa (in P. heldreichii and particularly in P. halepensis), and serotiny was maintained or increased across Mediterranean lineages, except in P. pinea and P. heldreichii. The evolutionary pathways of these two fire-adapted traits reflect the adaptive evolution of Mediterranean pines to ‘fiery’ environments, as recently found in a study examining 101 pine species (He et al., 2012) showing that taxa of this group have evolved towards thick fire-protective bark, as well as high cone serotiny. Finally, juvenile sprouting exemplifies a trait that was maintained in most Mediterranean taxa except in Canary Island pine (P. canariensis) where it increased (and exists at both juvenile and adult stages) due to intense volcanism (Climent et al., 2004), and in Bosnian pine (P. heldreichii) where it was lost. Juvenile sprouting could have been selected by distinct selective pressures beside fire, such as freezing, strong winds, drought or browsing (Keely et al., 2011). Our study indicates that its ancestral state in Mediterranean pines is towards intermediate sprouting capacity (see details in Table S2). This result contrasts with that

of He et al. (2012) who looked at resprouting capacity at the adult stage and considered only P. canariensis and P. roxburghii as resprouting taxa, a situation that would lead to an ancestral state with very low or no resprouting capacity in the Mediterranean pine clade. Apart from the cases reported above, traits under examination did not show any historical or evolutionary adaptive convergence in Mediterranean lineages. One interesting case is the atypical phenotypic values displayed by P. halepensis with regards to its sister species P. brutia: although these taxa have high taxonomic affinities, they seem to have diverged quite distinctly for several key adaptive traits (9 out of 13), P. halepensis displaying divergent phenotypes with regard to other Mediterranean pines. At this point, it is still unclear whether the divergence in trait evolution between P. halepensis and P. brutia could be attributed to distinct demographic histories (extreme bottlenecks and long-range colonization in P. halepensis, see Grivet et al., 2009, 2011; distinct roles of the two taxa in vegetation dynamic, and human influences, see Quézel, 2000), to differences in ecological behavior (in particular edaphic features and bioclimatic requirements, Quézel, 2000), or to both. Several genetic correlations were found for the Mediterranean pines, some involving known or expected relationships such as the inverse correlation between dispersal ability by wind and seed mass (Benkman, 1995). Other less obvious correlations seem evolutionarily relevant, such as the negative correlations between seed dispersal ability by wind and bark thickness as well as between serotiny and maximum life span, and finally the positive correlation between serotiny and seed dispersal. All three can be interpreted as genetic trade-offs between maximizing adult endurance and reproducing early in size terms, keeping an abundant aerial seed bank and dispersing seeds more efficiently after disturbances. This pattern is highly consistent with the postulated divergences in life history strategies in pines and other woody perennials (e.g. Keeley and Zedler, 1998; Richardson, 1998; Tapias et al., 2004). In addition our results confirm the evolutionary tradeoff between cone serotiny and bark thickness in Mediterranean pines as previously mentioned by other authors (Keeley and Zedler, 1998; Schwilk and Ackerly, 2001). Finally, although no specific evolutionary pattern was detected for genome size in terms of character evolution or coevolution, one result worth mentioning is the evolution of the genome towards larger size in Mediterranean pines (average of 32.68 pg compared to 27.02 pg in subgenus Pinus; Grotkopp et al., 2004). A study comparing 85 pine taxa showed that over evolutionary time, genome size has remained stable or increased in soft pines (subgenus Strobus), while it has decreased in most hard pine (subgenus Pinus) subsections (Grotkopp et al., 2004). Mediterranean pines belong to subgenus Pinus, subsection Pinaster, and constitute thus a noteworthy exception. Moreover, we did not find any correlations between genome size and the other traits examined. This result contrasts with other studies that found significant correlations with variables such as minimum generation time (in 18 North American pines, Wakamiya et al., 1993; in the genus and subgenus Pinus, Grotkopp et al., 2004), seed mass (Grotkopp et al., 2004), as well as growth indices and seed dimensions (Wakamiya et al., 1993). The absence of correlations in our study may be related to the low taxonomic level we are looking at (there is little variation in the mean genome size across the studied species, as this is generally the case for gymnosperms, see Leitch and Leitch, 2012), and/or to the rather low number of species examined (i.e. seven taxa) that may limit our ability to detect significant correlations. By investigating ecological patterns in a phylogenetic framework we free ourselves from covariation due to historical constraint, but we do not integrate historical effects (long-term dynamic nature of species assemblage) and lineage sorting (diversification rates of lineages), two other processes than can affect

D. Grivet et al. / Molecular Phylogenetics and Evolution 68 (2013) 555–566

covariation of life-history traits besides adaptive changes (Herrera, 1992). 4.4. Conclusion The examination of a common set of 21 low-copy loci allowed us to infer phylogenetic relationships in Mediterranean pines and to construct a well-supported phylogeny. It is expected that, as sequencing capacity increases (Neale and Kremer, 2011), unraveling the evolutionary history of pines at all taxonomic levels will be within reach. Based on the Mediterranean pine phylogeny, we detected two genes under selection that may have played a role in Mediterranean pines’ defense and stress responses. These genes constitute good candidates for specific gene expression experiments to further examine their functional role in pine adaptation strategies. At the phenotypic levels, with few exceptions, most life-history traits were found to have followed independent evolution pathways in Mediterranean pines, leading to a heterogeneous assembly of species. Finally, pairwise trait correlations taking into account phylogenetic constraints revealed various evolutionary trade-offs linking growth-development, reproduction and fire-related traits. Identifying key genomic and phenotypic elements involved in Mediterranean pines’ adaptation is central to understanding their past evolutionary success and their potential response in the face of a changing climate. Acknowledgments We are grateful to C. García-Barriga and S. Torre for their help in producing or editing some of the sequences. We thank M.R. Chambel and E. Hermoso for their contribution in producing the phenotypic data, as well as A. Eckert, J. Wegrzyn, O. Savolainen and K. Avia for facilitating some of the molecular data. Thanks are extended to G. Nieto-Feliner, J. Fuertes, J.G. Pausas, R. Alía, and Z. Lorenzo for their comments on early drafts, and to P. Alizoti, Y. Kurt, J.J. Robledo-Arnuncio and S. Mutke for checking and providing valuable suggestions on the phenotypic data. This work was supported by the EU EVOLTREE Network of Excellence, project CGL2011-30182-C02-01 (ADAPCON) from the Spanish Ministry of Economy and Competiveness, as well as project AT07-002 of INIA, and the Italian MIUR project ‘‘Bioversitalia’’ (RBAP10A2T4). D.G. acknowledges the support of the Spanish Ministry of Economy and Competiveness through a ‘Ramón y Cajal’ fellowship. Appendix A. Supplementary material Supplementary data associated with this article can be found, in the online version, at http://dx.doi.org/10.1016/ j.ympev.2013.03.032. References Ané, C., Larget, B., Baum, D.A., Smith, S.D., Rokas, A., 2007. Bayesian estimation of concordance among gene trees. Mol. Biol. Evol. 24, 412–426. Bakker, E.G., Traw, M.B., Toomajian, C., Kreitman, M., Bergelson, J., 2008. Low levels of polymorphism in genes that control the activation of defense response in Arabidopsis thaliana. Genetics 178, 2031–2043. Barbéro, M., Loisel, R., Quézel, P., Richardson, D.M., Romane, F., 1998. Pines of the Mediterranean Basin. In: Richardson, D.M. (Ed.), Ecology and Biogeography of Pinus. Cambridge University Press, Cambridge, pp. 153–170. Benkman, C.W., 1995. Wind dispersal capacity of pine seeds and the evolution of different seed dispersal modes in pines. Oikos 73, 221–224. Benkman, C.W., Miller, R.E., 1996. Morphological evolution in response to fluctuating selection. Evolution 50, 2499–2504. Berrigan, D., Charnov, E.L., 1994. Reaction norms for age and size at maturity in response to temperature: a puzzle for life histories. Oikos 70, 474–478. Bonser, S.P., Aarssen, L.W., 2009. Interpreting reproductive allometry: individual strategies of allocation explain size-dependent reproduction in plant populations. Perspect. Plant Ecol. Evol. Syst. 11, 31–40.

565

Charnov, E.L., Berrigan, D., 1990. Dimensionless numbers and life history evolution: age of maturity versus the adult lifespan. Evol. Ecol. 4, 273–275. Climent, J., Tapias, R., Pardos, J.A., Gil, L., 2004. Fire adaptations in Canary Islands pine (Pinus canariensis). Plant Ecol. 171, 185–196. Climent, J., San-Martín, R., Chambel, M.R., Mutke, S., 2011. Ontogenetic differentiation between Mediterranean and Eurasian pines sect. Pinus. at the seedling stage. Trees-Struct. Funct. 25, 175–186. Coll, N.S., Vercammen, D., Smidler, A., Clover, C., Breusegem, F.V., Dangl, J.L., 2010. Arabidopsis Type I metacaspases control cell death. Science 330, 1393–1397. Cornelissen, J.H.C., 1999. A triangular relationship between leaf size and seed size among woody species: allometry, ontogeny, ecology and taxonomy. Oecologia 118, 248–255. Creevey, C.J., McInerney, J.O., 2002. An algorithm for detecting directional and nondirectional positive selection, neutrality and negative selection in protein coding DNA sequences. Gene 300, 43–51. Creevey, C.J., McInerney, J.O., 2003. CRANN: detecting adaptive evolution in proteincoding DNA sequences. Bioinformatics 19, 1726. De Mita, S., Santoni, S., Ronfort, J., Bataillon, T., 2007. Adaptive evolution of the symbiotic gene NORK is not correlated with shifts of rhizobial specificity in the genus Medicago. BMC Evol. Biol. 7, 210–218. Drummond, A., Ashton, B., Buxton, S., et al. (15 coauthors) 2011. Geneious URL: . Eckert, A.J., Hall, B.D., 2006. Phylogeny, historical biogeography, and patterns of diversification for Pinus (Pinaceae): phylogenetic tests of fossil-based hypotheses. Mol. Phylogenet. Evol. 40, 166–182. Edwards, S.V., 2009a. Is a new and general theory of molecular systematics emerging? Evolution 63, 1–19. Edwards, S.V., 2009b. Natural selection and phylogenetic analysis. Proc. Natl. Acad. Sci. USA 106, 8799–8800. Eveno, E., Collada, C., Guevara, M.A., 2008. Contrasting patterns of selection at Pinus pinaster Ait. drought stress candidate genes as revealed by genetic differentiation analyses. Mol. Biol. Evol. 25, 417–437 (11 authors).. Fady, B., 2012. Biogeography of neutral genes and recent evolutionary history of pines in the Mediterranean Basin. Ann. Forest Sci. 69, 421–428. Felsenstein, J., 1985. Phylogenies and the comparative method. Am. Nat. 125, 1–15. Frankel, N., Carrari, F., Hasson, E., Iusem, N.D., 2006. Evolutionary history of the Asr gene family. Gene 378, 74–83. Gernandt, D.S., Lopez, G.G., Garcia, S.O., Liston, A., 2005. Phylogeny and classification of Pinus. Taxon 54, 29–42. Gernandt, D.S., Magallon, S., Lopez, G.G., Flores, O.Z., Willyard, A., Liston, A., 2008. Use of simultaneous analyses to guide fossil-based calibrations of Pinaceae phylogeny. Int. J. Plant Sci. 169, 1086–1099. González-Martínez, S.C., Wheeler, N.C., Ersoz, E., Nelson, C.D., Neale, D.B., 2007. Association genetics in Pinus taeda L.I. wood property traits. Genetics 175, 1–11. Grivet, D., Sebastiani, F., González-Martínez, S.C., Vendramin, G.G., 2009. Patterns of polymorphism resulting from long-range colonization in the Mediterranean conifer Aleppo pine. New Phytol. 184, 1016–1028. Grivet, D., Sebastiani, F., Alía, R., Bataillon, T., Torre, S., Zabal-Aguirre, M., Vendramin, G.G., González-Martínez, S.C., 2011. Molecular footprints of local adaptation in two Mediterranean conifers. Mol. Biol. Evol. 28, 101–116. Grotkopp, E., Rejmanek, M., Rost, T.L., 2002. Toward a causal explanation of plant invasiveness: seedling growth and life-history strategies of 29 pine Pinus species. Am. Nat. 159, 396–419. Grotkopp, E., Rejmánek, M., Sanderson, M.J., Rost, T.L., 2004. Evolution of Genome Size in Pines (Pinus) and its life-history correlates: supertree analyses. Evolution 58, 1705–1729. Guzmán, B., Lledó, M.D., Vargas, P., 2009. Adaptive radiation in Mediterranean Cistus (Cistaceae). PloS One 4, e6362. He, T., Pausas, J.G., Belcher, C.M., Schwilk, D.W., Lamont, B.B., 2012. Fire-adapted traits of Pinus arose in the fiery Cretaceous. New Phytol. 194, 751–759. Herrera, C.M., 1992. Historical effects and sorting processes as explanations for contemporary ecological patterns: character syndromes in Mediterranean woody plants. Am. Nat. 140, 421–446. Jones, C.S., 1999. An essay on juvenility, phase change, and heteroblasty in seed plants. Int. J. Plant Sci. 160, S105–S111. Jordan, G.J., Potts, B.M., Chalmers, P., Wiltshire, R.J.E., 2000. Quantitative genetic evidence that the timing of vegetative phase change in Eucalyptus globulus ssp globulus is an adaptive trait. Aust. J. Bot. 48, 561–567. Keeley, J., Zedler, P., 1998. Evolution of life history in pines. In: Richardson, D.M. (Ed.), Ecology and Biogeography of Pinus. Cambridge University Press, Cambridge, pp. 219–242. Keely, J., Bond, W.J., Bradstock, R.A., Pausas, J.G., Rundel, P.W., 2011. In: Fire in Mediterranean Ecosystems. Cambridge University Press, Cambridge, pp. 58–80. Klaus, W., 1989. Mediterranean pines and their history. Plant Syst. Evol. 162, 133– 163. Larget, B.R., Kotha, S.K., Dewey, C.N., Ané, C., 2010. BUCKy: Gene tree/species tree reconciliation with Bayesian concordance analysis. Bioinformatics 26, 2910– 2911. Leishman, M.R., Wright, I.J., Moles, A.T., Westoby, M., 2000. The evolutionary ecology of seed size. In: Fenner, M. (Ed.), Seeds: the Ecology of Regeneration in Plant Communities. UK: CAB International, Trowbridge, pp. 31–59. Leitch, A.R., Leitch, I.J., 2012. Ecological and genetic factors linked to contrasting genome dynamics in seed plants. New Phytol. 194, 629–646. Liston, A., Robinson, W.A., Piñero, D., Alvarez-Buylla, E.R., 1999. Phylogenetics of Pinus (Pinaceae) based on nuclear ribosomal DNA internal transcribed spacer region sequences. Mol. Phylogenet. Evol. 11, 95–109.

566

D. Grivet et al. / Molecular Phylogenetics and Evolution 68 (2013) 555–566

Liu, L., Pearl, D.K., 2007. Species trees from gene trees: reconstructing Bayesian posterior distributions of a species phylogeny using estimated gene tree distributions. Syst. Biol. 56, 504–514. López, G.G., Kamiya, K., Harada, K., 2002. Phylogenetic relationships of diploxylon pines (Subgenus Pinus) based on plastid sequence data. Int. J. Plant Sci. 163, 737–747. Maddison, W.P., Maddison, D.R., 2011. Mesquite: a Modular System for Evolutionary Analysis URL: . Martínez-Vilalta, J., Sala, A., Piño, J., 2004. The hydraulic architecture of Pinaceae – a review. Plant Ecol. 171, 3–13. Mirov, N.T., 1967. The Genus Pinus. Ronald Press Company, New York, USA. Neale, D.B., Kremer, A., 2011. Forest tree genomics: growing resources and applications. Nat. Rev. Genet. 12, 111–122. Ne’eman, G., Goubitz, S., Werger, M.J., Shmida, A., 2011. Relationships between tree size, crown shape, gender segregation and sex allocation in Pinus halepensis, a Mediterranean pine tree. Ann. Bot. 108, 197–206. Nielsen, R., Yang, Z., 1998. Likelihood models for detecting positively selected amino acid sites and applications to the HIV-1 envelope gene. Genetics 148, 929–936. Nylander, J.A.A., 2004. MrModeltest URL: . Page, C., 1974. Morphology and affinities of Pinus canariensis. Notes Roy. Bot. Gard. Edinburgh 33, 317–323. Pagel, M., Meade, A., 2006. Bayesian analysis of correlated evolution of discrete characters by reversible-jump Markov Chain Monte Carlo. Am. Nat. 167, 808– 825. Palmé, A.E., Pyhäjärvi, T., Wachowiak, W., Savolainen, O., 2009. Selection on nuclear genes in a Pinus phylogeny. Mol. Biol. Evol. 26, 893–905. Parks, M., Cronn, R., Liston, A., 2009. Increasing phylogenetic resolution at low taxonomic levels using massively parallel sequencing of chloroplast genomes. BMC Biol. 7, 84–100. Pausas, J., Bradstock, R., Keith, D., Keeley, J., 2004. Plant functional traits in relation to fire in crown-fire ecosystems. Ecology 85, 1085–1100. Poorter, L., McDonald, I., Alarcón, A., Fichtler, E., Licona, J.-C., Peña-Claros, M., Sterck, F., Villegas, Z., Sass-Klaassen, U., 2010. The importance of wood traits and hydraulic conductance for the performance and life history strategies of 42 rainforest tree species. New Phytol. 185, 481–492. Pot, D., McMillan, L., Echt, C., Le Provost, G., Garnier-Géré, P., Cato, S., Plomion, C., 2005. Nucleotide variation in genes involved in wood formation in two pine species. New Phytol. 167, 101–112. Quézel, P., 2000. Taxonomy and biogeography of Mediterranean pines. In: Ne’eman, G., Trabaud, L. (Eds.), EcologyBiogeography and Management of Pinus Halepensis and P. brutia Forest Ecosystems in the Mediterranean Basin. Backhuys Publishers, Leiden, The Netherlands, pp. 1–12. Richardson, D.R., 1998. Ecology and Biogeography of Pinus. Cambridge University Press, Cambridge, UK. Rochman, M.V., 2012. The QTN program and the alleles that matter for evolution: all that’s gold does not glitter. Evolution 66, 1–17. Schiller, G., 2000. Inter- and intra-specific genetic diversity of Pinus halepensis Mill. and P. brutia Ten. In: Ne’eman, G., Trabaud L. (Eds), Ecology, Biogeography and

Managment of Pinus Halepensis and P. brutia Forest Ecosystems in the Mediterranean basin. Backhuys Publishers, Leiden, The Netherlands, pp. 13–35. Schwilk, D.W., Ackerly, D.D., 2001. Flammability and serotiny as strategies: correlated evolution in pines. Oikos 94, 326–336. Smith, K.K., 2003. Time’s arrow: heterochrony and the evolution of development. Int. J. Dev. Biol. 47, 613–621. Stearns, S.C., 1977. The Evolution of life histories traits: a critique of the theory and a review of the data. Ann. Rev. Ecol. Evol. Syst. 8, 145–171. Strauss, S.H., Ledig, F.T., 1985. Seedling architecture and life-history evolution in pines. Am. Nat. 125, 702–715. Swofford, D.L., 1991. PAUP: Phylogenetic Analysis Using Parsimony URL: . Syring, J., Willyard, A., Cronn, R., Liston, A., 2005. Evolutionary relationship among Pinus (Pinaceae) subsection inferred from multiple low-copy nuclear loci. Am. J. Bot. 92, 2086–2100. Tapias, R., Gil, L., Fuentes-Utrilla, P., Pardos, J.A., 2001. Canopy seed banks in Mediterranean pines of South-eastern Spain: a comparison between Pinus halepensis Mill., P. pinaster Ait., P. nigra Arn. and P. pinea L. J. Ecol. 89, 629–638. Tapias, R., Climent, J., Pardos, J.A., Gil, L., 2004. Life histories of Mediterranean pines. Plant Ecol. 171, 53–68. van Gelder, H.A., Poorter, L., Sterck, F.J., 2006. Wood mechanics, allometry, and lifehistory variation in a tropical rain forest tree community. New Phytol. 171, 367–378. Verdú, M., 2002. Age at maturity and diversification in woody angiosperms. Evolution 56, 1352–1361. Wachowiak, W., Salmela, M.J., Ennos, R.A., Iason, G., Cavers, S., 2011. High genetic diversity at the extreme range edge: nucleotide variation at nuclear loci in Scots pine Pinus sylvestris L. in Scotland. Heredity 106, 775–787. Wakamiya, I., Newton, R.J., Johnston, J.S., Price, H.J., 1993. Genome size and environmental factors in the genus Pinus. Am. J. Bot. 80, 1235–1241. Wang, X.-R., Tsumura, Y., Yoshimaru, H., Nagasaka, K., Szmidt, A.E., 1999. Phylogenetic relationships of Eurasian pines Pinus (Pinaceae) based on chloroplast rbcL, matK, rpl20-rps18 spacer, and trnV intron sequences. Am. J. Bot. 86, 1742–1753. Wegrzyn, J.L., Lee, J.M., Tearse, B.R., Neale, D.B., 2008. TreeGenes: A Forest tree Genome Database. Int. J. Plant Genom. Article ID 412875. Willyard, A., Cronn, R., Liston, A., 2009. Reticulate evolution and incomplete lineage sorting among the ponderosa pines. Mol. Phylogenet. Evol. 52, 498–511. Yang, Z., 1998. Likelihood ratio tests for detecting positive selection and application to primate lysozyme evolution. Mol. Biol. Evol. 15, 568–573. Yang, Z., 2007. PAML 4: phylogenetic analysis by maximum likelihood. Mol. Biol. Evol. 24, 1586–1591. Yang, Z., Bielawski, J.P., 2000. Statistical methods for detecting molecular adaptation. Trends Ecol. Evol. 15, 496–503. Yang, Z., Nielsen, R., Goldman, N., Pedersen, A.M., 2000. Codon-substitution models for heterogeneous selection pressure at amino acid sites. Genetics 155, 431– 449. Yang, Z., Wong, W.S.W., Nielsen, R., 2005. Bayes empirical Bayes inference of amino acid sites under positive selection. Mol. Biol. Evol. 24, 1107–1118.