# Microgeographic adaptation and the effect of pollen flow on the adaptive potential of a temperate tree species

Julie Gauzere, Etienne K Klein, Oliver Brendel, Hendrik Davi, Sylvie

Oddou-muratorio

► To cite this version:

Julie Gauzere, Etienne K Klein, Oliver Brendel, Hendrik Davi, Sylvie Oddou-muratorio. Microgeographic adaptation and the effect of pollen flow on the adaptive potential of a temperate tree species. New Phytologist, 2020, 227 (2), pp.641 - 653. 10.1111/nph.16537. hal-04168537

HAL Id: hal-04168537

https://hal.inrae.fr/hal-04168537

Submitted on 21 Jul 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

<figure><img src="image_3.png" /></figure>

# Microgeographic adaptation and the effect of pollen flow on the adaptive potential of a temperate tree species

**Julie Gauzere**¹,²,³ ID, **Etienne K. Klein**² ID, **Oliver Brendel**⁴ ID, **Hendrik Davi**¹ ID and **Sylvie Oddou-Muratorio**¹ ID

¹INRAE, URFM, Avignon 84000, France; ²INRAE, BioSP, Avignon 84000, France; ³Institute of Evolutionary Biology, School of Biological Sciences, University of Edinburgh, Edinburgh, EH9 3JT, UK; ⁴Université de Lorraine, AgroParisTech, INRAE, UMR Silva, Nancy 54000, France

Author for correspondence:
**Julie Gauzere**
Tel: +44 131 650 7702
Email: julie.gauzere@ed.ac.uk

Received: 8 November 2019
Accepted: 29 February 2020

New Phytologist (2020) 227: 641–653
doi: 10.1111/nph.16537

**Key words:** adaptive divergence, climate change, common garden, elevation gradient, *Fagus sylvatica*, genetic diversity, natural selection, pollen dispersal.

## Summary

* In species with long-distance dispersal capacities and inhabiting a large ecological niche, local selection and gene flow are expected to be major evolutionary forces affecting the genetic adaptation of natural populations. Yet, in species such as trees, evidence of microgeographic adaptation and the quantitative assessment of the impact of gene flow on adaptive genetic variation are still limited.

* Here, we used extensive genetic and phenotypic data from European beech seedlings collected along an elevation gradient, and grown in a common garden, to study the signature of selection on the divergence of eleven potentially adaptive traits, and to assess the role of gene flow in resupplying adaptive genetic variation.

* We found a significant signal of adaptive differentiation among plots separated by < 1 km, with selection acting on growth and phenological traits. Consistent with theoretical expectations, our results suggest that pollen dispersal contributes to increase genetic diversity for these locally differentiated traits.

* Our results thus highlight that local selection is an important evolutionary force in natural tree populations and suggest that management interventions to facilitate movement of game-tes along short ecological gradients would boost genetic diversity of individual tree popula-tions, and enhance their adaptive potential to rapidly changing environments.

## Introduction

Recent interest in microgeographic adaptation – adaptation at spatial scales compatible with substantial amounts of gene disper-sal – suggests that we must reconsider the scale at which evolution occurs (Richardson et al., 2014). At this fine spatial scale, a major question that remains to be answered is whether gene flow is con-straining or facilitating local adaptation. Too much gene flow would overwhelm the effects of natural selection that drive adap-tive genetic divergence among populations (Lenormand, 2002) and decrease the adaptation of populations to their local environ-ment (migration load). Conversely, gene flow also can resupply the genetic diversity eroded by genetic drift and selection, and thus facilitate future evolutionary responses (Pease et al., 1989; Barton, 2001; Lenormand, 2002; Bridle et al., 2010; Polechova & Barton, 2015). In many taxa, and especially long-lived organ-isms, gene flow is expected to be the main process maintaining the standing genetic variation within natural populations (Bar-ton, 2001; Le Corre & Kremer, 2003; Bridle et al., 2010). In the particular case of environments varying through both space and time, gene flow also could spread alleles pre-adapted to future environmental conditions (Davis & Shaw, 2001; Aitken et al., 2008; Kremer et al., 2012; Fitzpatrick & Reid, 2019). Although theoretical developments on this topic have been rich and

provided interesting hypotheses to test in nature, empirical stud-ies on the role and importance of gene flow in adaptation are still lacking (but see Fitzpatrick et al., 2015, 2016).

Increasing evidence of microgeographic adaptation in nature provides strong support for the idea that gene flow rarely prevents the emergence of local adaptation or adaptive divergence (Ander-son et al., 2015; Eckert et al., 2015; Fitzpatrick et al., 2015; Moody et al., 2015; Peterson et al., 2016; for a review, see Richardson et al., 2014). This evidence strengthens the conclu-sions of studies at larger spatial scales, showing that geographic or genetic distances among populations are not correlated with the level of local adaptation or adaptive genetic divergence (McKay & Latta, 2002; Leimu & Fischer, 2008; Hereford, 2009). These results suggest that selective pressures can be strong enough to counterbalance the effects of gene flow, and to shape genetic vari-ation at a fine spatial scale. However, evidence of microgeo-graphic adaptation is reported mainly for certain taxa (e.g. fish, Fitzpatrick et al., 2015; Moody et al., 2015; and annual or peren-nial herbs, Anderson et al., 2015; Peterson et al., 2016), whereas others, such as temperate tree species, are under-represented. In trees, there is a long tradition of large-scale adaptation studies, through provenance trials, following a quantitative genetic approach (Savolainen et al., 2007; Alberto et al., 2013). The development of genome-wide data and statistical methods has
considerably changed this research field, and evidence of adaptations at regional and landscape scales are now flourishing (e.g. Csillery et al., 2014; Pluess et al., 2016). Yet, the adaptive response of phenotypic traits at the local scale is still poorly understood (but see Brousseau et al., 2013; Eckert et al., 2015). Methodological limitations that previously hampered the robust testing of local phenotypic adaptation have been largely overcome (Ovaskainen et al., 2011), permitting more accurate assessments of the scale at which phenotypic adaptations occur.

Empirical evidence of beneficial effects of gene flow on adaptive potential in nature are still rare, and often are assessed in laboratory settings (e.g. Swindell & Bouzat, 2006). Some experimental studies on inbred plant populations have demonstrated the beneficial effect of gene flow on the reduction of homozygosity and reintroduction of variation for fixed deleterious alleles (e.g. Costa e Silva et al. 2014; Bontrager & Angert, 2019). However, for outbred species, the role of gene flow on the maintenance of genetic diversity has been difficult to demonstrate, especially when relying on natural crosses. This effect of gene flow is thought to be particularly important when the populations exchanging genes are strongly genetically differentiated (Pease et al., 1989; Barton, 2001). Long-distance dispersal events could therefore play a major role in resupplying the genetic variation of populations and favour the spread of beneficial alleles (Savolainen et al., 2007; Kremer et al., 2012). In wind-pollinated species (anemophilous), these long-distance dispersal events are common, with mean pollen dispersal distances that can easily reach hundreds of metres, with rare events of spread over tens of kilometres in some species (Austerlitz et al., 2004; Petit & Hampe, 2006; Ashley, 2010; Kremer et al., 2012). Pollen dispersal may thus be highly efficient for spreading beneficial alleles

between populations (Kremer *et al.*, 2012). Even though a pollen grain brings half as many immigrating alleles as does a seed, gene flow by pollen or seed is expected to have similar consequences on a population's adaptation and migration load (Lopez *et al.*, 2008).

Wind-pollinated tree species are therefore good study systems to look for empirical evidence of the beneficial effects of long-distance gene flow on adaptive capacities in nature. However, to our knowledge, only two studies have looked for such evidence. In *Pinus contorta*, Yeaman & Jarvis (2006) showed a positive correlation between the level of genetic variance within populations for growth and the spatial heterogeneity of selection. This pattern is consistent with the theoretical expectation that gene flow among differentiated populations is an important source of genetic variation. Along a latitudinal gradient of genetically differentiated populations of *Pinus sylvestris*, Nilsson (1995) found that offspring harvested after one episode of reproduction showed a shift in average genetic values as compared to their mothers. This pattern is expected in the case of directional gene flow among locally differentiated populations (Fig. 1). This effect of directional pollen dispersal on an offspring's average genetic values for adaptive traits could be particularly beneficial if the environmental conditions of the recipient population change toward the one of the donor population (Davis & Shaw, 2001; Aitken *et al.*, 2008).

The aim of the present study is to investigate the effects of natural selection and long-distance gene flow on the adaptive potential of a major temperate tree species, the European beech (*Fagus sylvatica* L.). This anemophilous species can disperse pollen at long distance (Piotti *et al.*, 2011), a feature thought to explain the low level of neutral genetic differentiation among populations

<figure><img src="image_4.png" /><figcaption>**Fig. 1** Expected effects of random and directional gene flow between locally differentiated populations on individual genetic values. (a) In the case of random gene flow, migrant gene flow (orange arrows) tends to homogenize the frequencies of the red, grey and blue genotypes in each environment. Hence, in the core population (e.g. grey habitat), individuals and gametes from migrant gene flow should exhibit greater genetic variance for the traits under local selection ($v_1 > v_2$) than individuals and gametes from local gene flow (i.e. gene flow within population, represented by the green arrows), whereas the mean genetic value (dotted line) should be identical (assuming a gradual shift of the optimum). (b) In the case of directional gene flow (e.g. only the red genotypes migrate to the grey habitat), migrant gene flow should introduce individuals and gametes with shifted mean genetic values (dotted lines; $m_1 
eq m_2$), whereas the variance in genetic values between migrant and local gene flow should be similar (assuming that stabilizing selection is the same in each habitat).</figcaption></figure>
observed from the local to the European scale (Comps et al., 2001; Buiteveld et al., 2007). Several phenotypic traits show sig-
nificant genetic differentiation among beech provenances at the
European scale (Robson et al., 2018), in particular budburst date
(Gomory & Paule, 2011; Kramer et al., 2017), height (Rose
et al., 2009) and, to a lesser extent, wood anatomy, hydraulic
traits and physiological foliar traits (Knutzen et al., 2015; Hajek
et al., 2016; Stojnic et al., 2018). These patterns are usually inter-
preted as a long-term response to selection, conferring better
adaptation of beech populations to local climate and, most
importantly, to summer drought. Our study aimed to test
whether such patterns of phenotypic adaptation also exist at a fine
spatial scale.

Here, we analyzed the phenotypic and genotypic data of
~2300 beech seedlings from 60 maternal families grown in a
common garden. These open-pollinated families were collected
from three natural beech plots, spreading along an elevation gra-
dient at fine spatial scale (c. 1.5 km in-between the most distant
plots) at the southern, warm margin of this species' distribution.
We focused on 11 potentially adaptive phenotypic traits, for
which we previously found a significant genetic component in
phenotypic variation (Gauzere et al., 2016a). The first ques-
tion investigated was: (1) Is there genetic divergence across
elevations for the quantitative traits of interest? We expected
trees from the low-elevation, warmer and drier plot to have
evolved traits favouring adaptation to drought, and trees from
the high-elevation, colder plot to have evolved traits favouring
adaptation to low temperatures. Our second question was: (2)
Is there a signature of local selection in the differentiation
pattern of the trait studied. We tested this by comparing
traits differentiation to neutral differentiation shown by
microsatellite markers. Our final question was: (3) Does gene
flow contribute to the additive variance of traits under selec-
tion at the local scale? To that purpose, we used parentage

<figure><img src="image_5.png" /><figcaption>**Fig. 2** Climatic space explored by the three study plots in comparison to the whole species range of *Fagus sylvatica*. This climatic space is represented as the variation of temperature and summer precipitation over the beech distribution area, extracted from the SAFRAN database for the period 1958–2015 (collected on a 8 km² grid represented by grey dots). Each black triangle represents the average climate recorded since 2007 with HOBO weather stations for plots N1 (1020 m above sea level (asl)), N2 (1140 m asl) and N4 (1340 m asl).</figcaption></figure>

analysis to identify the offspring likely originating from local
or distant pollen immigration events, and we tested how
long-distance pollen dispersal affected the genetic values of
offspring after one episode of reproduction.

**Materials and Methods**

Species, sampling plots and experimental design

*Fagus sylvatica* L. (European beech) is a monoecious, anemophilous and predominantly outcrossed tree species (Merzeau *et al.*, 1994; Wang, 2003). It is a shade-tolerant species requiring well-drained and moderately deep soils and relatively high humidity rates (Breda *et al.*, 2006; Jump *et al.*, 2006). Its distribution ranges from the northern Mediterranean regions to the south of Scandinavia. On Mont-Ventoux, a mountain in southeastern France, beech forests are located at the southern limit of their ecological range (Fig. 2), in a mountainous Mediterranean ecosystem. On the northern face of Mont-Ventoux, the beech forest ranges almost continuously from 750 to 1700 m above sea level. This steep elevation gradient provides almost linear variation in mean temperature and humidity with elevation (Davi *et al.*, 2011). Three plots were defined along this climatic gradient, at 1020 m above sea level (asl) (N1), 1140 m asl (N2) and 1340 m asl (N4) (Supporting Information Table S1), over a total distance of *c.* 1.5 km. A previous study showed that these plots have high pollen immigration rates (on average 53.3%), moderate distances of pollen dispersal (on average 45.7 m), and non-negligible long-distance pollen dispersal events (Gauzere *et al.*, 2013b).

In August 2009, 20 open-pollinated families were sampled in each of the three plots (60 maternal progenies in total). Mother-trees were chosen for their high fertility and random location within the plot. In April 2010, the seeds were germinated and a subset of 46.4 seedlings on average per family were randomly planted in 25 blocks in a common garden experiment (with 1.8 seedlings per family in each block) at the State nursery of Aix-en-Provence (lat. 43°23′30″N, long. 5°23′24″E). All seedlings were planted in independent pots of 1.2 l with sand substrate and fertilizer, and irrigated regularly. The experiment ran for 3 years (from April 2010 to September 2013). Details on the seed collection, germination protocol and progeny test design can be found in Gauzere *et al.* (2016a).

All of the potentially reproductive adults within the three
plots, including the mother-trees, and an average of 23.9 off-
spring per family were genotyped at 13 microsatellites markers
(690 adults and 1437 offspring in total; see Gauzere *et al*., 2013b,
for the genotyping details).

Phenotypic measurements

A set of 11 different phenotypic traits were measured on the off-
spring grown in the common garden. These traits, already intro-
duced in Gauzere *et al*. (2016a), are presented briefly below (see
also Table 1). Height and diameter growths were measured
between August 2010 and November 2011 (ΔH, ΔD). For the
<table><caption>**Table 1** Description and descriptive statistics for the traits measured on seedlings of *Fagus sylvatica*.</caption><thead><tr><th>Abbreviation</th><th>Trait</th><th>Category</th><th>N<sub>w</sub></th><th>Unit</th><th>Transformation</th></tr></thead><tbody><tr><td>ΔH<sub>log</sub></td><td>Total height increase</td><td>Growth</td><td>2380</td><td>Rate</td><td>log</td></tr><tr><td>ΔD</td><td>Total diameter increase</td><td>Growth</td><td>2380</td><td>Rate</td><td>—</td></tr><tr><td>t<sub>b2→3(2011)</sub></td><td>Budburst date in 2011</td><td>Phenology</td><td>2334</td><td>DOY</td><td>—</td></tr><tr><td>t<sub>b2→3(2012)</sub></td><td>Budburst date in 2012</td><td>Phenology</td><td>1835</td><td>DOY</td><td>—</td></tr><tr><td>t<sub>s1→2</sub></td><td>Senescence date in 2011</td><td>Phenology</td><td>2386</td><td>DOY</td><td>—</td></tr><tr><td>VD<sub>2011</sub></td><td>Vegetation season duration in 2011</td><td>Phenology</td><td>2193</td><td>DOY</td><td>—</td></tr><tr><td>LM<sub>√</sub></td><td>Dry leaf mass</td><td>Morphology</td><td>2491</td><td>g</td><td>Square-root</td></tr><tr><td>LA<sub>√</sub></td><td>Leaf area</td><td>Morphology</td><td>2491</td><td>cm<sup>2</sup></td><td>Square-root</td></tr><tr><td>LMA</td><td>Leaf mass area</td><td>Morphology</td><td>2491</td><td>g m<sup>-2</sup></td><td>—</td></tr><tr><td>% N</td><td>Nitrogen content in leaves</td><td>Physiology</td><td>1031</td><td>%</td><td>—</td></tr><tr><td>δ<sup>13</sup>C</td><td>Water use efficiency</td><td>Physiology</td><td>1031</td><td></td><td>—</td></tr></tbody></table>

DOY, day of the yr. N<sub>w</sub>, number of phenotyped individuals.

whole trial, phenological events also were recorded: the date of budburst in 2011 and 2012 (t<sub>b2→3(2011)</sub>, t<sub>b2→3(2012)</sub>), the date of leaf senescence in 2011 (t<sub>s1→2(2011)</sub>) and the duration of the growing season in 2011 (VD<sub>2011</sub> = t<sub>s1→2(2011)</sub> - t<sub>b2→3(2011)</sub>). At the end of the growing season in 2011, three light-exposed leaves were collected on the stem of each seedling to measure morphological and physiological traits. On fresh leaves, we first measured the leaf area (LA), and put them to dry at 60°C to then measure the leaf dry mass (LA) of each seedling. The leaf mass area, related to the photosynthetic capacity and stomatal conductance of the plant (Reich et al., 2003), was calculated as LMA = LM/LA. Finally, for a subset of 1031 individuals (evenly sampled across families), we measured the leaf carbon isotope composition (δ<sup>13</sup>C), as a surrogate for water use efficiency (Farquhar et al., 1989). The nitrogen content in the leaves (%N) also was estimated for this subset of individuals, as leaf nitrogen density is linked to the photosynthetic capacity of leaves (Kattge et al., 2009). Note that three of these traits were transformed before analyses to limit departure from linear model assumptions: ΔH with a logarithm transformation, and LM and LA with a square-root transformation (following Gauzere et al., 2016a).

## Test for population genetic differentiation

In order to assess the genetic differentiation of quantitative traits among plots we ran a specific linear mixed model (called ‘animal model’) that uses relatedness information from a pedigree to dissociate the contribution of genetic, population and environmental effects on the total phenotypic variance of traits. Thus, we analyzed each trait Y measured in common garden using the following model, and focused particularly on the significance of the fixed plot effect P:

$$Y_{i,f,p,k} = \mu + P_i + B_f + M_p + a_k + \varepsilon_{i,f,p,k} \quad \text{Eqn 1}$$

with, as fixed effects, μ the intercept, P<sub>i</sub> the effect of the plot of origin, B<sub>f</sub> the effect of the block and M<sub>p</sub> the observer effect (only included for the analysis of the phenological traits). The random terms are the additive genetic values of the individuals k,

$\{a_k\} : N(0, AV_A)$ with A the matrix of pairwise relatedness derived from a pedigree, and V<sub>A</sub> is the additive genetic variance, and the residuals, $\{\varepsilon_{i,f,p,k}\} : N(0, I_{de}V_R)$ with I<sub>de</sub> the identity matrix and V<sub>R</sub> the residual variance.

The (co)variance structure of the additive genetic effects was defined using a one-generation pedigree reconstructed in a previous study (Gauzere et al., 2016a). Here, we did not model material effects because of the lack of power to accurately dissociate maternal and additive genetic effects (Gauzere et al., 2016b). Model Eqn 1 was fitted in R/ASREML (Gilmour et al., 2006). Using a similar model, Gauzere et al. (2016a) previously showed that all traits presented moderate narrow-sense heritability, with h<sup>2</sup> ranging from 0.12 to 0.35.

## Signature of natural selection accounting for genetic drift

We used the method developed by Ovaskainen et al. (2011) to detect signatures of natural selection in trait differentiation among populations. This method requires genotypic and phenotypic data, that are first analysed under a neutral model which assumes that the genes coding for the trait are as divergent as the neutral markers (Q<sub>st</sub> = F<sub>st</sub> case). This first model thus considers the population mean genetic additive values as normally distributed, with the covariance between pairs of population means being proportional to the pairwise average coancestry between these populations (i.e. assuming neutrality for the trait):

$$A^P : N(0, 2V_A \times \mathbf{0}^P) \quad \text{Eqn 2}$$

with A<sup>P</sup> the vector of the population-level effects, V<sub>A</sub> the additive genetic variance of the trait, 0<sup>P</sup> the matrix of pairwise population-level coancestry coefficients. Note that the 0<sup>P</sup> matrix provides an estimation of F<sub>st</sub>.

The method then considers an animal model for each quantitative trait Y, that decomposes the additive value into a population and an individual effect. The population-level additive effect is determined by the mean genetic additive values in the parent population and the individual-level additive effects a<sub>k</sub> are defined the same way as in model Eqn 1. The presence of selection on the
trait is tested with the S statistics, which evaluate whether the realized pattern of population mean genetic additive values is likely under the neutral model Eqn 2.

In order to apply this method to our dataset, we first used the genotyping of the adult trees and the admixture F-model for neutral divergence among populations implemented in R/RAFM, to estimate 0^P and the index of neutral genetic differentiation, F_st (Karhunen & Ovaskainen, 2012). We then analyzed the traits measured in the common garden using the same fixed effects as in model Eqn 1 and the one-generation pedigree using R/DRIFT-SEL, to estimate V_A, A^P and S (Karhunen et al., 2013). These Bayesian analyses were repeated three times per trait to ensure convergence of the Markov chain Monte Carlo (MCMC) chains. We provide the average posterior S estimates (the chains provided very similar S values; Table S2). We also measured a low autocorrelation within the chains, indicating that the space of the parameter values is explored efficiently.

## Immigrant characterization based on molecular markers

**Categorical assignments of paternity** We first tried to identify the father of the genotyped seedlings to characterize their immigration status. To that purpose, we used the likelihood-based method of categorical paternity assignment implemented in CERVUS (Marshall et al., 1998) with allelic frequencies estimated from the genotypes of all adult trees. All of the reproductive adults from the three sampled plots were considered as candidate fathers, but only fathers within the mother's plot were assigned (see below). We considered 0% typing error, 100% sampling of candidate fathers, a confidence level of 95%, and allowed selfing (following Gauzere et al., 2016a). These choices were made to favour assignments to genotyped fathers, despite the risk of type I errors (i.e. a wrong tree is assigned while the true father is not sampled), knowing that quantitative genetic methods are robust to pedigree errors (Charmantier & Réale, 2005; Gauzere et al., 2016b). We successfully assigned the paternity for 45.8% of the genotyped offspring. We considered that individuals assigned to a father located within their mother's plot come from 'local' pollen pool, whereas individuals with no compatible father were considered 'immigrant'. Noncategorized individuals include those genotyped at fewer than six markers, offspring not assigned because of compatibility with multiple fathers, and offspring assigned to a father outside their mother's plot. Indeed, given the large number of trees potentially contributing to median- and long-distance mating events, we considered that our ability to retrieve the true father outside the mother's plot was low and presented large expected type I error rate.

**Probability of long-distance pollination** In order to refine the 'immigrant' status, we also used the whole genetic dataset and a noncategorical paternity analysis to assess the likelihood for each genotyped seedlings to belong to long-distance immigrant pollen pool. To that purpose, we estimated the Mendelian transition probabilities (T) that each offspring has the observed genotype g_0 while originating from cross-fertilization of female j with an ungenotyped male from a population with the allelic frequencies

BAF: T(g_0|g_j, BAF). We estimated this transition probability using (1) the global allelic frequencies estimated from a larger-scale study on the whole north-face of Mont-Ventoux (AF_global; Lander et al., 2011), and (2) the local allelic frequencies estimated from the genotyping of all the reproductive trees in the three study plots (AF_local). This allowed us to compute the probabilities that each offspring has the observed genotype g_0 given that it originates from long-distance and local pollination events, T(g_0|g_j, AF_global) and T(g_0|g_j, AF_local), respectively. From these two probabilities, we derived for each offspring the logarithm of the odds (LOD)-score of long-distance immigration vs local immigration:

$$v_{g_0} = \log \left( \frac{T(g_0|g_j, AF_{\text{global}})}{T(g_0|g_j, AF_{\text{local}})} \right) \qquad \text{Eqn 3}$$

$v_{g_0} > 0$ indicating that offspring g_0 is more likely originating from a long-distance pollination event and $v_{g_0} < 0$ that it is more likely originating from a short-distance pollination event.

These transitions probabilities were estimated using the MEMMi model (Gauzere et al., 2013b).

## Impact of pollen immigration on quantitative traits

We explicitly tested whether (1) offspring originating from distant fathers presented more genetic variance for quantitative traits than offspring originating from local fathers, which is expected in the case of migration from genetically differentiated populations in an island model, and whether (2) offspring originating from distant fathers had shifted average genetic values for quantitative traits as compared to the ones originating from local fathers, which is expected in the case of directional pollen immigration from genetically differentiated populations along the gradient (Fig. 1). In both cases (1) and (2), pollen dispersal would result in increasing genetic variance and adaptive capacities. We ran two different models depending on how the immigration status was characterized (categorical or quantitative variable).

First, for each phenotypic trait Y, we looked for differences in the mean and variance of additive genetic values ($a_k$) between 'local' and 'immigrant' groups, using the following animal model:

$$Y_{i,f,p,k} = \mu + P_i + B_f + M_p + a_k + \varepsilon_{i,f,p,k} \qquad \text{Eqn 4}$$

with the same fixed effects than in model Eqn 1, but here the (co)variance structure of the additive genetic effects depend on the immigration classes, with $l$='local' and $m$='immigrant':
$$ \begin{bmatrix} V_{Al} & 0 \\ 0 & V_{Am} \end{bmatrix} $$
Note that the effect of 'local' and 'immigrant' categories on the average trait values was tested and found to be never significant. This fixed effect is thus absent in the final model.

From model Eqn 4, we thus estimated two variance components: $V_{Al}$ and $V_{Am}$.

We then performed a random regression model to investigate the effect of long-distance pollen dispersal, v, on the variation in
additive genetic variance. Random regression models are commonly used in evolutionary biology and breeding to allow additive genetic effects (or any random effect) to vary with a covariate (e.g. environment, age; Nussey *et al.*, 2008). This model defines a random intercept and slope at the additive genetic level, as deviation from the mean linear effect of *v* on *Y*, to model the variation of *V*<sub>A</sub> with *v* as: $V_{A}(v) = V_{Aintercept} + v^{2}$. $V_{Aslope} + 2 \cdot v \cdot cov_{A}$ (intercept, slope). The effects of *v* on $V_{A}$ was tested using the following models:

$$Y_{i,f,p,k} = \mu + P_i + B_f + M_p + \alpha v_k + a_k + \varepsilon_{i,f,p,k} \quad \text{Eqn 5}$$

$$Y_{i,f,p,k} = \mu + P_i + B_f + M_p + \alpha v_k + \tilde{a}_k + \varepsilon_{i,f,p,k} \quad \text{Eqn 6}$$

where $\tilde{a}_k = a_k + \beta_k v_k$ describes the random intercept and slope of *v* on *Y* at the additive genetic level in model Eqn 6. In both models Eqn 5 and 6 the term $\alpha v_k$ is the fixed effect of the LOD-score of long-distance vs local migration on the trait.

From the model Eqn 6, we estimated the variances of the random intercept and slope, and the covariance between them. Note that for each individual *k*, we only have one value of *v* and *Y*, but the model estimates these (co)variance components by comparing phenotypic similarities and the effects of *v* on *Y* between relatives.

The significance of the qualitative and quantitative effects of pollen immigration on $V_{A}$ were tested by comparing the two hierarchical models Eqn 1 and 4, and the two models Eqn 5 and 6, with a log-likelihood ratio test. We performed these analyses with R/ASREML, considering only the information about maternal relatednesses, as we (often) ignore the father's identity of 'immigrant' individuals.

## Relatedness and Nep estimates

We used the genetic dataset to estimate the average paternal relatedness ($\bar{\rho}_p$) and effective number of pollen donors ($Nep = 1/\bar{\rho}_p$) in the 'local' and 'immigrant' pollen clouds. We first extracted the haplotype corresponding to the paternal contribution to each diploid genotype following Gauzere *et al.* (2016a). These paternal genotypes were then used to estimate the pairwise paternal relatednesses among each pair of offspring *k* and *k'* ($\rho_{k,k'}$) using the kinship coefficient of Loiselle *et al.* (1995) implemented in SPAGEDI (Hardy & Vekemans, 2002), and the local allelic frequencies ($AF_{local}$). These coefficients were estimated removing all the selfed individuals to only compare the outcrossed pollen pools.

## Results

### Neutral genetic differentiation among plots and pollen immigration

The admixture *F*-model of Karhunen & Ovaskainen (2012) estimated a posterior $F_{st}$ value of 0.026 [0.024; 0.029], consistent with previous estimates of neutral genetic differentiation among populations measured with a comparable method but from 51

populations over a much wider area of Mont-Ventoux ($F_{st}$ = 0.026; Lander *et al.*, 2011). The coancestry estimates within plots (diagonal elements of the $\theta^P$ matrix) indicated that plot N1 tended to have lower effective population size than plot N2 and N4 (Table 2). Plots N1 and N2 exhibited slightly higher between-plots coancestry coefficient, suggesting that the two closest plots exchanged more gene flow ($\theta_{N1N2} > \theta_{N1N4} \sim \theta_{N2N4}$). Overall, the coancestry estimates were low (maximum $\theta_{N1}$ = 0.057; Table 2) indicating that effective population sizes were globally large. This result is consistent with the large effective number of pollen donors Nep estimated for these sites using the same dataset (Table S1; Gauzere *et al.*, 2013b).

Using paternity assignment, we categorized *n* = 658 offspring as 'local' (45.5%) and *n* = 492 as 'immigrant' (34%). A total of 295 (20%) offspring were outside of these two categories (i.e. ambiguous origin). Note that this 'immigrant' status includes offspring with nonsampled fathers which could be located either at a short or long distance to the mother-tree's plot. We checked that the LOD-score of long-distance vs local immigration was indeed higher on average for the 'immigrant' than 'local' individuals, with $\bar{v}$ = -0.33 and 0.18 for local and immigrant offspring, respectively ($P$-value < 0.001, ANOVA test; Fig. S1).

The analysis of the pairwise coefficient of paternal relatedness showed that within the immigrant pollen pool fathers were significantly less related than within the local pollen pool, with $\bar{\rho}_p$ = 0.011 and $\bar{\rho}_p$ = 0.007 within the local and immigrant pollen pools, respectively ($P$ < 0.001, ANOVA). Consequently, we estimated a higher effective number of pollen donors (Nep) in the immigrant than in the local pollen pool (Nep = 26.4 and 31.3 for local and immigrant pollen pools, respectively).

### Patterns of genetic trait differentiation and test for selection

Significant genetic differentiation across elevations was found for $\Delta H_{log}$, $\Delta D$, $t_{b2\to3(2011)}$ and $t_{b2\to3(2012)}$, whereas marginal differentiation was found for $t_{s1\to2(2011)}$ (as revealed by the *P* effect in model Eqn 1; Fig. 3). Using the method of Karhunen & Ovaskainen (2012), we detected a signal of natural selection for three of these four traits: $\Delta H_{log}$, $\Delta D$ and $t_{b2\to3(2011)}$ (i.e. $S$ > 0.90; Fig. 3). For growth traits, plots at intermediate and high elevations presented the highest height and diameter increments

<table><caption>**Table 2** Matrix of the coancestry coefficients estimated by the admixture *F*-model of Karhunen & Ovaskainen (2012) for the three plots of *Fagus sylvatica*.</caption><thead><tr><th></th><th>N1</th><th>N2</th><th>N4</th></tr></thead><tbody><tr><td>N1</td><td>0.057 [0.033; 0.085]</td><td></td><td></td></tr><tr><td>N2</td><td>0.018 [0.005; 0.035]</td><td>0.029 [0.017; 0.045]</td><td></td></tr><tr><td>N4</td><td>0.012 [0.00; 0.024]</td><td>0.013 [0.002; 0.023]</td><td>0.035 [0.02; 0.052]</td></tr></tbody></table>

The diagonals correspond to the within-plot coancestry coefficient, with its median value and confidence intervals (i.e. values of the estimates that bound 95% of the posterior distribution). Lower diagonal elements are the median coancestry coefficients among-plots, with confidence intervals in brackets.
(ΔH<sub>log N2</sub> = 0.88, ΔH<sub>log N4</sub> = 1.03 and ΔD<sub>N2</sub> = 0.82, ΔD<sub>N4</sub> = 0.72), whereas the site at the lowest elevation showed the lowest growth (ΔH<sub>log N1</sub> = 0.56 and ΔD<sub>N1</sub> = 0.50; Fig. 3). For the phenological traits, buds flushed earlier for plots from high elevations as compared to low elevation plots (t<sub>b2→3(2011)N1</sub> = 101.9; t<sub>b2→3(2011)N2</sub> = 100.3; t<sub>b2→3(2011)N4</sub> = 99.8; Fig. 3). Interestingly, contrary to growth and phenological traits, morphological and physiological traits (except LMA) tended to be slightly less differentiated than expected under neutral processes (with S ∈ [0.24; 0.39]), with a strong conservation of the genetic values between plots (e.g. for the leaf mass LM<sub>√N1</sub> = 0.44, LM<sub>√N2</sub> = 0.45, LM<sub>√N4</sub> = 0.44; Fig. 3).

Effect of pollen immigration on within-plot genetic variation

Using the categorical assignment and the model Eqn 4, we only
estimated marginally significant difference in V_A between ‘local’
and ‘immigrant’ individuals for t_{b2→3(2012)} (P=0.088; for all of
the other traits P>0.1). In that case, we predicted higher V_A
within the ‘immigrant’ than the ‘local’ group (V_{Aimmig} = 16.57
and V_{Alocal} = 3.16; Table 3). Using the random regression model
Eqn 6, we found significant variation in additive genetic vari-
ance with the likelihood of provenance from long-distance pol-
len pool for ΔH_{log} (P=0.005), and a marginally significant
effect for t_{b2→3(2012)} (P=0.078; for all the other traits P>0.1;
Table 3). The estimated variance-covariance terms showed
increasing V_A with higher likelihood of long-distance pollina-
tion event for both ΔH_{log} and t_{b2→3(2012)} (Fig. 4). Nevertheless,
for these traits, the fixed effect of v was not significant, which
indicated no global trend in the variation of breeding values
with the likelihood that an offspring originates from long-dis-
tance pollen dispersal. Interestingly, these two different analyses
highlighted the same trend towards a higher genetic variance for
offspring originating from long-distance dispersal than from
local dispersal.

Discussion

The scale of phenotypic adaptation

This study provides evidence that divergent selection can act on
phenotypic traits at fine spatial scale in temperate tree species.
The study beech plots showed important dispersal capacities,
with pollen dispersal that can easily reach 100 m, a fat-tailed
pollen dispersal kernel and high rates of immigration (Gauzere
et al., 2013b; Oddou-Muratorio et al., 2018), similar to other
anemophilous tree species (see Ashley, 2010; Kremer et al.,
2012, for reviews). Because dispersal distances overlap the spa-
tial scale of the selective environment, our results provide
strong evidence that adaptive genetic differentiation can occur
despite high gene flow (see Fitzpatrick et al., 2015; Peterson
et al., 2016, for recent evidence in other taxa). This also sug-
gests that selective pressures on seedling growth (in diameter
and height) and on budburst timing were strong enough to
counteract the homogenizing effect of gene flow on structuring

the genetic variation at the scale of a few hundred metres up to
one kilometre. We found that the seedlings derived from the
plots at high elevation, inhabiting the lowest temperature con-
ditions, flushed earlier and had a higher height and diameter
growth in the common garden than the plot at low elevation.
For the budburst date, this pattern of genetic differentiation is
opposed to the pattern of phenotypic variation measured
*in situ*, where seedlings from plots at low elevations flushed ear-
lier than plots at high elevations (Davi *et al.*, 2011). Beech pop-
ulations originating from higher longitude or elevation also
have been shown to be genetically earlier in provenance tests
(Teissier du Cros *et al.*, 1988; Gomory & Paule, 2011; Kramer
*et al.*, 2017), suggesting that these populations evolved pheno-
logical traits promoting a longer growing season length. In
forestry, growth is considered a good approximation of the per-
formance of trees in a given condition (e.g. Rehfeldt *et al.*,
1999, 2002). *In situ* and *ex situ*, the highest diameter growth
values were measured at intermediate to high elevations
(Cailleret & Davi, 2011, and present study), suggesting that
the genetic and environmental influences on growth traits are
in the same direction. At this southern edge of the ecological
conditions tolerated by beech (Fig. 2), the reduced allocation to
stem growth at the low-elevation plot is likely an adaptive
response to drought, which has been described previously by
comparing marginal vs central beech populations (Rose *et al.*,
2009). These studies, including ours, shed light on the pheno-
typic traits under selection and underlying adaptation at differ-
ent spatial scales. However, they do not provide a direct
assessment of the strength of natural selection or the ecological
bases of selection, nor do they quantify phenotypic mismatch
and maladaptation, which would require selection analyses
(MacColl, 2011).

The weak adaptive divergence of morphological and physio-
logical traits as compared to growth and phenological traits along
ecological gradients remains an open question. Meta-analyses
highlighted that growth and phenological traits are the most
genetically differentiated phenotypic traits at large spatial scale in
temperate and boreal tree species (Savolainen et al., 2007; Alberto
et al., 2013). Here, our study suggests that this conclusion also
holds at a narrower spatial scale. The existing literature also tend
to indicate that morphological and physiological traits generally
exhibit a lower fraction of genetically controlled phenotypic dif-
ferentiation, and a higher within-population genetic variation
than other traits (Hajek et al., 2016). A putative explanation is
that latitudinal or elevation gradients do not fully capture the
spatial scale of variation for drought stress, even if they capture a
strong variation in temperature. In particular in our study, the
important microheterogeneity in soil water capacity on (Nourtier
et al., 2013) may generate micro-local differentiation of morpho-
logical and physiological traits, blurring the signal of divergent
selection for morphological and physiological traits at a larger
spatial scale (Cubry et al., unpublished). Another possible expla-
nation is that morphological and physiological traits could be
under stabilizing selection (e.g. Lamy et al., 2011). Unfortu-
nately, the power to detect the signal of stabilizing selection is
low for species with low neutral genetic differentiation, such as
<figure><img src="image_6.png" /><figcaption>**Fig. 3** Level of trait adaptive divergence between the three plots located at 1020 m above sea level (asl) (N1), 1140 m asl (N2) and 1340 m asl (N4) for the growth, phenological, physiological and morphological traits measured on seedlings of *Fagus sylvatica*. Histograms represent the *a posteriori* distributions of additive genetic values estimated for each plot. The dotted lines represent the median of these distributions. We also provide the test for population differentiation (*P*) and signal of selection (*S*) for these quantitative traits. An *S* coefficient close to 1 and 0 indicates signal of divergent and stabilizing selection, respectively. ***, *P* < 0.001; **, 0.001 < *P* < 0.01; *, 0.01 < *P* < 0.05; ,, 0.05 < *P* < 0.09; ns, *P* > 0.09.</figcaption></figure>
<table><caption>**Table 3** Effect of pollen immigration on the variation in additive genetic variance for traits measured in seedlings of *Fagus sylvatica*.</caption><thead><tr><th colspan="3">Model Eqn 4</th></tr></thead><tbody><tr><td>Trait</td><td>V<sub>Aimmig</sub></td><td>V<sub>Alocal</sub></td></tr><tr><td>t<sub>b2→3(2012)</sub></td><td>16.57 (8.2)</td><td>3.16 (5.6)</td></tr><tr><th colspan="3">Model Eqn 6</th></tr><tr><td colspan="3">Additive genetic effects for ΔH<sub>log</sub></td></tr><tr><td>ΔH<sub>logintercept</sub></td><td>ΔH<sub>logintercept</sub></td><td>ΔH<sub>logslope</sub></td></tr><tr><td>ΔH<sub>logslope</sub></td><td>0.029 (0.010)</td><td>0.0092 (0.008)</td></tr><tr><td colspan="3">Additive genetic effects for t<sub>b2→3(2012)</sub></td></tr><tr><td>t<sub>b2→3(2012)</sub><sub>intercept</sub></td><td>t<sub>b2→3(2012)</sub><sub>intercept</sub></td><td>t<sub>b2→3(2012)</sub><sub>slope</sub></td></tr><tr><td>t<sub>b2→3(2012)</sub><sub>slope</sub></td><td>7.36 (4.8)</td><td>2.47 (2.18)</td></tr></tbody></table>

This effect was found significant or marginally significant for two genetically differentiated traits: using model Eqn 4 for $t_{b2→3(2012)}$ with $P = 0.088$; using model Eqn 6 for $ΔH_{log}$ with $P = 0.005$ and for $t_{b2→3(2012)}$ with $P = 0.078$. Model Eqn 4 estimates an additive genetic variance ($V_A$) for each 'immigrant' and 'local' group (SE in brackets). Model Eqn 6 estimates the additive genetic variance-covariance terms for the intercept and slope of the trait with the logarithm of the odds (LOD)-score of long-distance immigration vs local immigration $v$ (SE in brackets). The diagonal elements thus provide $V_{Aintercept}$ and $V_{Aslope}$, and the off-diagonal element $cov_A$ (intercept, slope). See Supporting Information Tables S3 and S4 for the model outputs for all phenotypic traits.

trees (Savolainen *et al.*, 2007), and patterns of weak quantitative trait differentiation are difficult to interpret (Whitlock, 2008). With the development of genomic resources, this question of microgeographic adaptation is now addressed mainly using genome-wide association studies (e.g. Eckert *et al.*, 2015). However, associations between genotype, phenotype and environment often provide weaker evidence of adaptive differentiation at the phenotypic levels (Le Corre & Kremer, 2012), probably because of the complex multilocus determinism of the traits underlying phenotypic adaptation. Using a more precise characterization of ecological gradients (e.g. by considering the spatial-scale of water

availability) and new methods to compare the similarity of habitats and phenotypes offer promising prospects to better understand the spatial scale of phenotypic adaptation and the ecological drivers of adaptive divergence (Karhunen *et al.*, 2014; Csillery *et al.*, 2020).

## Effects of pollen dispersal on quantitative trait variation at a local scale

For populations locally adapted and experiencing strong selective pressures that vary through space, theory suggests that moderate gene flow can maintain quantitative genetic variation within these populations (Barton, 2001; Le Corre & Kremer, 2003; Bridle *et al.*, 2010). Although adaptive genetic differentiation is widespread, only a few studies have provided evidence that gene flow is an important process to maintain the high levels of genetic variation usually detected in nature (e.g. Yeaman & Jarvis, 2006). Here, consistent with this theoretical expectation, we detected an effect of pollen dispersal on increased genetic variance for two locally differentiated traits: height growth and budburst date. However, we likely lacked the power to detect this effect for the third heritable trait under selection – diameter growth. Because female flowers are receptive before the emission of local pollen (protogyny) and lower populations flower earlier than higher populations, we hypothesized a directional pollen flow across elevations in Gauzere *et al.* (2013a), which would have resulted in shifted genetic values between immigrant and local individuals (Fig. 1; Nilsson, 1995). Our results did not support this expectation.

We used two different approaches to characterize the pollen immigration status. Both analyses suggested that immigrant pollen brings more genetic variance for locally differentiated traits than local pollen. These effects were quantitatively high, with genetic variance for immigrants more than twice as high as for local offspring. However, this result had large standard errors

<figure><img src="image_7.png" /><figcaption>**Fig. 4** Change in additive genetic variance (VA) with the logarithm of the odds (LOD)-score of long-distance immigration vs local immigration ($v$), estimated using model Eqn 6 in *Fagus sylvatica*. This relationship is represented only for two of 11 traits tested, for which we found that this effect was significant or marginally significant: (left) the log-transformed height growth ($ΔH_{log}$) and (right) the budburst date in 2012 ($t_{b2→3(2012)}$). The solid line represents the estimate, and the dashed lines the 95% confidence intervals around this estimate. From the random regression model Eqn 6, we can estimate the change in $V_A$ with a covariate $x$ as: $V_A(x) = var(\text{trait}_{\text{intercept}} + \text{trait}_{\text{slope}} \cdot x) = V_{\text{Aintercept}} + x^2 \cdot V_{\text{Aslope}} + 2 \cdot x \cdot \text{cov}_A(\text{intercept, slope})$, using the estimates provided in Table 3.</figcaption></figure>
around the estimates, and the significance changed depending on
how the immigration status was characterized. We checked that
'immigrant' individuals originated from more diverse fathers
(higher effective number of pollen donors Nep), and had paternal
genotypes more similar to global than local allelic frequencies
(higher LOD-score of long-distance vs local immigration, v).
Nevertheless, our immigrant status still likely suffers some inac-
curacy. In particular, the low neutral genetic differentiation
among the study beech plots (F<sub>st</sub> = 0.026) and, overall, on the
northern face of Mont-Ventoux (Lander et al., 2011), may
reduce our ability to retrieve the origin of a pollen grain using
molecular markers. Moreover, the low levels of quantitative
genetic differentiation measured along the elevation gradient sug-
gest that the expected effect of pollen flow on genetic values after
one episode of reproduction is necessarily moderate. Here, this
effect may thus be difficult to detect, even knowing the true
immigration status of all offspring.

The approach used here provides a direct quantification of the
impact of contemporary pollen immigration on the genetic vari-
ance of quantitative traits. Although using molecular markers to
distinguish migrants and local offspring could have been done on
numerous existing datasets, to our knowledge no study before
ours had reiterated the approach first used by Nilsson (1995).
However, we cannot estimate the immigration load associated
with this increase in genetic variance for the locally differentiated
traits. Although this genetic load is expected to be beneficial
under rapidly changing environmental conditions (Yeaman &
Jarvis, 2006; Savolainen et al., 2007; Kremer et al., 2012), this
remains to be tested. Applying this approach to other species and
populations may improve our ability to both understand and
quantify the impact of gene flow on the adaptive potential of
quantitative traits under natural conditions. This effect should be
more detectable for populations departing from migration–selec-
tion equilibrium, such as isolated (or marginal) populations, or
transplanted plant populations (see Fitzpatrick et al., 2015, 2016,
for examples in guppies).

The present results about trait differentiation and genetic vari-
ance of migrants and nonmigrants extend previous studies about
pollen dispersal kernels and migration rates (Gauzere et al.,
2013a), trait heritabilities (Gauzere et al., 2016a) and variance of
fecundities (Oddou-Muratorio et al., 2018). Pulling all effects
together suggests building a predictive mechanistic model to both
(1) check whether the observed levels of genetic variance in
migrant pollen are consistent with the amounts of migrant
genetic variance that would be predicted mechanistically and (2)
further investigate the adaptive dynamics at the scale of few gen-
erations (following, e.g., Kuparinen et al., 2010). This work is
beyond the scope of the present study because it still requires a
map of beech density over the whole northern face of Mont-Ven-
toux and an extrapolation model to propose a spatially explicit
distribution of adaptive diversity. Point (2) would require infor-
mation about the variation in selection gradients or adaptive
landscapes with elevation. Unfortunately, empirical estimates of
selection need specific approaches in natural populations of long-
lived species and, hence, remain scarce in forest trees (but see
Bontemps et al., 2017).

Gene flow in the context of climate change

One of the most documented consequences of climate change on tree ecosystem functioning is probably the impact of climate warming on the advancement of spring phenology (Menzel *et al.*, 2006; Fu *et al.*, 2015). Understanding the adaptive or maladaptive value of the response of phenology to climate has become a key focus of evolutionary studies (Tansey *et al.*, 2017; Gauzere *et al.*, 2020). An underrated issue is that these phenological changes also are likely to affect gene flow among populations, and the efficiency of long-distance pollen dispersal. For many temperate tree species, we now have strong evidence that climate warming reduces spring phenology synchrony among individuals experiencing the same environmental conditions (Zohner *et al.*, 2018), but also leads to more uniform spring phenology across climatic gradients (Vitasse *et al.*, 2018). The evolutionary consequences of pollen dispersal may therefore change with increasing climate warming. In this study, we detected an effect of local selection and pollen dispersal on the adaptive potential of traits that are known to affect dispersal or reproduction in anemophilous tree species (tree height and phenology; Klinkhamer *et al.*, 1997; Gauzere *et al.*, 2013a; Oddou-Muratorio *et al.*, 2018). In particular, long-distance pollen dispersal events are only effective if the donor and recipient trees have synchronous phenology (Aitken *et al.*, 2008). Our results thus suggest that refining the models of phenotypic adaptation by incorporating the feedbacks between traits under divergent selection and dispersal–reproductive capacities may be important to thoroughly understand the evolutionary consequences of gene flow in plant species (e.g. Soularue & Kremer, 2012).

With the concern that populations could be maladaptated and
lag behind with increasing climate change, new conservation and
management strategies are needed to mitigate these effects.
Assisted gene flow – the managed movement of individuals or
gametes between populations within species ranges – is one of
the tools proposed to help populations adapt to new climatic
conditions (Aitken & Whitlock, 2013; Aitken & Bemmels,
2016). However, assisted gene flow among distant populations,
at different latitudes, may have unexpected consequences for the
response of individuals or new hybrids, because of outbreeding
depression (Aitken & Whitlock, 2013), but also because environ-
mental cues other than temperature can govern physiological
responses in trees (see, for instance, the role of photoperiod in
bud development; Way & Montgomery, 2015). Here, our results
suggest some effect of local pollen dispersal on increasing genetic
variation for locally differentiated traits after solely one episode of
reproduction. Therefore, programmes considering the move-
ments of gametes along short ecological gradients could be effi-
cient to boost genetic diversity and introduce genetic innovation
within natural or managed tree populations. Our results also
encourage the development of evolution-orientated forest man-
agement and silviculture, whereby thinning could be done to
favour the spread of long-distance gene flow (e.g. by reducing
local densities; Lefèvre *et al.*, 2014). Nonetheless, unlike assisted
gene flow, natural gene flow also can lead to increased maladapta-
tion through migration load (Aitken & Whitlock, 2013).
## Conclusion

By combining population genetics and quantitative genetics, we provided a comprehensive analysis of the effects of gene flow on the adaptive potential of a tree species in nature. We found that gene flow among trees at different elevations did not prevent adaptive differentiation among plots separated by < 1 km, and that divergent selection at this fine spatial scale acted mainly on growth and phenological traits, probably as an adaptation to drought and temperature variation. For two genetically differentiated traits, our analyses suggested an effect of pollen dispersal on increasing genetic diversity after one episode of reproduction. However, our power to detect and estimate these effects probably was quite low. Our results thus suggest that local selection is an important evolutionary force in natural tree populations and that, under rapidly changing environmental conditions, gene flow can increase the adaptive capacities of natural populations.

## Acknowledgements

We are grateful to Ian White for help with the implementation and interpretation of the random regression models in ASREML, and to Tonya Lander for useful comments on the manuscript. We thank Olivier Gilg, Frank Rei, Norbert Turion, Frédéric Jean and Mehdi Pringarbe, INRAE UEFM, Avignon for sample collection, field measurements and seed germination. We also thank Anne Roig and Matthieu Lingrand for genotyping and managing the laboratory work and the employees of the nursery of Aix-Les-Milles for management and measurement work on the progeny trial. Isotopic measurements were performed by C. Hossann at Silvatech (Silvatech, INRAE, 2018. Structural and functional analysis of tree and wood facility, https://doi.org/10.15454/1.5572400113627854E12) from UMR 1434 Silva at the research centre INRAE Nancy-Lorraine. Silvatech facility is supported by the French National Research Agency through the Laboratory of Excellence ARBRE (ANR-11-LABX-0002-01). This study was funded by the EC-supported Network of Excellence Evoltree (GOCE-016322, BEECH initiative), the ERA-Net BiodivERsA, LINKTREE (ANR-08-Biodiversa-006-06), the TipTree (ANR-12-EBID-0003) and MeCC (ANR-13-ADAP-0006) projects.

## Author contributions

JG, EKK and SO-M conceived and designed the study; JG collected the genetic and phenotypic data, performed the analyses, wrote the first draft of the text and contributed the revisions; SO-M conceived and coordinated the beech experiment; JG, HD, OB, EKK and SO-M discussed and interpreted the findings; and EKK, HD, OB and SO-M contributed to the analyses and the revisions of preliminary drafts.

## Data availability statement

Data available upon request from the authors.

## ORCID

Oliver Brendel https://orcid.org/0000-0003-3252-0273

Hendrik Davi https://orcid.org/0000-0001-8828-3145

Julie Gauzere https://orcid.org/0000-0003-4420-6185

Etienne K. Klein https://orcid.org/0000-0003-4677-0775

Sylvie Oddou-Muratorio https://orcid.org/0000-0003-2374-8313

## References

* Aitken SN, Bemmels JB. 2016. Time to get moving: assisted gene flow of forest trees. *Evolutionary Applications* 9: 271–290.

* Aitken SN, Whitlock MC. 2013. Assisted gene flow to facilitate local adaptation to climate change. *Annual Review of Ecology, Evolution and Systematics* 44: 367–388.

* Aitken SN, Yeaman S, Holliday JA, Wang T, Curtis-McLane S. 2008. Adaptation, migration or extirpation: climate change outcomes for tree populations. *Evolutionary Applications* 58: 95–111.

* Alberto FJ, Aitken SN, Alia R, Gonzalez-Martinez SC, Hanninen H, Kremer A, Lefevre F, Lenormand T, Yeaman S, Whetten R et al. 2013. Potential for evolutionary responses to climate change evidence from tree populations. *Global Change Biology* 19: 1645–1661.

* Anderson JT, Perera N, Chowdhury B, Mitchell-Olds T. 2015. Microgeographic patterns of genetic divergence and adaptation across environmental gradients in *Boechera stricta* (Brassicaceae). *American Naturalist* 186: S60–S73.

* Ashley MV. 2010. Plant parentage, pollination, and dispersal: How DNA microsatellites have altered the landscape. *Critical Reviews in Plant Sciences* 29: 148–161.

* Austerlitz F, Dick C, Dutech C, Klein E, Oddou-Muratorio S, Smouse P, Sork V. 2004. Using genetic markers to estimate the pollen dispersal curve. *Molecular Ecology* 13: 937–954.

* Barton N. 2001. Adaptation at the edge of a species range. In: Silvertown J, Antonovics J, eds. *Integrating Ecology and Evolution in a Spatial Context*. New York: Blackwell Science, 365–392.

* Bontemps A, Davi H, Lefevre F, Rozenberg P, Oddou-Muratorio S. 2017. How do functional traits syndromes covary with growth and reproductive performance in a water-stressed population of *Fagus sylvatica*? *Oikos* 126: 1472–1483.

* Bontrager M, Angert AL. 2019. Gene flow improves fitness at a range edge under climate change. *Evolution Letters* 3: 55–68.

* Breda N, Huc R, Granier A, Dreyer E. 2006. Temperate forest trees and stands under severe drought: a review of ecophysiological responses, adaptation processes and long-term consequences. *Annals of Forest Science* 63: 625–644.

* Bridle JR, Polechova J, Kawata M, Butlin RK. 2010. Why is adaptation prevented at ecological margins? New insights from individual-based simulations. *Ecology Letters* 13: 485–494.

* Brousseau L, Bonal D, Cigna J, Scotti I. 2013. Highly local environmental variability promotes intrapopulation divergence of quantitative traits: an example from tropical rain forest trees. *Annals of Botany* 112: 1169–1179.

* Buiteveld J, Vendramin GG, Leonardi S, Kamer K, Geburek T. 2007. Genetic diversity and differentiation in European beech (*Fagus sylvatica* L.) stands varying in management history. *Forest Ecology and Management* 247: 98–106.

* Cailleret M, Davi H. 2011. Effects of climate on diameter growth of co-occurring *Fagus sylvatica* and *Abies alba* along an altitudinal gradient. *Trees* 25: 265–276.

* Charmantier A, Réale D. 2005. How do misassigned paternities affect the estimation of heritability in the wild? *Molecular Ecology* 14: 2839–2850.

* Comps B, Gomory D, Letouzey J, Thiebaut B, Petit R. 2001. Diverging trends between heterozygosity and allelic richness during postglacial colonization in the European beech. *Genetics* 157: 389–397.
Costa e Silva J, Potts BM, Lopez GA. 2014. Heterosis may result in selection favouring the products of long-distance pollen dispersal in Eucalyptus. *PLoS ONE* 9: e93811.

Csillery K, Lalaguee H, Vendramin GG, Gonzalez-Martinez SC, Fady B, **Oddou-Muratorio** S. 2014. Detecting short spatial scale local adaptation and epistatic selection in climate-related candidate genes in European beech (Fagus sylvatica) populations. *Molecular Ecology* 23: 4696-4708.

Csillery K, Ovaskainen O, Sperisen C, Buchmann N, Widmer A, Gugerli F. 2020. Adaptation to local climate in multi-trait space: evidence from silver fir (Abies alba mill.) populations across a heterogeneous environment. *Heredity* 124: 77-92.

Davi H, Gillmann M, Ibanez T, Cailleret M, Bontemps A, Fady B, Lefevre F. 2011. Diversity of leaf unfolding dynamics among tree species: New insights from a study along an altitudinal gradient. *Agricultural and Forest Meteorology* 151: 1504-1513.

Davis M, Shaw R. 2001. Range shifts and adaptive responses to Quaternary climate change. *Science* 292: 673-679.

Eckert AJ, Maloney PE, Vogler DR, Jensen CE, Mix AD, Neale DB. 2015. Local adaptation at fine spatial scales: an example from sugar pine (Pinus lambertiana, Pinaceae). *Tree Genetics & Genomes* 11: 42.

Farquhar G, Ehleringer J, Hubick K. 1989. Carbon isotope discrimination and photosynthesis. *Annual Review of Plant Physiology and Plant Molecular Biology* 40: 503-537.

Fitzpatrick SW, Gerberich JC, Angeloni LM, Bailey LL, Broder ED, Torres-Dowdall J, Handelsman CA, Lopez-Sepulcre A, Reznick DN, Ghalambor CK *et al*. 2016. Gene flow from an adaptively divergent source causes rescue through genetic and demographic factors in two wild populations of Trinidadian guppies. *Evolutionary Applications* 9: 879-891.

Fitzpatrick SW, Gerberich JC, Kronenberger JA, Angeloni LM, Funk WC. 2015. Locally adapted traits maintained in the face of high gene flow. *Ecology Letters* 18: 37-47.

Fitzpatrick SW, Reid BN. 2019. Does gene flow aggravate or alleviate maladaptation to environmental stress in small populations? *Evolutionary Applications* 12: 1402-1416.

Fu YH, Piao S, Vissaste Y, Zhao H, De Boeck HJ, Liu Q, Yang H, Weber U, Hanninen H, Janssens IA. 2015. Increased heat requirement for leaf flushing in temperate woody species over 1980-2012: effects of chilling, precipitation and insolation. *Global Change Biology* 21: 2687-2697.

Gauzere J, Klein E, Davi H, Brendel O, Oddou-Muratorio S. 2016a. Using partial genotyping to estimate the genetic and maternal determinants of adaptive traits in a progeny trial of *Fagus sylvatica*. *Tree Genetics and Genomes* 12: 115.

Gauzere J, Klein EK, Oddou-Muratorio S. 2013a. Ecological determinants of mating system within and between three *Fagus sylvatica* populations along an elevational gradient. *Molecular Ecology* 22: 5001-5015.

Gauzere J, Oddou-Muratorio S, Gay L, Klein E. 2016b. Partial genotyping at polymorphic markers can improve heritability estimates in sibling groups. *Molecular Ecology Resources* 16: 1340-1352.

Gauzere J, Oddou-Muratorio S, Pichot C, Lefevre F, Klein E. 2013b. Biases in quantitative genetic analyses using open-pollinated progeny tests from natural tree populations. *Acta Botanica Gallica* 160: 227-238.

Gauzere J, Teuf B, Davi H, Chevin L.-M., Caignard T, Delzon S, Ronce O, Chuine I. 2020. Where is the optimum? Predicting the variation of selection along climatic gradients and the adaptive value of plasticity. A case study on tree phenology. *Evolution Letters* 4: 109-123.

Gilmour A, Gogel B, Cullis B, Thompson R. 2006. Asreml user guide release 2.0. Hemel Hempstedt, UK: VSN International.

Gomory D, Paule L. 2011. Trade-off between height growth and spring flushing in common beech (Fagus sylvatica L.). *Annals of Forest Science* 68: 975-984.

Hajek P, Kurjak D, vonWuehlisch G, Delzon S, Schuldt B. 2016. Intraspecific variation in wood anatomical, hydraulic, and foliar traits in ten European beech provenances differing in growth yield. *Frontiers in Plant Science* 7: 791.

Hardy O, Vekemans X. 2002. SPAGEDi: a versatile computer program to analyse spatial genetic structure at the individual or population levels. *Molecular Ecology Notes* 2: 618-682, 620.

Hereford J. 2009. A quantitative survey of local adaptation and fitness trade-offs. *American Naturalist* 173: 579-588.

Jump AS, Hunt JM, Martinez-Izquierdo JA, Penuelas J. 2006. Natural selection and climate change: temperature-linked spatial and temporal trends in gene frequency in *Fagus sylvatica*. *Molecular Ecology* 15: 3469-3480.

Karhunen M, Merila J, Leinonen T, Cano JM, Ovaskainen O. 2013. driftsel: an R package for detecting signals of natural selection in quantitative traits. *Molecular Ecology Resources* 690: 746-754.

Karhunen M, Ovaskainen O, Herczeg G, Merila J. 2014. Bringing habitat information into statistical tests of local adaptation in quantitative traits: a case study of nine-spined sticklebacks. *Evolution* 68: 559-568.

Karhunen M, Ovaskainen O. 2012. Estimating population-level coancestry coefficients by an admixture F model. *Genetics* 192: 609.

Kattge J, Knorr W, Raddatz T, Wirth C. 2009. Quantifying photosynthetic capacity and its relationship to leaf nitrogen content for global-scale terrestrial biosphere models. *Global Change Biology* 15: 976-991.

Klinkhamer P, de Jong T, Metz H. 1997. Sex and size in cosexual plants. *Trends in Ecology & Evolution* 12: 260-265.

Knutzen F, Meier IC, Leuschner C. 2015. Does reduced precipitation trigger physiological and morphological drought adaptations in European beech (Fagus sylvatica L.)? Comparing provenances across a precipitation gradient. *Tree Physiology* 35: 949-963.

Kramer K, Ducoussou A, Gomory D, Hansen JK, Ionita L, Liesebach M, Lorent A, Schueler S, Sulkowska M, de Vries S *et al*. 2017. Chilling and forcing requirements for foliage bud burst of European beech (Fagus sylvatica L.) differ between provenances and are phenotypically plastic. *Agricultural and Forest Meteorology* 234: 172-181.

Kremer A, Ronce O, Robledo-Arnuncio JJ, Guillaume F, Bohrer G, Nathan R, Bridle JR, Gomulkiewicz R, Klein EK, Ritland K *et al*. 2012. Long-distance gene flow and adaptation of forest trees to rapid climate change. *Ecology Letters* 15: 378-392.

Kuparinen A, Savolainen O, Schurr FM. 2010. Increased mortality can promote evolutionary adaptation of forest trees to climate change. *Forest Ecology and Management* 259: 1003-1008.

Lamy JB, Bouffier L, Burlett R, Plomion C, Cochard H, Delzon S. 2011. Uniform selection as a primary force reducing population genetic differentiation of cavitation resistance across a species range. *PLoS ONE* 6: e23476.

Lander TA, Oddou-Muratorio S, Prouillet-Leplat H, Klein EK. 2011. Reconstruction of a beech population bottleneck using archival demographic information and bayesian analysis of genetic data. *Molecular Ecology* 20: 5182-5196.

Le Corre V, Kremer A. 2003. Genetic variability at neutral markers, quantitative trait loci and trait in a subdivided population under selection. *Genetics* 164: 1205-1219.

Le Corre V, Kremer A. 2012. The genetic differentiation at quantitative trait loci under local adaptation. *Molecular Ecology* 21: 1548-1566.

Lefevre F, Boivin T, Bontemps A, Courbet F, Davi H, Durand-Gillmann M, Fady B, Gauzere J, Gidoin C, Karam M *et al*. 2014. Considering evolutionary processes in adaptive forestry. *Annals of Forest Science* 71: 723-739.

Leimu R, Fischer M. 2008. A meta-analysis of local adaptation in plants. *PLoS ONE* 3: e4010.

Lenormand T. 2002. Gene flow and the limits to natural selection. *Trends in Ecology and Evolution* 17: 183-189.

Loiselle B, Sork VL, Nason JD, Graham C. 1995. Spatial genetic structure of a tropical understory shrub, *Psychotria officinalis* (rubiaceae). *American Journal of Botany* 82: 1420-1425.

Lopez S, Rousset F, Shaw FH, Shaw RG, Ronce O. 2008. Migration load in plants: role of pollen and seed dispersal in heterogeneous landscapes. *Journal of Evolutionary Biology* 21: 294-309.

MacColl ADC. 2011. The ecological causes of evolution. *Trends in Ecology & Evolution* 26: 514-522.

Marshall T, Slate J, Kruuk L, Pemberton J. 1998. Statistical confidence for likelihood based paternity inference in natural populations. *Molecular Ecology* 7: 639-655.

McKay JK, Latta R. 2002. Adaptive population divergence: markers, QTL and traits. *Trends in Ecology and Evolution* 17: 285-291.

Menzel A, Sparks TH, Estrella N, Koch E, Aasa A, Ahas R, Alm-Kuebler K, Bissolli P, Braslavska O, Briede A *et al*. 2006. European phenological response
to climate change matches the warming pattern. *Global Change Biology* 12: 1969-1976.

Merzeau D, Comps B, Thiébaut B, Letouzey J. 1994. Estimation of fagus sylvatica mating system parameters in natural populations. *Annals of Forest Science* 51: 163-173.

Moody KN, Hunter SN, Childress MJ, Blob RW, Schoenfuss HL, Blum MJ, Ptacek MB. 2015. Local adaptation despite high gene flow in the waterfall-climbing Hawaiian goby, *Sicyopterus stimpsoni*. *Molecular Ecology* 24: 545-563.

Nilsson J. 1995. Genetic-variation in the natural pollen cloud of Pinus-Sylvestris - a study based on progeny testing. *Scandinavian Journal of Forest Research* 10: 140-148.

Nourtier M, Chanzy A, Cailleret M, Yingge X, Huc R, Davi H. 2013. Transpiration of silver fir (Abies alba mill.) during and after drought in relation to soil properties in a Mediterranean mountain area. *Annals of Forest Science* 71: 683-695.

Nussey DH, Wilson AJ, Morris A, Pemberton J, Clutton-Brock T, Kruuk LEB. 2008. Testing for genetic trade-offs between early- and late-life reproduction in a wild red deer population. *Proceedings of the Royal Society of London. Series B: Biological Sciences* 275: 745-750.

Oddou-Muratorio S, Gauzere J, Bontemps A, Rey JF, Klein EK. 2018. Tree, sex and size: Ecological determinants of male vs. female fecundity in three *Fagus sylvatica* stands. *Molecular Ecology* 27: 3131-3145.

Ovaskainen O, Karhunen M, Zheng C, Arias JMC, Merila J. 2011. A new method to uncover signatures of divergent and stabilizing selection in quantitative traits. *Genetics* 189: 621-U729.

Pease C, Lande R, Bull J. 1989. A model of population-growth, dispersal and evolution in a changing environment. *Ecology* 70: 1657-1664.

Peterson ML, Kay KM, Angert AL. 2016. The scale of local adaptation in Mimulus guttatus: comparing life history races, ecotypes, and populations. *New Phytologist* 211: 345-356.

Petit R, Hampe A. 2006. Some evolutionary consequences of being a tree. *Annual Review of Ecology Evolution and Systematics* 37: 187-214.

Piotti A, Leonardi S, Buiteveld J, Geburek T, Gerber S, Kramer K, Vettori C, Vendramin GG. 2011. Comparison of pollen gene flow among four European beech (*Fagus sylvatica* L.) populations characterized by different management regimes. *Heredity* 108: 322-331.

Pluess AR, Frank A, Heiri C, Lalaguee H, Vendramin GG, Oddou-Muratorio S. 2016. Genome-environment association study suggests local adaptation to climate at the regional scale in *Fagus sylvatica*. *New Phytologist* 210: 589-601.

Polechova J, Barton NH. 2015. Limits to adaptation along environmental gradients. *Proceedings of the National Academy of Sciences, USA* 112: 6401-6406.

Rehfeldt G, Tchebakova N, Parfenova Y, Wykoff W, Kuzmina N, Milyutin L. 2002. Intraspecific responses to climate in *Pinus sylvestris*. *Global Change Biology* 8: 912-929.

Rehfeldt G, Ying C, Spittlehouse D, Hamilton D. 1999. Genetic responses to climate in *Pinus contorta*: Niche breadth, climate change, and reforestation. *Ecological Monographs* 69: 375-407.

Reich P, Wright I, Cavender-Bares J, Craine J, Oleksyn J, Westoby M, Walters M. 2003. The evolution of plant functional variation: Traits, spectra, and strategies. *International Journal of Plant Sciences* 164: S143-S164.

Richardson JL, Urban MC, Bolnick DI, Skelly DK. 2014. Microgeographic adaptation and the spatial scale of evolution. *Trends in Ecology & Evolution* 29: 165-176.

Robson TM, Garzon MB, BeechCOSTe52 database consortium Consortium. 2018. Phenotypic trait variation measured on European genetic trials of *Fagus sylvatica* L. *Scientific Data* 5: 180149.

Rose L, Leuschner C, Koeckemann B, Buschmann H. 2009. Are marginal beech (*Fagus sylvatica* L.) provenances a source for drought tolerant ecotypes? *European Journal of Forest Research* 128: 335-343.

Savolainen O, Pyhajarvi T, Knurr T. 2007. Gene flow and local adaptation in trees. *Annual Review of Ecology Evolution and Systematics* 38: 595-619.

Soularue JP, Kremer A. 2012. Assortative mating and gene flow generate clinal phonological variation in trees. *BMC Evolutionary Biology* 12: 79.

Stojnic S, Suchocka M, Benito-Garzon M, Torres-Ruiz JM, Cochard H, Bolte A, Cocozza C, Cvjetkovic B, de Luis M, Martinez-Vilalta J *et al*. 2018. Variation in xylem vulnerability to embolism in European beech from geographically marginal populations. *Tree Physiology* 38: 173-185.

Swindell W, Bouzat J. 2006. Gene flow and adaptive potential in Drosophila melanogaster. *Conservation Genetics* 7: 79-89.

Tansey CJ, Hadfield JD, Phillimore AB. 2017. Estimating the ability of plants to plastically track temperature-mediated shifts in the spring phenological optimum. *Global Change Biology* 23: 3321-3334.

Teissier du Cros E, Thiebaut B, Duval H. 1988. Variability in beech - budding, height growth and tree form. *Annals of Forest Science* 45: 383-398.

Vitasse Y, Signarbieux C, Fu YH. 2018. Global warming leads to more uniform spring phenology across elevations. *Proceedings of the National Academy of Sciences, USA* 115: 1004-1008.

Wang KS. 2003. Mating system in isolated stands of European beech (*Fagus sylvatica* L.). *Forest Genetics* 10: 159-164.

Way DA, Montgomery RA. 2015. Photoperiod constraints on tree phenology, performance and migration in a warming world. *Plant, Cell & Environment* 38: 1725-1736.

Whitlock MC. 2008. Evolutionary inference from Q(ST). *Molecular Ecology* 17: 1885-1896.

Yeaman S, Jarvis A. 2006. Regional heterogeneity and gene flow maintain variance in a quantitative trait within populations of lodgepole pine. *Proceedings of the Royal Society of London. Series B: Biological Sciences* 273: 1587-1593.

Zohner CM, Mo L, Renner SS. 2018. Global warming reduces leaf-out an flowering synchrony among individuals. *eLife* 7: e40214.

## Supporting Information

Additional Supporting Information may be found online in the Supporting Information section at the end of the article.

**Fig. S1** Distribution of the LOD-score of long-distance vs local immigration between 'immigrant' and 'local' groups.

**Table S1** Location, climate and ecological characteristics of the study plots.

**Table S2** Test for plot differentiation and signal of selection in quantitative traits.

**Table S3** Differences in additive genetic variance between 'immigrant' and 'local' groups estimated from model Eqn 4.

**Table S4** Effect of pollen immigration on the variation in additive genetic variance estimated from model Eqn 6.

Please note: Wiley Blackwell are not responsible for the content or functionality of any Supporting Information supplied by the authors. Any queries (other than missing material) should be directed to the *New Phytologist* Central Office.