# Genomic Signatures of Climate-Driven (Mal)Adaptation in an Iconic Conifer, the English Yew (Taxus baccata L.)

Thomas Francisco, Maria Mayol, Elia Vajana, Miquel Riba, Marjana Westergren, Stephen Cavers, Sara Pinosio, Francesca Bagnoli, Maurizio Marchi, Filippos A Aravanopoulos, et al.

► To cite this version:

Thomas Francisco, Maria Mayol, Elia Vajana, Miquel Riba, Marjana Westergren, et al. "Genomic Signatures of Climate-Driven (Mal)Adaptation in an Iconic Conifer, the English Yew (Taxus baccata L.)". Evolutionary Applications, 2025, 18 (10), pp.e70160. 10.1111/eva.70160. hal-05502553

HAL Id: hal-05502553

https://hal.inrae.fr/hal-05502553v1

Submitted on 10 Feb 2026

**HAL** is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published 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 CC BY 4.0 - Attribution - International License
<figure><img src="image_2.png" /></figure>

ORIGINAL ARTICLE

OPEN ACCESS

# Genomic Signatures of Climate-Driven (Mal)Adaptation in an Iconic Conifer, the English Yew (Taxus baccata L.)

Thomas Francisco¹ ID | Maria Mayol² ID | Elia Vajana³ ID | Miquel Riba²,⁴ ID | Marjana Westergren⁵ ID |
Stephen Cavers⁶ ID | Sara Pinosio³ ID | Francesca Bagnoli³ ID | Maurizio Marchi³ ID | Filippos A. Aravanopoulos⁷ ID |
Anna-Maria Farsakoglou⁷ ID | Ivan Scotti⁸ ID | Bruno Fady⁸ ID | Giovanni G. Vendramin³ ID | Juliette Archambeau¹ ID |
Andrea Piotti³ ID | Santiago C. González-Martínez¹ ID

¹National Research Institute for Agriculture, Food and the Environment (INRAE), University of Bordeaux, BIOGECO, Cestas, France | ²Centro de Investigación Ecológica y Aplicaciones Forestales (CREAF), Cerdanyola del Vallès, Catalonia, Spain | ³Institute of Biosciences and BioResources (IBBR), National Research Council (CNR), Sesto Fiorentino, Italy | ⁴Universitat Autònoma Barcelona, Cerdanyola del Vallès, Spain | ⁵Slovenian Forestry Institute, Ljubljana, Slovenia | ⁶UK Centre for Ecology & Hydrology (UKCEH), Edinburgh, UK | ⁷Laboratory of Forest Genetics, Faculty of Agriculture, Forestry and Natural Environment, Aristotle University of Thessaloniki, Thessaloniki, Greece | ⁸National Research Institute for Agriculture, Food and the Environment (INRAE), URFM, Ecology of Mediterranean Forests, Avignon, France

**Correspondence:** Santiago C. González-Martínez (santiago.gonzalez-martinez@inrae.fr)

**Received:** 29 May 2025

**Revised:** 6 September 2025

**Accepted:** 8 September 2025

**Funding:** T.F. was funded by the INRAE Scientific Department ECODIV. E.V. was supported by the MSCA European fellowship MedForAct (GA 101107604). F.A.A. was supported by project Taxol-GR (MIS5004922) co-financed by Greece and the EU. Sampling in Slovenia and M.W. were funded by the Slovenian Research Agency (research core funding no. P4-0107). This research was also funded by EU H2020 GENTREE (grant agreement No. 676876), TAXUS (CGL2007-63107/BOS), ADAPCON (CGL2011-30182-C02-02), READAPT (PID2020-112738GB-I00), and EU H2020 FORGENIUS (grant agreement No. 862221) projects. Views and opinions expressed are, however, those of the authors only and do not necessarily reflect those of the European Union. Neither the European Union nor the granting authority can be held responsible for them.

**Keywords:** climate change | genomic offset | genotype-environment association | local adaptation | *Taxus baccata*

## ABSTRACT

The risk of climate maladaptation is increasing for numerous species, including trees. Developing robust methods to assess population maladaptation remains a critical challenge. Genomic offset approaches aim to predict climate maladaptation by characterizing the genomic changes required for populations to maintain their fitness under changing climates. In this study, we assessed the risk of climate maladaptation in European populations of English yew (*Taxus baccata*), a long-lived tree with a patchy distribution across Europe, the Atlas Mountains, and the Near East, where many populations are small or threatened. We found evidence suggesting local climate adaptation by analyzing 8616 SNPs in 475 trees from 29 European *T. baccata* populations, with climate explaining 18.1% of genetic variance and 100 unlinked climate-associated loci identified via genotype-environment association (GEA). Then, we evaluated the deviation of populations from the overall gene-climate association to assess variability in local adaptation or different adaptation trajectories across populations and found the highest deviations in low latitude populations. Moreover, we predicted genomic offsets and successfully validated these predictions using phenotypic traits assessed in plants from 26 populations grown in a comparative experiment. Finally, we integrated information from current local adaptation, genomic offset, historical genetic differentiation, and effective migration rates to show that Mediterranean and high-elevation *T. baccata* populations face higher vulnerability to climate change than low-elevation Atlantic and continental populations. Our study demonstrates the practical use of the genomic offset framework in conservation genetics, offers insights for its further development, and highlights the need for a population-centered approach that incorporates additional statistics and data sources to credibly assess climate vulnerability in wild plant populations.

This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited.

© 2025 The Author(s). *Evolutionary Applications* published by John Wiley & Sons Ltd.
Genetic adaptation is one of the main mechanisms for long-
term persistence of organisms in a changing environment
(Aitken et al. 2008). For genetic adaptation to occur, individ-
uals must exhibit heritable phenotypic differences, and these
differences must affect their survival and/or reproduction
(i.e., be related to their fitness; Darwin 1859; Mayr 1982).
While evidence for local adaptation in plants (Leimu and
Fischer 2008) and animals (Hereford 2009) is common, it can
only occur or be maintained when the selective pressures are
greater than other evolutionary forces such as genetic drift or
gene flow (Savolainen et al. 2013). Random genetic drift can,
purely by chance, promote mutations that are not beneficial
but rather neutral or even deleterious, potentially under-
mining local adaptation (Yeaman and Otto 2011). Gene flow
can hinder adaptation by introducing maladaptive alleles or
promote it by spreading beneficial alleles across populations
and increasing overall genetic diversity (Joron et al. 2011;
Blanquart et al. 2013; Savolainen et al. 2013; Deacon and
Cavender-Bares 2015).

The predicted intensification of climate change in the near
future poses serious threats to all biodiversity levels (Díaz
et al. 2019). Rapid environmental changes expose species
and populations to novel conditions leading to a phenotypic
mismatch, which may result in a reduction in fitness (Chen
et al. 2011; Brady et al. 2019). A growing body of evidence
shows range contractions and shifts as serious threats to bio-
diversity conservation, particularly affecting sessile organ-
isms that will most likely not be able to adapt and disperse at
the same rate as the predicted climate change (Zhu et al. 2012;
Dyderski et al. 2018; Zu et al. 2021). Additionally, fragmenta-
tion often leads to smaller population sizes, which increase the
likelihood of genetic drift and inbreeding, while limiting gene
flow (Young et al. 1996; Cheptou et al. 2017). Therefore, iso-
lated populations may be more vulnerable to climate change
due to reduced adaptive capacity and migration potential
(Leimu et al. 2010). Thus, the development of reliable meth-
ods to study and predict maladaptation is becoming urgent, in
particular for long-lived, sessile keystone species with highly
fragmented ranges, such as some trees (e.g., the English yew,
*Taxus baccata* L.).

To quantify levels of maladaptation in populations under
changing climate, Fitzpatrick and Keller (2015) proposed an
approach based on genomic data and climate projections,
known as genomic (or genetic) offset. This increasingly pop-
ular approach relies on calculating the magnitude of genetic
change required to maintain the current gene–environment
relationships under novel climates, assuming that populations
are locally adapted (Bernatchez et al. 2024). Several methods
have been used to compute genomic offsets: for example latent
factor mixed models (LFMMs; Gain and François 2021), gra-
dient forest (GF; Fitzpatrick and Keller 2015) or redundancy
analysis (RDA; Capblancq and Forester 2021), and have been
applied to a variety of study systems: crops (Aguirre-Liguori
et al. 2019; Rhoné et al. 2020), trees (Capblancq, Morin,
et al. 2020; Theraroz et al. 2024; Archambeau et al. 2024;
Bonnier et al. 2025), invertebrates (Adam et al. 2022; Bourret
et al. 2024), amphibians (Wu et al. 2025), birds (Bay et al. 2018;

Smith et al. 2021), and, recently, mammals (Hoste et al. 2024;
McLennan et al. 2025).

Recent studies have shown that genomic offset predictions
can be validated using simulation (Láruson et al. 2022) or em-
pirical (Bay et al. 2018; Fitzpatrick et al. 2021; Archambeau
et al. 2024) approaches. However, to date, most validations
have relied on the same data used to train the model. For
this reason, some authors have highlighted the need to use
independent datasets to train and test genomic offset models
in order to avoid model overfitting (e.g., Rellstab et al. 2021;
Lotterhos 2024b). Furthermore, the majority of the evalua-
tion studies using empirical data produced contrasting results
depending on the fitness proxy used (Fitzpatrick et al. 2021;
Lind et al. 2024; Archambeau et al. 2024). For example, in
trees, Lind et al. (2024) found that height was more associ-
ated with genomic offset predictions than mortality for Jack
pine (Pinus banksiana Lamb.) and Douglas fir (Pseudotsuga
menziesii Mirb.), while Archambeau et al. (2024) found the op-
posite pattern for maritime pine (Pinus pinaster Aiton). These
contrasting results suggest that the most appropriate fitness
proxies may vary across species, even for those with similar
life-history traits. Taken together, these findings highlight the
need for more case studies across species and ecological sce-
narios, as well as evaluation data that use a variety of fitness
proxies (e.g., variables related to reproductive, phenological,
or ecophysiological traits), before adopting the genomic offset
approach to support population conservation or management
guidelines.

Although the genomic offset framework is promising and
offers a variety of methods, several assumptions, common
to all methods, need to be carefully considered (see Ahrens
et al. 2023). These include: (i) all populations are equally lo-
cally adapted and (ii) the genotype-environment relationships
can be extrapolated across space and will remain valid in the
future, thus neglecting the adaptive capacity of populations
(Rellstab et al. 2021; Ahrens et al. 2023). Populations devi-
ating from overall gene-climate relationships may be partic-
ularly common in species with fragmented distributions as
some populations would not or only partially exchange genes
with other populations (Deacon and Cavender-Bares 2015). To
date, these key assumptions remain largely untested (but see
the RONA framework; Rellstab et al. 2016). Moreover, a num-
ber of additional processes can exacerbate climate malapda-
tion (e.g., inbreeding depression) or mitigate it (e.g., beneficial
gene flow; Aguirre-Liguori et al. 2021; Rellstab et al. 2021;
Chen et al. 2022; but see Lachmuth et al. 2024). In particular
in trees, the overall effect of gene flow is expected to facili-
tate the evolutionary change required by populations to adapt
to new climatic conditions, as they typically exhibit long-
distance gene flow (Kremer et al. 2012). Several recent stud-
ies have interpreted together (Capblancq, Morin, et al. 2020;
Lazic et al. 2024) or combined (Barratt et al. 2024) various
sources of information with genomic offset to provide more in-
tegrative climate maladaptation or vulnerability predictions,
but the integration of key processes such as gene flow is still
largely missing.

English yew (*Taxus baccata* L.) is a long-lived, dioecious, wind-pollinated and animal-dispersed conifer tree
species, with a fragmented distribution across Europe, the Atlas Mountains, and the Near East, from sea level up to 2,200 m a.s.l. Many populations of *T. baccata* are small and isolated, and they are threatened by fires, browsing, and, in some cases, current or past overexploitation and long-term climate variability (Thomas and Polwart 2003; Mysterud and Østbye 2004; Ruprecht et al. 2010; Benham et al. 2016). The species is found in a wide range of climatic conditions, including high-elevation, Atlantic, continental, and Mediterranean climates (Thomas and Polwart 2003). Cold temperature and drought appear to be among the main limiting factors explaining, respectively, its northern (Thomas and Polwart 2003) and southern (Gegechkori 2018; but see Linares 2013) distribution. Several studies have described overlaying levels of population genetic structure in *T. baccata* at different geographical scales. Indeed, while Mayol et al. (2015) found two to three main gene pools at the rangewide geographical scale, studies at local and regional spatial scales found strong population genetic structure and marked genetic differentiation even between populations in close spatial proximity (e.g., less than 500 m in Dubreuil et al. 2010; see also González-Martínez et al. 2010; Chybicki et al. 2011). Different lines of evidence suggest local adaptation to climate in this species, including phenotypic differences observed under common garden conditions (Mayol et al. 2020), as well as genomic signatures of natural selection in defense-related genes (Burgarella et al. 2012) and climate-related genes (Mayol et al. 2020). Demographic decline is thought to have begun a long time ago in several parts of the *T. baccata* range (e.g., around 300,000 years ago in Iberian populations; Burgarella et al. 2012) and may have intensified more recently (particularly over the past 4,000 years; Thomas and Polwart 2003). Climate change is expected to exacerbate this decline, particularly for Mediterranean populations at the warm and arid range limit, and for those in other sensitive environments, such as high-elevation mountain regions, where climate-induced stress is already being observed and/

or is predicted to be particularly strong in the near future (Mendoza et al. 2009; Knight 2022).

The main objectives of this study are (i) to investigate patterns of local adaptation to climate in *T. baccata*, (ii) to evaluate the potential of genomic offset to predict population maladaptation under novel climatic conditions for a species with a highly fragmented distribution, and (iii) to gain insight into the potential vulnerability of *T. baccata* populations to future climate by linking information on the current degree of local adaptation, the historical capacity of gene flow, and the predicted future climate maladaptation. To meet these objectives, we first used climatic and genomic data from 29 populations across the species' European range to identify gene-climate relationships, as well as candidate climate-associated loci; second, we computed the distance between the observed and predicted genomic composition to identify populations deviating from average gene-climate relationships (subsequently called 'genomic discrepancy index'); and third, we calculated historical effective migration rates and genetic differentiation between populations to determine potential trends in gene flow capacity in the near future, assuming conservatism of historical patterns. Finally, we calculated genomic offsets using two methods and evaluated their predictions using phenotypic traits related to growth, growth phenology, reproductive phenology, and drought/temperature tolerance, measured in plants from 26 populations grown in a comparative experiment.

## 2 | Materials and Methods

### 2.1 | Sampling and Genotyping

Genomic data were obtained for 501 trees from 29 populations covering nearly the entire range of *T. baccata* in Europe (Figure 1 and Table S1). A first dataset involving 25,726 single

<figure><img src="image_3.png" /><figcaption>**FIGURE 1** | Map showing the two sets of *T. baccata* populations used for genomic offset (and related indexes) estimation and validation. In green, the species' European range (EUFORGEN; Caudullo et al. 2017). The black dots represent the 29 SNP-genotyped populations, while the pink triangles indicate the 26 populations from which material was sourced to establish the comparative experiment in Central Spain (clonal bank) used for empirical evaluation of genomic offset models. Population labels correspond to those in Table S1.</figcaption></figure>
Nucleotide Polymorphisms (SNPs) from 120 trees sampled in 12 natural populations was obtained from Mayol et al. (2020). A subset of these SNPs (7,405 SNPs) was used as target SNPs to genotype 11 additional populations (202 individuals) using Single Enrichment Primer Technology (SPET) by Farsakoglou (2025), data available at Data INRAE public repository (doi: 10.15454/7COC1A). In addition, 179 trees from six populations in Slovenia and Italy were genotyped using the same SPET technology. For the latter, DNA was extracted using the NucleoSpin Plant II kit (Macherey Nagel). SNP calling followed GATK best practices (version 3, with -filterExpression 'MQ <40 || MQrank > 12.5 || DP <7.0 || Q <50.0 || QD <1.5 || FS > 60.0 || ReadPosRanksum > -8'; Auwera and der Auwera and O'Connor 2020). Additionally, we filtered for allele balance (between 0.2 and 0.8) using the vcfR R package (Knaus and Grünwald 2017). A subset of three individuals that were genotyped by both Mayol et al. (2020) and Farsakoglou (2025) was used to ensure reproducibility across experiments. Low-quality SNPs and those with inconsistent calls across experiments were discarded, resulting in a final genomic dataset of 11,367 SNPs for 490 trees.

From these data, two datasets were produced for specific analyses (see Table S2). For the first dataset, we filtered for missing data > 15% (for both individuals and SNPs) without filtering for minor allele count (i.e., we retained low-frequency alleles), resulting in 8252 SNPs for 452 individuals. This set was used for population genetic structure analyses, as they are sensitive to high missing data rates and removal of low-frequency SNPs (Linck and Battey 2019). Genotyping errors resulting from the lack of filtering for minor allele count were considered unlikely due to the robust study design and filtering procedures (see above). For the second dataset, used for all other analyses (except for estimating historical gene flow, see below), we filtered for missing data > 30% (for both individuals and SNPs) and removed all SNPs with a minor allele count lower than 20 (i.e., twice the number of individuals in the smallest population before filtering, which corresponded to a minor allele frequency of 2.11%), as low-frequency alleles can lead to high false-positive rates, resulting in 8616 SNPs for 475 individuals (ranging from 8 to 31 individuals per population). Then, we imputed missing data in this dataset by using the most common genotype within each gene pool identified by a STRUCTURE clustering analysis (Pritchard et al. 2000; see below).

Finally, a dataset consisting of 4993 individuals from 238 populations genotyped with seven nuclear microsatellites (nuSSRs) was retrieved from Mayol et al. (2015). This dataset was used to estimate historical gene flow using EEMS (see below and Table S2). Analyses were also performed based only on populations with at least 10 trees sampled (4,658 individuals, 176 populations).

## 2.2 | Climatic Data

Climatic data at high resolution (30 arc-seconds, approximately 1 km at the equator) were obtained from the Climate Downscaling tool (ClimateDT, available at www.ibbr.cnr.it/climate-dt; Marchi et al. 2024). The 1901–1950 time interval, which likely reflects the long-term climatic conditions to which

current *T. baccata* populations are adapted, was defined as the reference period. To predict maladaptation to future climate, we used climatic projections for the 2041–2070 time interval from five global climate models (GCMs; Table S3) under the Shared Socio-economic Pathway SSP3-7.0 (AR6) associated with a severe climate forcing scenario (Fujimori et al. 2017).

To identify climatic variables potentially involved in *T. baccata* adaptation, we applied four selection steps. First, we preselected climatic variables based on their relevance according to previous studies and/or the species' ecology (Melzack and Watts 1982; Moir 1999; Thomas and Polwart 2003; Mayol et al. 2015, 2020; Gegechkori 2018; Sanchez-Martinez et al. 2021; Cedro 2023). Second, we performed forward selection on the preselected variables to identify those explaining most of the genomic variation among populations. For this, we used the ordiR2step function from the vegan R package (Oksanen et al. 2025), ranging from a null model without explanatory variables to a full model with all the explanatory variables. The stopping criterion was the significance of the variable in increasing the model's R², using a *p*-value threshold of 0.05 and 1000 permutations. Twenty runs were performed to ensure robustness (Table S4). Third, we removed climate variables with Pearson's correlation greater than |0.75|, keeping, when possible, the most informative climate variables identified during the forward selection. Fourth, we removed additional climatic variables to achieve variance inflation factors (VIF) below 10, following James et al. (2021). The final six selected climatic variables were: mean annual temperature, mean diurnal range temperature, temperature seasonality, mean temperature of the driest quarter, mean annual precipitation, and precipitation seasonality (Table 1).

To establish how well our SNP data depicted rangewide climatic pressures, we visually compared the climatic envelope of the 29 SNP-genotyped populations with the rangewide one based on all 298 locations with confirmed presence of *T. baccata* used in this study (i.e., including the locations genotyped with nuSSRs in Mayol et al. 2015; Figure S1). The climatic envelopes were built by extracting, for each location, the values for the six climate predictors identified above and combining them using PCA (FactoMineR; Lê et al. 2008). In addition, we compared the climatic envelope of these 29 populations with that of the 26 populations from the comparative experiment used to validate the genomic offset models.

## 2.3 | Phenotypic Traits

Phenotypic data were obtained from plants from 26 natural populations grown in a comparative experiment (a clonal bank established in 1992) under common garden conditions. Only three of the populations (RAS, CAR, and BUJA) in the comparative experiment were also in the SNP dataset but probably sampled different individual trees; therefore, we considered both sets of populations to be nearly independent. The comparative experiment was located in Valsaín, in the foothills of the Sierra de Guadarrama mountains in central Spain (40.9106°N, 4.0125°W, 1,140 m a.s.l). Phenotypic traits were measured several times over the years for over 250 trees (see details in Mayol et al. 2020). We selected growth, growth phenology, reproductive phenology, and drought/temperature tolerance as relevant
<table><caption>**TABLE 1** | Description of the final set of six climatic predictors used in genotype-environment association (GEA) and genomic offset analyses. “OrdiR2step model” corresponds to climatic predictors identified by the forward selection as explaining most of the genomic variation among populations.</caption><thead><tr><th>Variable</th><th>Abbreviation</th><th>Relevance</th></tr></thead><tbody><tr><td>Bio1: Mean Annual Temperature (°C)</td><td>Annual_Tc</td><td>Thomas and Polwart (2003); Gegechkori (2018); Mayol et al. (2020); Ahmadi et al. (2020)</td></tr><tr><td>Bio2: Mean Diurnal Range (°C)</td><td>Diurnal_range_Tc</td><td>OrdiR2step model</td></tr><tr><td>Bio4: Temperature Seasonality (standard deviation × 100)</td><td>Tc_Seasonality</td><td>Cedro and Iszkuło (2011); Ahmadi et al. (2020); OrdiR2step model</td></tr><tr><td>Bio9: Mean Temperature Driest Quarter (°C)</td><td>Tc_Driest_quarter</td><td>Moir (1999); Cedro and Iszkuło (2011); OrdiR2step model</td></tr><tr><td>Bio12: Mean Annual Precipitation (mm)</td><td>Annual_P</td><td>Moir (1999); Mayol et al. (2020)</td></tr><tr><td>Bio15: Precipitation Seasonality (coefficient of variation)</td><td>P_Seasonality</td><td>Ahmadi et al. (2020)</td></tr></tbody></table>

fitness components for climate adaptation, based on previous studies and the ecology of the species, as explained below (see also Table S5 and Mayol et al. 2020).

Growth is commonly considered a proxy for fitness in trees (Fitzpatrick et al. 2021; Capblancq et al. 2023). Individual growth was approximated using shoot volume measured as the volume of the leading current-year shoot from the two longest stems, averaged between the stems and over four years (2009–2012).

Phenological traits often follow latitudinal/altitudinal clines, suggesting climate-related selective pressures that may lead to adaptation (Alberto et al. 2013; Mayol et al. 2020; Silvestro et al. 2023). This is especially relevant for *T. baccata* because some populations occur at high elevation or latitude, where they are exposed to frost events, while others occur at low elevation or latitude, and are exposed to drought during the growing season or benefited from milder climate (coastal Atlantic populations). Growth phenology was estimated from the mean ratio of shoot length before summer (May–June) to the total shoot length at the end of the growing season (October–November), over three consecutive years (2010–2012; referred to as shoot elongation). Reproductive phenology was estimated over the same period by estimating the ratio of open male strobili to the total number of male strobili in several branches during the main period of pollen dispersal (see details in Mayol et al. 2020).

Finally, drought/temperature tolerance traits have garnered increasing interest for evaluating the potential of populations to adapt to warming climates (Keller et al. 2011). Previous studies in *T. baccata* have shown that drought may be one of the limiting factors in the southern part of the species' range (Thomas and Polwart 2003; Gegechkori 2018). We used the mean leaf thickness, an ecophysiological trait associated with drought and thermal tolerance, measured with a digital calliper (to the nearest 0.01 mm) on the two to five largest leaves of the leading shoot at the end of the growing year, over four years (2009–2011, 2021).

For all traits, we calculated best linear unbiased predictors (BLUPs) to estimate population-level phenotypic values,

accounting for covariates influencing trait variation when needed. As the traits showed nearly normal distributions, we used Gaussian mixed-effects models as follows:

$$Y_{ijk} = \mu + pop_{i} + geno_{j(i)} + \varepsilon_{ijk}$$

where $Y_{ijk}$ is the average trait value across years of measurement for individual *k* of genotype *j* in population *i*, $\mu$ is the overall phenotypic mean, $pop_i$ and $geno_{j(i)}$ are the population and genotype (nested in population) random intercepts, respectively, and $\varepsilon_{ijk}$ are the residuals of the model. For the growth model, stem length was used as a covariate (fixed effect) to account for differences in age between trees. BLUPs used in this study correspond to the mean of the posterior distribution of $pop_i$.

All mixed-effects models were fitted within a Bayesian framework using Markov chain Monte Carlo (MCMC) methods implemented in the MCMCglmm R package (Hadfield 2010). Because of the lack of prior knowledge, we used weakly informative priors (inverse-Wishart with $V=1$ and $\nu=0.002$) for all random effects and residuals (Hadfield 2010). Each model ran for 1,500,000 iterations, with a 50,000-iteration burn-in to stabilize the chains and a thinning interval of 500 iterations to reduce autocorrelation. Four independent MCMC chains were run for each model, and chain convergence was assessed using the Gelman-Rubin criteria (Gelman and Rubin 1992) (Table S5).

## 2.4 | Population Genetic Structure and Historical Gene Flow

Population genetic structure analyses were conducted using principal component analysis (PCA), as implemented in the vegan R package (Oksanen et al. 2025), and the STRUCTURE Bayesian clustering approach (version 2.3.4; Pritchard et al. 2000). STRUCTURE models were run for a number of clusters (*K*) ranging from 1 to 10, with 10 independent MCMC runs for each *K*, using 500,000 iterations each, including a burn-in period of 100,000. The most likely number of genetic clusters was estimated from averaged values for each *K*, using the $\Delta K$ method (Evanno et al. 2005).
<table><caption>**TABLE 2** | Description of the six genotype–environment association (GEA) analyses used to detect climate-associated loci (see additional details in Supporting Information).</caption><thead><tr><th>Method</th><th>Functional shape</th><th>Population structure correction</th><th>Procedure for candidate-SNP selection</th><th>Sources</th></tr></thead><tbody><tr><td>Redundancy analysis (RDA)</td><td>Linear</td><td>No</td><td>Extreme Mahalanobis distance with false discovery rate of 10%</td><td>Legendre and Legendre (2012); Capblancq and Forester (2021)</td></tr><tr><td>Partial redundancy analysis (pRDA)</td><td></td><td>Yes, with first two PCs</td><td></td><td></td></tr><tr><td>Latent factor mixed model (LFMM)</td><td>Linear</td><td>Yes, with two latent factors</td><td>False discovery rate of 10%</td><td>Caye et al. (2019)</td></tr><tr><td>BayPass</td><td>Linear</td><td>Yes, with covariance matrix</td><td>Mean Bayes factor across five runs > 8 (substantial evidence on the Jeffreys' scale; Kass and Raftery 1995)</td><td>Gautier (2015)</td></tr><tr><td>Gradient forest (GF-raw)</td><td>Non-linear</td><td>No</td><td>Top 5% overlapping SNPs across five runs</td><td>Ellis et al. (2012); Fitzpatrick et al. (2021)</td></tr><tr><td>Gradient forest corrected (GF-corrected)</td><td></td><td>Yes, with LFMM corrected matrix</td><td></td><td></td></tr></tbody></table>

Historical gene flow patterns were investigated based on the genetic differentiation ($F_{ST}$) of each population relative to a common predicted ancestral gene pool, calculated using BayeScan software version 2.1 (Foll and Gaggiotti 2008), as well as Jost's D (Jost 2008) estimated using the mmod R package (Winter 2012). Potential areas of historical permeability or barriers to gene flow across the range of *T. baccata* were further identified by estimating effective migration surfaces (EEMS, Petkova et al. 2016). This analysis evaluates the consistency of molecular data with expectations under an isolation-by-distance (IBD) model (see details in Table S7 and in Supporting Information).

## 2.5 | Patterns of Local Adaptation to Climate

### 2.5.1 | Genetic Variance Partitioning

To disentangle the influence of demographic history, geography, and climate on population allele frequencies in *T. baccata*, we used a combination of redundancy analysis (RDA) and partial redundancy analysis (pRDA), following Capblancq and Forester (2021). We first ran a full RDA model with all predictors to estimate the total proportion of genetic variance they explain. Second, we ran three pRDA models to estimate the proportion of genetic variance explained separately by demographic history, geography, and climate, while partialing out the effects of the others. Population demographic history was accounted for by using the first two PCs (explaining 24.5% of the variance) of the PCA performed on the SNP data not filtered for minor allele count. As proxies for geography, we used all positive correlation axes from a spatial autocorrelation analysis (distance-based Moran's eigenvector maps, dbMEMs) computed with the adespatial R package (Dray 2016). Finally, the selected climate variables (Table 1) were used as climate predictors. All the RDA and pRDA models, and the PCA, were performed using the vegan R package (Oksanen et al. 2025).

### 2.5.2 | Climate-Related Candidate SNPs

SNPs potentially involved in climate adaptation were identified using six genotype-environment association (GEA) methods: RDA, pRDA, LFMM, BayPass, and Gradient Forest with or without correction for population structure (see Table 2, Figure 2, and additional Methods in Supporting Information). These methods have contrasting assumptions and performance under different scenarios, and combining their results helps reduce false positive rates (Forester et al. 2018).

Candidate SNPs identified by at least two GEA methods were retained in an *outlier set*. This set was filtered to include only unlinked SNPs within contigs, using a linkage disequilibrium threshold of $R^2 < 0.7$, calculated using the genetics R package (Warnes et al. 2022). A set of putatively neutral SNPs (hereafter *random set*) was generated by randomly selecting the same number of SNPs as in the *outlier set* from the SNPs that were not identified as candidates by any of the GEA methods, matching the allele frequencies of the *outlier set*. Finally, we also performed the downstream analyses using all the available SNPs (hereafter the *all set*).

### 2.5.3 | Identifying Populations Deviating From GEA Patterns

We identified populations that deviate from the global GEA trends provided by the RDA model with the best correlation with fitness proxies (i.e., the RDA model based on the *random set* of SNPs; see Figure 3). To achieve this, we extracted the observed and predicted population scores (using the reference time period 1901–1950) based on the linear combination of the climatic variables in the constrained RDA ordination space. These scores were extracted for the most explanatory RDA axes only ($K=2$ in our case) and were weighted by axis importance. Euclidean distances between the observed
<figure><img src="image_4.png" /><figcaption>**FIGURE 2** | UpSet plot for climate-related SNPs by method. The horizontal bars show the number of candidate SNPs identified by each GEA method, while the vertical bars indicate the number of candidate SNPs unique to a single method and those shared across multiple methods, with the corresponding method(s) represented by the dots directly below the vertical bars. Colors refer to whether the method corrects (in green) or not (in orange) for population structure.</figcaption></figure>

and predicted population scores were calculated and subsequently standardized using a min-max normalization to range from zero to one. The resulting index is hereafter referred to as *genomic discrepancy index* (GDI-RDA). Populations with higher GDI-RDA values deviate more from the mean GEA relationship.

**2.6 | Potential Maladaptation Under Future Climate**

**2.6.1 | Genomic Offset Predictions**

The potential maladaptation of *T. baccata* populations to future climate was estimated using the genomic offset approach (Fitzpatrick and Keller 2015). Although several methods exist to compute this index, the general framework is similar. The first step estimates gene-climate relationships using genomic data and climate conditions at the population locations for a reference period (1901–1950 in our study). In the second step, these relationships are used to project the genomic composition of the populations for both the reference and future periods (2041–2070 under SSP3-7.0 in our study). The final step calculates the genetic distance between the reference and future projected genomic composition, that is the genomic offset. Genomic offsets were calculated using two methods,

redundancy analysis (RDA; Capblancq and Forester 2021) and gradient forest (GF; Fitzpatrick and Keller 2015), whose main difference is the shape of the gene-climate relationship (i.e., linear for RDA and linear as well as nonlinear for GF). For both methods, we predicted the genomic offset using the three SNP sets described above (i.e., the *outlier set*, the *random set*, and the *all set*). We also used two additional SNP sets to assess the consistency of our predictions for the *random* and *outlier sets*: a second random SNP set, referred to as *random_2*, and a set including all 935 outliers identified by at least one GEA method, referred to as *all_outlier*. All genomic offset predictions were made using five distinct GCMs (see Table S3), and the mean genomic offset value across GCMs was standardized using min-max normalization. Genomic offset predictions were projected for a relatively near-future interval given the long generation time of *T. baccata* (reported to be as long as 75 years; Thomas and Polwart 2003), to address the assumption that genotype-climate associations based on current genomic data will remain valid for the prediction time period.

**2.6.2 | Evaluation of Genomic Offset Predictions Using Phenotypic Traits**

To assess the relevance of genomic offset predictions, we estimated their Pearson's correlations with the BLUPs of the four
phenotypic traits (growth, growth phenology, reproductive phe- nology, and drought/temperature tolerance) described above. Only populations with at least three individual plants measured for a given trait were included in the analyses. Single traits are often integrated into complex phenotypes, with different compo- nents of lifetime fitness being approximated by different traits, especially in long-lived organisms such as trees (see review in Climent et al. 2024). Therefore, we also correlated genomic off- set predictions with a multitrait index (hereafter Composite trait index), computed as follows:

$$
\text{Composite trait index}_i = \frac{\sum_{m=1}^{n} (\text{Std}(\text{trait}_{im}))}{n}
$$

where *i* is the population, trait<sub>im</sub> denotes the population inter- cept (or BLUPs) of population *i* and trait *m*, *n* is the total num- ber of traits, and Std corresponds to a min–max standardization (applied to each trait to constrain the values between 0 and 1). While simple, given the lack of better information on trait–fit- ness correlation, this index combines all available traits and gives them equal weights.

Genomic offsets for populations from the comparative experi-
ment were calculated using the GEA models built with the 29
SNP-genotyped populations and the Euclidean distance between
the predicted optimal genomic composition at the population's
reference climate of origin (1901–1950 period) and the experi-
mental site's climate from the planting to the last measurement
date (1992–2012 for growth, growth phenology, and reproduc-
tive phenology, and 1992–2021 for drought/temperature toler-
ance). Relationships between trait BLUPs and genomic offset
predictions for each trait and the composite trait index were in-
vestigated using both Pearson's correlation coefficients, as sug-
gested by Lotterhos (2024b), and linear models, as suggested by
Fitzpatrick et al. (2021). For the latter, analyses were conducted
under both a frequentist framework using the stats R package (R
Core Team 2023) and a Bayesian framework using the brms R
package (Bürkner 2021; see Supporting Information for details).
Because of the particular climatic conditions of the compar-
ative experiment, we expected a negative correlation between
genomic offset predictions and phenotypic traits, assuming that
later phenology, higher growth, and larger leaf thickness would
confer higher fitness in this environment.

3 | Results

3.1 | Climatic Variation and Phenotypic Traits

Visual exploration of the climatic envelope of the 29 SNP-
genotyped populations suggests that it adequately covers the cli-
matic niche of *T. baccata*, as determined from the 298 sites with
known presence of the species across its range (Figure S2a). The
climatic envelope of the populations planted in the compara-
tive experiment matched the climatic envelope of the 29 SNP-
genotyped populations (Figure S2b). Finally, the climates of the
source populations of the plants in the comparative experiment
(reference period 1901–1950) were milder than the current cli-
mate of the experimental site (time periods 1992–2012 and
1992–2021), which is characterized by higher temperature, es-
pecially during the dry season, higher temperature seasonality,

and lower annual precipitation (Figure S3). This showed that the
experimental location was suitable for validating genomic offset
predictions associated with climate change.

Additionally, we found that single phenotypic traits were not
highly correlated; therefore, supporting the idea that they con-
tribute to different fitness components (Figure S4).

3.2 | Population Genetic Structure and Historical Gene Flow

Bayesian clustering analysis identified two main gene pools corresponding to the Western and Eastern parts of the species' range, while centrally located populations, such as those from Switzerland and Germany, exhibited admixed genetic compositions. PCA revealed a similar geographic pattern, consistent with an east-west genetic structure (Figures S5-S7). According to ΔK, the best number of clusters was K=2 (Figure S7c), as largely confirmed by the topology of both individuals and populations along the first two PC axes (Figures S5 and S6). The Western gene pool includes the populations from Spain, the United Kingdom, Norway, France, and southern Italy, while the Eastern one includes those from Sweden, Germany, Switzerland, northern Italy, Slovenia, Slovakia, Bosnia-Herzegovina, and Greece. The two population genetic differentiation indices, namely Bayesian's population-specific FST and Jost's D, were highly correlated (Pearson's correlation = 0.93). Therefore, we only used Bayesian's population-specific FST in the following sections. Population-specific FST ranged from 0.05 to 0.29, with most of the highest values found in the Mediterranean area, but with no region consistently exhibiting the highest genetic differentiation (Figure 4a and Table S6). Additionally, EEMS models enabled the identification of multiple populations with higher (e.g., Sainte-Baume in France, Omberg in Sweden, Gabrovo in Bulgaria, or Tosande in Spain) or lower (e.g., Skarzynsko-Kamienna in Poland, Yenice in Turkey, Las Hurdes in Spain, or the Middle Atlas in Morocco) migration than under an isolation by distance (IBD) model (Figure S8). Levels of estimated historical migration did not seem to follow any particular geographical or environmental pattern. In addition, the correlation between population-specific FST and historical migration estimates from EEMS for the 29 SNP-genotyped populations was very low and non-significant (Pearson's correlation = -0.21, p-value = 0.27).

3.3 | Patterns of Local Adaptation to Climate

3.3.1 | Genetic Variance Partitioning

The full RDA model explained 58.6% of the total genetic vari-
ance across the 29 SNP-genotyped populations and provided
strong evidence against the null hypothesis (i.e., no association
between the model predictors and genetic variation; Table 3).
Partial RDA models showed that 30.9% of the genetic variance
was attributed solely to climate, 17.9% to the demographic his-
tory, and 22.2% to geography, while 29% of the genetic vari-
ance remained confounded between the predictors (see also
Figures S5b and S9).
<table><caption>**TABLE 3** | Partition of the total amount of genetic variance showing proportions explained by climate (clim.), demographic history (pop.struc.) and geography (geog.), as well as by the confounding effect among these factors.</caption><thead><tr><th>Partial RDA models</th><th>Proportion of total explained variance (R<sup>2</sup>)</th><th>Proportion of relative explained variance</th><th>p</th></tr></thead><tbody><tr><td>Full model: F~clim. + pop.struc. + geog.</td><td>0.586</td><td>1</td><td>0.001</td></tr><tr><td>Pure climate: F~clim. | (pop.struc. + geog.)</td><td>0.181</td><td>0.309</td><td>0.041</td></tr><tr><td>Pure population structure: F~pop.struc. | (clim. + geog.)</td><td>0.105</td><td>0.179</td><td>0.001</td></tr><tr><td>Pure geography: F~geog. | (clim. + pop.struc.)</td><td>0.130</td><td>0.222</td><td>0.030</td></tr><tr><td>Confounded climate/population structure/geography</td><td>0.170</td><td>0.290</td><td></td></tr><tr><td>Total unexplained</td><td>0.414</td><td></td><td></td></tr></tbody></table>

<figure><img src="image_5.png" /><figcaption>**FIGURE 3** | Pearson's correlations between the genomic offset predictions and the phenotypic traits assessed under common garden conditions. Genomic offset was calculated for each population as the Euclidean distance between the predicted optimal genomic composition in the climate of origin of the population for the reference period (1901–1950) and the predicted optimal genomic composition in the climate of the comparative experiment from the planting to the last measurement dates (1992–2012 for growth, growth phenology, and reproductive phenology, and 1992–2021 for drought tolerance and the composite trait index), for each method and SNP set considered. The y-axis is inverted to display negative correlations in the upper part, reflecting our expectation that genomic offset predictions and fitness components are negatively correlated under the specific climatic conditions of the comparative experiment (see main text and Figure S3).</figcaption></figure>

### 3.3.2 | Climate-Related Candidate SNPs

Genotype-environment association (GEA) methods identified different numbers of climate-related candidate SNPs (Figure 2 and Figures S10–S14 in Supporting Information). RDA and pRDA identified the highest number of candidates (373 and

397 SNPs, respectively), and also had the largest number of candidates in common (57 SNPs). GF uncorrected for population genetic structure (GF-raw) and BayPass identified a similar number of candidate SNPs (110 and 133, respectively), while GF corrected for population genetic structure (GF-corrected) and LFMM identified the smallest number of candidate SNPs (59 and 24, respectively). The number of overlapping candidate SNPs between at least two methods was small compared to the overall number of candidate SNPs (136 vs. 935 SNPs). Furthermore, no candidate SNP was identified by more than three methods. After pruning SNPs for LD within contigs, 100 candidate SNPs were retained as outliers potentially involved in climate adaptation across populations (i.e., the *outlier set*).

### 3.3.3 | Identifying Populations Deviating From GEA Patterns

We found a general pattern of higher GDI-RDA at lower latitudes, supported by a negative Pearson's correlation of −0.51 between GDI-RDA and latitude (p-value=0.005), with Greek, Bosnian-Herzegovinian, Spanish, and Italian populations having, in most cases, higher GDI-RDA than Swedish, Norwegian, British, and Swiss ones (Figure 4b).

## 3.4 | Genomic Offset Predictions

### 3.4.1 | Variability Across Methods and SNP Sets

Variability in genomic offset predictions was observed across GCMs, with RDA predictions showing greater variability than those of GF (Figures S15 and S16). We used the mean genomic offset predictions across the five GCMs for subsequent comparisons and analyses (Table S8). For each method, genomic offset predictions across the main SNP sets were highly correlated (i.e., all Pearson's correlations above 0.90 for RDA and 0.70 for GF; Figure S17) and thus population ranks based on genomic offset were largely consistent (Figure S18). Moreover, the *random_2* SNP set provided predictions consistent with the initial *random set's*, and the *all_outlier* SNP set also provided predictions very similar to the *outlier set's* (all Pearson's correlations above 0.90; Figure S19). In contrast, we observed some variability in genomic offset predictions across methods (albeit substantial correlation remained: lowest Pearson's correlation = 0.30; see Figure S17). For instance, the Pearson's correlation between
<figure><img src="image_6.png" /><figcaption>**FIGURE 4** | Potential maladaptation of *T. baccata* to future climate. (A) Population-specific *F*<sub>ST</sub>, an estimate of historically realised gene flow, on the European-scale environmental zone map of Metzger (2018). The size of the 29 red dots is proportional to their *F*<sub>ST</sub> value. Population labels correspond to those in Table S1. (B) Genomic discrepancy index (GDI-RDA) based on the *random set* of SNPs. A higher GDI-RDA indicates a greater difference between the predicted and observed genomic composition at the reference period. (C) Genomic offset predictions for the RDA (lower half-circle) and GF (upper half-circle) models using the *random set* of SNPs and the mean climate across five general circulation models for the 2041–2070 time interval under the SSP3-7.0. Red stars represent the populations where the predictions of the two genomic offset models differed by more than two classes. In (B, C), GDI-RDA and genomic offset predictions were standardized to the range 0–1 and categorised into five equal-size classes from low (green) to high (red) values. The species’ European range is presented in dark gray (EUFORGEN; Caudullo et al. 2017).</figcaption></figure>

RDA and GF using the *outlier set* was 0.41, resulting in discrepancies in population rankings (e.g., the southern French population of Sainte-Baume (BAU) had the highest genomic offset using RDA, but the second lowest using GF; Figure S18).

### 3.4.2 | Evaluation Using Phenotypic Traits

Evaluation of genomic offset predictions using phenotypic data provided consistent results for almost all genomic offset models. Indeed, according to Pearson's correlations, all genomic offset

predictions were negatively associated with the phenotypic traits (except for the RDA model based on *all set* and reproductive phenology), as expected, with significant *p*-values for several correlations (Figure 3; Figure S20 and Table S9). Several of these associations were also supported by frequentist and Bayesian linear models. For instance, the negative associations observed between RDA models (with *random* and *outlier sets*) and GF models (with *random set*) and the composite trait index were highly unlikely to be obtained by chance (see Tables S10 and S11). Additionally, we found that the negative Pearson's correlations between the genomic offset predictions and the composite
trait index were greater than most correlations for single traits (Figure 3). Finally, we also found that several genomic offset models appeared similarly associated with the composite trait index, and thus no prediction model can be confidently identified as displaying the best genomic offset predictions overall. Therefore, we retained the two models most strongly associated with fitness proxies and displayed both predictions on the same map for each population by organizing predictions into five equal-size classes (Figure 4c). Predictions were considered to be congruent (and therefore robust) only when the two models differed by less than two of these classes.

### 3.4.3 | Geographical Trends

Based on robust predictions (i.e., those supported by both GF and RDA models with a *random set*), we found a negative correlation between genomic offset and latitude (Pearson's correlation of -0.71, *p*-value = 0.0001), with most of the populations that had higher genomic offset being located at low latitudes in Mediterranean environments (VOU, OLY, and CHO in Greece, CAR in Spain, and ANC and VAL in Italy). Additionally, we also found a positive correlation between genomic offset and elevation (Pearson's correlation of 0.65, *p*-value = 0.0005). However, this correlation was substantially reduced when Atlantic low-elevation populations (CED, RWM, and YEW in the UK, and BRA in Norway) were removed (Pearson's correlation of 0.39, *p*-value = 0.08), suggesting that the correlation is mostly driven by the very low genomic offset of Atlantic coastal populations (but notice the high genomic offset of some high-elevated populations: BLI in the Swiss Alps, and BUJA and BUJG in the Spanish Pyrenees). Finally, populations located in continental environments were also predicted to have low genomic offset (e.g., PAT in Germany or HAR in Slovakia). GF and RDA genomic offset models using other SNP sets (i.e., *outlier*, *random_2*, and *all_outlier*) showed rather similar geographical patterns and comparable numbers of populations with inconsistent predictions (except for the models using the *outlier set*; see Figures S21–S23).

## 4 | Discussion

In this study, we found evidence of local adaptation to climate in *T. baccata*, by unveiling substantial genetic variation among populations explained by climate (18.1% of the total genetic variance) and identifying a sizeable number of climate-associated loci using several GEA methods. We found variability in genomic offset predictions across methods (GF and RDA), but not across SNP sets (*all*, *random*, *outlier*). Despite this, we found negative correlations between nearly all the genomic offset models and four phenotypic traits related to fitness (growth, growth phenology, reproductive phenology, and drought tolerance) measured in plants from 26 nearly independent populations planted in a comparative experiment under common garden conditions. By combining information from estimates of historical genetic differentiation and effective migration, local adaptation, and validated future climate maladaptation predictions, we provided valuable and credible insights into population vulnerability to future climate in *T. baccata*. Overall, our results highlight the vulnerability of Mediterranean and high-elevation mountain populations to ongoing climate change.

### 4.1 | Range-Wide Population Genetic Structure and Historical Gene Flow

At the European scale, SNP markers revealed two main gene pools (western and eastern parts of the species' range), with central populations showing admixture, consistent with patterns reported by Mayol et al. (2015) based on nuclear microsatellites at a similar spatial scale. This study also detected a weak but significant isolation-by-distance (IBD) pattern across the entire species distribution. Nevertheless, despite these large-scale patterns, many *T. baccata* populations appear to be strongly isolated. This is supported by the high population-specific $F_{ST}$ values observed in our study, and a large body of literature reporting strong $F_{ST}$ values at local or regional geographical scales (e.g., Dubreuil et al. 2010; González-Martínez et al. 2010; Maroso et al. 2021; Casier et al. 2024; Chybicki et al. 2024). Notably, population-specific $F_{ST}$ and EEMS provided divergent estimates of historical gene flow in many cases, likely reflecting the influence of multiple interacting processes beyond simple IBD (e.g., historical fragmentation, genetic drift, sporadic gene flow, or local adaptation). However, both methods converged in specific cases, for example VAL (Northern Italy) and CHO (Greece), which were consistently inferred to have experienced reduced historical gene flow, while SUE (Spain) and RWM (United Kingdom) were among the populations with signals of historically high connectivity (a pattern previously reported for United Kingdom populations in Gargiulo et al. 2019).

### 4.2 | Evidence Supporting Local Adaptation to Climate

The percentage of genetic variance explained by climate predictors in *T. baccata* (18.1%) was similar to or even higher than the values reported for other trees (*Fagus sylvatica* L.: 17%, Capblancq, Morin, et al. 2020; *Pseudotaxus chienii* W.C. Cheng: 8.2%, Liu et al. 2021; *Picea rubens* Sarg.: 14%, Capblancq et al. 2023). In addition, we identified a substantial number of climate-related candidate SNPs in this species, including two from genes previously reported to be potentially under positive selection for climate in *T. baccata* (Table S12; Mayol et al. 2020). Nevertheless, congruent genomic offset predictions and correlations with phenotypic traits across SNP sets (including random sets) suggest that the signals captured by GEA may result from linkage disequilibrium (LD) with genome-wide loci of small effect, and thus that climate adaptation along the identified climate gradients is polygenic rather than oligogenic in this species. This hypothesis is supported by recent literature on trees, which reports widespread polygenic adaptation for adaptive traits (e.g., De La Torre et al. 2019; de Miguel et al. 2022). Nevertheless, we cannot exclude that the signals captured by GEA could also reflect neutral demographic history rather than local adaptation, leading to inflated false positive rates among identified climate-associated SNPs, and thus potentially biased GEA models and genomic offset estimates (see Holliday et al. 2010; Rellstab et al. 2015).

One of the main assumptions underlying GEA and genomic offset methods is that populations are equally locally adapted to their environments (Rellstab et al. 2015; Ahrens et al. 2023; Lotterhos 2024a). However, numerous studies in diverse species
have shown that this assumption is often not met. For example, Marcora et al. (2021) found higher performance (e.g., germination success, sapling growth, and survival) under site-of-origin conditions in *Polylepis australis* (Bitter), an endemic tree species of central Argentina, but only in marginal populations. There is also evidence of adaptation lags to current climate affecting unequally different parts of a species' range, probably associated with contrasting demographic histories (Browne et al. 2019; Fréjaville et al. 2020). Furthermore, under strong population structure, or complex mechanisms such as conditional neutrality or genetic redundancy, GEA approaches may fail to detect candidate loci for climate adaptation (Rellstab et al. 2015; Lind et al. 2018). To assess variability in local adaptation across populations or potential differences in adaptive trajectories, we calculated a new *genomic discrepancy index* (GDI-RDA, see Methods). Interestingly, in this study, GDI-RDA followed the inverse pattern of historical effective migration detected by EEMS (Pearson's correlation of -0.39, *p*-value=0.03), with higher GDI-RDA associated with lower predicted historical migration rates. This pattern was also supported by population-specific *F*<sub>ST</sub> (Pearson's correlation of 0.52, *p*-value=0.004). This suggests that isolation and/or higher genetic drift may have hindered local adaptation or led to different adaptive trajectories by limiting the influx of beneficial alleles from other populations (e.g., Deacon and Cavender-Bares 2015) or, alternatively, may have exacerbated the effects of genetic drift (Lande 1993). However, inconsistencies between estimates derived from EEMS and population-specific *F*<sub>ST</sub> highlight the need for further research to clarify the role of gene flow in local adaptation in *T. baccata*. Here, we used GDI-RDA to assess the degree of local adaptation by identifying populations deviating from the overall GEA patterns identified by RDA. To test empirically the validity of this index, future work could use reciprocal transplant experiments (Kawecki and Ebert 2004). The genomic discrepancy index approach could also be extended beyond RDA to other GEA methods (e.g., gradient forest, which does not assume a linear relationship) or incorporated directly into the genomic offset framework.

## 4.3 | Toward a More Comprehensive Genomic Offset Framework

We found variability in genomic offset predictions across methods (GF and RDA). However, despite these differences, we still found substantial positive correlations among them, regardless of the SNP set used, and mostly consistent population rankings. These results are reassuring compared to those of Archambeau et al. (2024), who found almost no correlations between RDA and GF predictions in maritime pine (*Pinus pinaster* Aiton), another European conifer, and inconsistent population rankings. Similarly, Lind et al. (2024) reported inconsistent genomic offset predictions between GF and the Risk of Non-Adaptiveness (RONA, Rellstab et al. 2016) for Jack pine (*Pinus banksiana* Lamb.), but not for Douglas-fir (*Pseudotsuga menziesii* Mirb.). These stark differences between studies illustrate a likely species-specific effect on the consistency of genomic offset predictions across methods, possibly related to levels of population genetic structure, shape of climate adaptive clines, and/or the genomic sampling used (e.g., candidate genes vs. random markers). Taken together, these studies highlight that, in the current

framework, the use of genomic offset predictions without empirical validation should be avoided.

Empirical validation in our study was performed by using phenotypic traits measured under common garden conditions. We found that random SNP sets had predictive accuracy similar to that of candidate SNPs identified by GEA. This is in agreement with recent studies on genomic offset prediction based on GF (Fitzpatrick et al. 2021; Lind et al. 2024) and RDA (Archambeau et al. 2024), which have also reported comparable performance between candidate and random marker sets. In our study, these results might be explained by the fact that both the candidate and random SNPs were drawn from coding regions. As a result, even the random SNPs may carry biologically relevant signals due to linkage with functional variants, as also noted by Lind et al. (2024). Moreover, as the phenotypic traits used for validation are expected to be highly polygenic (see, e.g., de Miguel et al. 2022 for growth), predictive models may benefit from diffuse genome-wide LD patterns rather than individual large-effect SNPs.

Almost all the genomic offset models had moderate to substantial negative correlations with the phenotypic traits, consistent with current knowledge on the adaptation of trees to Mediterranean mountains (Liu et al. 2020; Arroyo et al. 2021; Quan et al. 2024). While a negative relationship between genomic offset predictions and fitness proxies has already been demonstrated using common gardens (Rhoné et al. 2020; Fitzpatrick et al. 2021; Lachmuth et al. 2023; Archambeau et al. 2024), to our knowledge, it has never been tested before using fitness proxies from a nearly independent set of populations. Therefore, our study addresses a pending question in the literature (Lotterhos 2024b) by suggesting that the observed negative correlation between genomic offset predictions and fitness under common garden conditions is not the result of model overfitting.

Finally, we found a stronger negative association between genomic offset and the index integrating the different phenotypic traits (i.e., the *composite trait index*) than for any single trait. This finding suggests (1) that phenology and drought/temperature tolerance traits may provide information on fitness components other than growth and, thus, (2) that combining multiple traits into a synthetic index may help to better approximate fitness, especially for long-lived organisms, such as trees, for which direct lifetime fitness cannot be measured (Climent et al. 2024). Nevertheless, although promising, the calculation of the *composite trait index* in this study may not properly account for the complex interaction of phenotypic traits that determines fitness, and further work on methods for phenotypic trait integration in trees is needed. Additionally, the methodology used in this study to approximate population phenotypic means can lead to anti-conservatism (i.e., increased false-positive rates), as estimate uncertainty was neglected (see Hadfield et al. 2010). This issue is particularly important for the *composite trait index*, where error propagation may accumulate across traits. Consequently, the observed correlations involving phenotypic data should be interpreted with caution.

To conclude, we argue that studies assessing population maladaptation to climate in a conservation context should restrict their estimates to the genotyped populations (i.e., to adopt a
more “population-centered” framework), and trusted only when predictions from different methods are in agreement. Indeed, extrapolating genomic offset models to unsampled areas can lead to misleading conclusions in the context of population conservation (see Lind and Lotterhos 2024). We have also shown in this study that conflicting genomic offset predictions can arise across models for some genotyped populations, even when the predictions are linked with fitness proxies using empirical data. Moreover, integrating additional information, such as the degree of local adaptation or connectivity, seems advisable for a more comprehensive assessment of population vulnerability in the face of climate change.

## 4.4 | The Fate of English Yew in the Face of Climate Change

Given the long-term demographic decline in *T. baccata*, due to various factors including long-term climate variability, one may expect climate change to be an additional threat for population survival. By combining information on the potential degree of local adaptation (GDI-RDA), historical genetic differentiation, effective migration surfaces, and potential maladaptation to future climate, we identified several *T. baccata* populations potentially at risk in the face of climate change.

Populations located along the Mediterranean coast and in high-elevation mountain environments (i.e., in Greece, northern and central Italy, the Pyrenees, Mediterranean Spain, and southern Switzerland) mostly share similar features regarding the above-mentioned factors, suggesting a higher vulnerability to future climate. Apart from high genomic offsets, it appears that most of these populations have evolved under strong historical isolation, so their access to beneficial alleles through gene flow may also be compromised in the near future if these patterns are maintained. In addition, most of these populations show substantial GDI-RDA, suggesting that they deviate from average climate adaptation patterns in the species. This may be because of pre-adaptation, distinct adaptive patterns not captured by the statistical framework, or a lower degree of adaptation due to adaptation lags or random genetic drift, which highlights the higher degree of uncertainty in patterns of adaptation to future climate for high GDI-RDA populations.

Overall, our results are in accordance with the current body of literature that shows the considerable impact of climate change on tree species inhabiting Mediterranean and high-elevation mountain regions (Capblancq, Morin, et al. 2020; Calama et al. 2024; Longo-Minnolo et al. 2025). They also suggest that *T. baccata* populations in Atlantic and continental environments may cope better with future climate. Indeed, populations from northern Spain, the United Kingdom, Norway, Germany, and Slovakia have very low genomic offsets, follow overall climate adaptation patterns (i.e., low GDI-RDA), and, in most cases, seem to have evolved under lower historical isolation.

Caution must be exercised while considering these interpretations, as there are several limitations to our analyses. First, population-specific $F_{ST}$ and EEMS provided inconsistent information on the historical gene flow patterns, suggesting a complex demographic history that is not fully accounted for by

either of the two approaches. Additionally, these analyses provided information on historical patterns but not on the current or future gene flow patterns, which remain unknown. At present, many *T. baccata* populations have low population sizes or suffer from habitat fragmentation, and thus may even be more isolated than in the past, as shown in other tree species with currently fragmented distributions (Jump and Peñuelas 2006; Sebbenn 2011). Second, the genomic offset framework relies on several assumptions that may limit its validity (see Rellstab et al. 2021; Ahrens et al. 2023, amongst others). For example, these approaches assume that the genotype-environment relationships will remain constant, thus neglecting the adaptive capacity of populations. We have addressed some of these assumptions by using common garden data to establish a relationship between genomic offset predictions and fitness, and by limiting our predictions to near-future time intervals, making standing genetic variation the most likely source of future adaptation (Barrett and Schlüter 2008). Nevertheless, we must remain aware of other untested assumptions that may impact future climate maladaptation predictions (e.g., overall lack of local adaptation). Furthermore, the insights into future climate vulnerability provided in this study do not take into account other processes that may mitigate the effects of climate change in the short term, such as phenotypic plasticity, or exacerbate its effects, such as inbreeding depression (Capblancq, Fitzpatrick, et al. 2020; Aguirre-Liguori et al. 2021). Extreme climatic events and projected changes in their frequency, duration, and intensity are not accounted for in our models and predictions. As a result, the patterns of climate adaptation and maladaptation presented in this study may differ from actual outcomes, particularly in certain parts of the species' range (e.g., Northern and Mediterranean regions). More broadly, biotic interactions, such as grazing intensification, pest outbreaks, and competition, may also be altered by climate change and have a strong impact on *T. baccata* populations (Thomas and Polwart 2003; Mysterud and Østbye 2004; Ruprecht et al. 2010).

## 5 | Conclusions

Our study demonstrates the complexity of applying the genomic offset framework for the genetic conservation of *T. baccata*, a declining tree species. It also highlights the benefits of combining different sources of information (e.g., genomic and phenotypic data) and analyses (e.g., genomic offsets, local adaptation estimates, and historical gene flow) to gain insight into population vulnerability to climate change. By using a nearly independent set of populations growing under common garden conditions, we provided a robust interpretation of genomic offset predictions and showed how integrating multiple phenotypic traits in simple indices can help to approximate fitness for long-lived organisms such as trees. However, we also found potential variability in the degree of local adaptation across populations and different degrees of historical isolation, which were not fully consistent across estimation methods. Finally, by combining insights into levels of local adaptation and historical effective gene flow with maladaptation predictions, we suggested that European Mediterranean and high-elevation mountain populations of *T. baccata* will suffer the most in the face of climate change, as opposed to Atlantic and continental ones.
Acknowledgments

We wish to thank Thibaut Capblancq for helpful insights and discussion, and Natalija Dovč, Tjaša Baloh, Rok Damjanić, Carlo Urbinati, Alessandro Vitali, Camilla Avanzi, Ilaria Spanu, Celeste Labriola, Daniele Di Santo, Paolo Piovani, Simone Barbarotti, Stefano Zanzucchi, Eleftheria Dalmaris, and the personnel of Reparto Carabinieri Biodiversità Foresta Umbra for field work. Thanks are extended to Dalibor Ballian for kindly providing yew samples from Bosnia-Herzegovina. Data on *T. baccata* distribution were provided by EUFORGEN (https://www.euforgen.org/). T.F. was funded by the INRAE Scientific Department ECODIV. E.V. was supported by the MSCA European fellowship MedForAct (GA 101107604). F.A.A. was supported by project Taxol-GR (MIS5004922) co-financed by Greece and the EU. Sampling in Slovenia and M.W. were funded by the Slovenian Research Agency (research core funding no. P4-0107). This research was also funded by EU H2020 GENTREE (grant agreement No. 676876), TAXUS (CGL2007-63107/BOS), ADAPCON (CGL2011-30182-C02-02), READAPT (PID2020-112738GB-I00), and EU H2020 FORGENIUS (grant agreement No. 862221) projects. Views and opinions expressed are, however, those of the authors only and do not necessarily reflect those of the European Union. Neither the European Union nor the granting authority can be held responsible for them.

Conflicts of Interest

The authors declare no conflicts of interest.

Data Availability Statement

Genomic data are openly available in Data INRAE at DOI: 10.57745/JJ0ZEI. Scripts for genotype-environment association (GEA), and genomic offset computation and validation can be found in a GitHub public repository at https://github.com/Thomas-Francisco/Genomic-signatures-of-mal-adaptation-to-climate-in-English-Yew.

References

Adam, A. A. S., L. Thomas, J. Underwood, J. Gilmour, and Z. T. Richards. 2022. "Population Connectivity and Genetic Offset in the Spawning Coral *Acropora digitifera* in Western Australia." *Molecular Ecology* 31: 3533-3547. https://doi.org/10.1111/mec.16498.

Aguirre-Liguori, J. A., S. Ramirez-Barahona, and B. S. Gaut. 2021. "The Evolutionary Genomics of Species' Responses to Climate Change." *Nature Ecology & Evolution* 5: 1350-1360. https://doi.org/10.1038/s41559-021-01526-9.

Aguirre-Liguori, J. A., S. Ramirez-Barahona, P. Tiffin, and L. E. Eguiarte. 2019. "Climate Change Is Predicted to Disrupt Patterns of Local Adaptation in Wild and Cultivated Maize." *Proceedings of the Biological Sciences* 286: 20190486. https://doi.org/10.1098/rspb.2019.0486.

Ahmadi, K., S. J. Alavi, G. Z. Amiri, S. M. Hosseini, J. M. Serra-Diaz, and J. C. Svenning. 2020. "The Potential Impact of Future Climate on the Distribution of European Yew (*Taxus baccata* L.) in the Hyrcanian Forest Region (Iran)." *International Journal of Biometeorology* 64: 1451-1462. https://doi.org/10.1007/s00484-020-01922-z.

Ahrens, C. W., P. D. Rymer, and A. D. Miller. 2023. "Genetic Offset and Vulnerability Modelling: Misinterpretations of Results and Violations of Evolutionary Principles." *Authorea*. https://doi.org/10.22541/au.16872971.18670759/v1.

Aitken, S. N., S. Yeaman, J. A. Holliday, T. Wang, and S. Curtis-McLane. 2008. "Adaptation, Migration or Extirpation: Climate Change Outcomes for Tree Populations." *Evolutionary Applications* 1: 95-111. https://doi.org/10.1111/j.1752-4571.2007.00013.x.

Alberto, F. J., S. N. Aitken, R. Alia, et al. 2013. "Potential for Evolutionary Responses to Climate Change—Evidence From Tree Populations." *Global Change Biology* 19: 1645-1661. https://doi.org/10.1111/gcb.12181.

Archambeau, J., M. Benito-Garzón, M. De-Miguel, et al. 2024. "Evaluating Genomic Offset Predictions in a Forest Tree with High Population Genetic Structure." *bioRxiv*. https://doi.org/10.1101/2024.05.17.594631.

Arroyo, M. T. K., Í. Tamburrino, P. Pliscoff, V. Robles, M. Colldecarrera, and P. C. Guerrero. 2021. "Flowering Phenology Adjustment and Flower Longevity in a South American Alpine Species." *Plants* 10: 461. https://doi.org/10.3390/plants10030461.

Barratt, C. D., R. E. Onstein, M. L. Pinsky, et al. 2024. "Life on the Edge: A New Toolbox for Population-Level Climate Change Vulnerability Assessments." *Methods in Ecology and Evolution* 15, no. 11: 2038-2058. https://doi.org/10.1111/2041-210x.14429.

Barrett, R. D. H., and D. Schlüter. 2008. "Adaptation From Standing Genetic Variation." *Trends in Ecology & Evolution* 23: 38-44. https://doi.org/10.1016/j.tree.2007.09.008.

Bay, R. A., and R. J. Harrigan, V. L. Underwood, H. L. Gibbs, T. B. Smith, and K. Ruegg. 2018. "Genomic Signals of Selection Predict Climate-Driven Population Declines in a Migratory Bird." *Science* 359: 83-86. https://doi.org/10.1126/science.aan4380.

Benham, S. E., T. Houston Durrant, G. Caudullo, and D. de Rigo. 2016. "*Taxus baccata* in Europe: Distribution, Habitat, Usage and Threats." *In European Atlas of Forest Tree Species*, edited by J. San-Miguel-Ayanz, D. de Rigo, G. Caudullo, T. Houston Durrant, and A. Mauri, e015921. Publ. Off. EU.

Bernatchez, L., A.-L. Ferchaud, C. S. Berger, C. J. Venney, and A. Xuereb. 2024. "Genomics for Monitoring and Understanding Species Responses to Global Climate Change." *Nature Reviews. Genetics* 25: 165-183. https://doi.org/10.1038/s41576-023-00657-y.

Blanquart, F., O. Kaltz, S. L. Nuismer, and S. Gandon. 2013. "A Practical Guide to Measuring Local Adaptation." *Ecology Letters* 16: 1195-1205. https://doi.org/10.1111/ele.12150.

Bonnier, J., E. Sáez Laguna, T. Francisco, et al. 2025. "Wet Season Environments Drive Local Adaptation in the Timber Tree *Dicorynia guianensis* in French Guiana." *Molecular Ecology* 8: e17759. https://doi.org/10.1111/mec.17759.

Bourret, A., C. Leung, G. N. Puncher, et al. 2024. "Diving Into Broad-Scale and High-Resolution Population Genomics to Decipher Drivers of Structure and Climatic Vulnerability in a Marine Invertebrate." *Molecular Ecology* 33: e17448. https://doi.org/10.1111/mec.17448.

Brady, S. P., D. I. Bolnick, A. L. Angert, et al. 2019. "Causes of Maladaptation." *Evolutionary Applications* 12: 1229-1242. https://doi.org/10.1111/eva.12844.

Browne, L., J. W. Wright, S. Fitz-Gibbon, P. F. Gugger, and V. L. Sork. 2019. "Adaptational Lag to Temperature in Valley Oak (*Quercus lo-bata*) Can Be Mitigated by Genome-Informed Assisted Gene Flow." *Proceedings of the National Academy of Sciences* 116: 25179-25185. https://doi.org/10.1073/pnas.1908771116.

Burgarella, C., M. Navascués, M. Zabal-Aguirre, et al. 2012. "Recent Population Decline and Selection Shape Diversity of Taxol-Related Genes." *Molecular Ecology* 21: 3006-3021. https://doi.org/10.1111/j.1365-294X.2012.05532.x.

Bürkner, P.-C. 2021. "Bayesian Item Response Modeling in R With Brms and Stan." *Journal of Statistical Software* 100, no. 1: 1-54. https://doi.org/10.18637/jss.v100.i05.

Calama, R., C. Martínez, J. Gordo, M. Del Río, M. Menéndez-Miguélez, and M. Pardos. 2024. "The Impact of Climate and Management on Recent Mortality in *Pinus pinaster* Resin-Tapped
Forests of Inland Spain." *Forestry* 97: 120-132. https://doi.org/10.1093/forestry/cpad023.

Capblancq, T., M. C. Fitzpatrick, R. A. Bay, M. Exposito-Alonso, and S. R. Keller. 2020. "Genomic Prediction of (Mal)Adaptation Across Current and Future Climatic Landscapes." *Annual Review of Ecology, Evolution, and Systematics* 51: 245-269. https://doi.org/10.1146/annur ev-ecolsys-020720-042553.

Capblancq, T., and B. R. Forester. 2021. "Redundancy Analysis: A Swiss Army Knife for Landscape Genomics." *Methods in Ecology and Evolution* 12: 2298-2309. https://doi.org/10.1111/2041-210X.13722.

Capblancq, T., S. Lachmuth, M. C. Fitzpatrick, and S. R. Keller. 2023. "From Common Gardens to Candidate Genes: Exploring Local Adaptation to Climate in Red Spruce." *New Phytologist* 237: 1590-1605. https://doi.org/10.1111/nph.18465.

Capblancq, T., X. Morin, M. Gueguen, J. Renaud, S. Lobreaux, and E. Bazin. 2020. "Climate-Associated Genetic Variation in *Fagus sylvatica* and Potential Responses to Climate Change in the French Alps." *Journal of Evolutionary Biology* 33: 783-796. https://doi.org/10.1111/jeb.13610.

Casier, M., A. Van Diest, R. Aerts, et al. 2024. "Genetic Diversity and Structure of Endangered Native Yew *Taxus baccata* in Remnant Populations in Belgium." *Forest Ecology and Management* 553: 121633. https://doi.org/10.1016/j.foreco.2023.121633.

Caudullo, G., E. Welk, and J. San-Miguel-Ayanz. 2017. "Chorological Maps for the Main European Woody Species." *Data in Brief* 12: 662-666. https://doi.org/10.1016/j.dib.2017.05.007.

Caye, K., B. Jumentier, J. Lepeule, and O. François. 2019. "LFMM 2: Fast and Accurate Inference of Gene-Environment Associations in Genome-Wide Studies." *Molecular Biology and Evolution* 36: 852-860. https://doi.org/10.1093/molbev/msz008.

Cedro, A. 2023. "Dendrochronology and Dendroclimatology of Yew in Poland." *Dendrochronologia* 78: 126068. https://doi.org/10.1016/j.dendro.2023.126068.

Cedro, A., and G. Iszkuło. 2011. "Do Females Differ From Males of European Yew (*Taxus baccata* L.) in Dendrochronological Analysis?" *Tree-Ring Research* 67: 3-11. https://doi.org/10.3959/2009-9.1.

Chen, I.-C., J. K. Hill, R. Ohlemüller, D. B. Roy, and C. D. Thomas. 2011. "Rapid Range Shifts of Species Associated With High Levels of Climate Warming." *Science* 333: 1024-1026. https://doi.org/10.1126/science.1206432.

Chen, Y., Z. Jiang, P. Fan, et al. 2022. "The Combination of Genomic Offset and Niche Modelling Provides Insights Into Climate Change-Driven Vulnerability." *Nature Communications* 13: 4821. https://doi.org/10.1038/s41467-022-32546-z.

Cheptou, P.-O., A. L. Hargreaves, D. Bonte, and H. Jacquemyn. 2017. "Adaptation to Fragmentation: Evolutionary Dynamics Driven by Human Influences." *Philosophical Transactions of the Royal Society, B: Biological Sciences* 372: 20160037. https://doi.org/10.1098/rstb.2016.0037.

Chybicki, I. J., A. Oleksa, and J. Burczyk. 2011. "Increased Inbreeding and Strong Kinship Structure in *Taxus baccata* Estimated From Both AFLP and SSR Data." *Heredity* 107: 589-600. https://doi.org/10.1038/hdy.2011.51.

Chybicki, I. J., J. J. Robledo-Arnuncio, J. Bodziarczyk, et al. 2024. "Disrupted Connectivity Within a Metapopulation of a Wind-Pollinated Declining Conifer, *Taxus baccata* L." *Forest Ecosystems* 11: 100240. https://doi.org/10.1016/j.fecs.2024.100240.

Climent, J., R. Alía, K. Karkkainen, et al. 2024. "Trade-Offs and Trait Integration in Tree Phenotypes: Consequences for the Sustainable Use of Genetic Resources." *Current Forestry Reports* 10: 196-222. https://doi.org/10.1007/s40725-024-00217-5.

Darwin, C. 1859. "Origin of the Species." In *British Politics and the Environment in the Long Nineteenth Century*. Routledge.

De La Torre, A. R., B. Wilhite, and D. B. Neale. 2019. "Environmental Genome-Wide Association Reveals Climate Adaptation Is Shaped by Subtle to Moderate Allele Frequency Shifts in Loblolly Pine." *Genome Biology and Evolution* 11: 2976-2989. https://doi.org/10.1093/gbe/evz220.

de Miguel, M., I. Rodríguez-Quilón, M. Heuertz, et al. 2022. "Polygenic Adaptation and Negative Selection Across Traits, Years and Environments in a Long-Lived Plant Species (*Pinus pinaster* Ait., Pinaceae)." *Molecular Ecology* 31: 2089-2105. https://doi.org/10.1111/mec.16367.

Deacon, N. J., and J. Cavender-Bares. 2015. "Limited Pollen Dispersal Contributes to Population Genetic Structure but Not Local Adaptation in *Quercus oleoides* Forests of Costa Rica." *PLoS One* 10: e0138783. https://doi.org/10.1371/journal.pone.0138783.

der Auwera, G. A. V., and B. D. O'Connor. 2020. *Genomics in the Cloud: Using Docker, GATK, and WDL in Terra*. O'Reilly Media, Inc.

Díaz, S. M., J. Settele, and E. Brondizio. 2019. *The Global Assessment Report on Biodiversity and Ecosystem Services: Summary for Policy Makers*. Intergovernmental Science-Policy Platform on Biodiversity and Ecosystem Services.

Dray, S. 2016. "adespatial: Multivariate Multiscale Spatial Analysis." 0.3-24.

Dubreuil, M., M. Riba, S. C. González-Martínez, G. G. Vendramin, F. Sebastiani, and M. Mayol. 2010. "Genetic Effects of Chronic Habitat Fragmentation Revisited: Strong Genetic Structure in a Temperate Tree, *Taxus baccata* (Taxaceae), With Great Dispersal Capability." *American Journal of Botany* 97: 303-310. https://doi.org/10.3732/ajb.0900148.

Dyderski, M. K., S. Paź, L. E. Frelich, and A. M. Jagodziński. 2018. "How Much Does Climate Change Threaten European Forest Tree Species Distributions?" *Global Change Biology* 24: 1150-1163. https://doi.org/10.1111/gcb.13925.

Ellis, N., S. J. Smith, and C. R. Pitcher. 2012. "Gradient Forests: Calculating Importance Gradients on Physical Predictors." *Ecology* 93: 156-168. https://doi.org/10.1890/11-0252.1.

Evanno, G., S. Regnaut, and J. Goudet. 2005. "Detecting the Number of Clusters of Individuals Using the Software STRUCTURE: A Simulation Study." *Molecular Ecology* 14: 2611-2620. https://doi.org/10.1111/j.1365-294X.2005.02553.x.

Farsakoglou, A.-M. 2025. "Using Genomic Approaches for the Conservation of Conifer Species Populations at Different Spatial Scales." PhD Dissertation. https://doi.org/10.12681/eadd/59169.

Fitzpatrick, M. C., V. E. Chhatre, R. Y. Soolanayakanahally, and S. R. Keller. 2021. "Experimental Support for Genomic Prediction of Climate Maladaptation Using the Machine Learning Approach Gradient Forests." *Molecular Ecology Resources* 21: 2749-2765. https://doi.org/10.1111/1755-0998.13374.

Fitzpatrick, M. C., and S. R. Keller. 2015. "Ecological Genomics Meets Community-Level Modelling of Biodiversity: Mapping the Genomic Landscape of Current and Future Environmental Adaptation." *Ecology Letters* 18: 1-16. https://doi.org/10.1111/ele.12376.

Foll, M., and O. Gaggiotti. 2008. "A Genome-Scan Method to Identify Selected Loci Appropriate for Both Dominant and Codominant Markers: A Bayesian Perspective." *Genetics* 180: 977-993. https://doi.org/10.1534/genetics.108.092221.

Forester, B. R., J. R. Lasky, H. H. Wagner, and D. L. Urban. 2018. "Comparing Methods for Detecting Multilocus Adaptation With Multivariate Genotype-Environment Associations." *Molecular Ecology* 27: 2215-2233. https://doi.org/10.1111/mec.14584.

Fréjaville, T., N. Vizcaíno-Palomar, B. Fady, A. Kremer, and M. Benito Garzón. 2020. "Range Margin Populations Show High Climate Adaptation Lags in European Trees." *Global Change Biology* 26: 484-495. https://doi.org/10.1111/gcb.14881.
Fujimori, S., T. Hasegawa, T. Masui, et al. 2017. "SSP3: AIM Implementation of Shared Socioeconomic Pathways." *Global Environmental Change* 42: 268-283. https://doi.org/10.1016/j.gloen vcha.2016.06.009.

Gain, C., and O. François. 2021. "LEA 3: Factor Models in Population Genetics and Ecological Genomics With R." *Molecular Ecology Resources* 21: 2738-2748. https://doi.org/10.1111/1755-0998.13366.

Gargiulo, R., M. Saubin, G. Rizzuto, et al. 2019. "Genetic Diversity in British Populations of *Taxus baccata* L.: Is the Seedbank Collection Representative of the Genetic Variation in the Wild?" *Biological Conservation* 233: 289-297. https://doi.org/10.1016/j.biocon.2019.01.014.

Gautier, M. 2015. "Genome-Wide Scan for Adaptive Divergence and Association With Population-Specific Covariates." *Genetics* 201: 1555-1579. https://doi.org/10.1534/genetics.115.181453.

Gegechkori, A. 2018. "Patterns of Distribution and Survival of European Yew (*Taxus baccata* L.) in an Alpine Tree Line Ecotone in the Greater Caucasus (Georgia)." *Annals of Agrarian Science* 16: 170-176. https://doi.org/10.1016/j.aasci.2018.02.006.

Gelman, A., and D. B. Rubin. 1992. "Inference From Iterative Simulation Using Multiple Sequences." *Statistical Science* 7: 457-472. https://doi.org/10.1214/ss/1177011136.

González-Martínez, S. C., M. Dubreuil, M. Riba, G. G. Vendramin, F. Sebastiani, and M. Mayol. 2010. "Spatial Genetic Structure of *Taxus baccata* L. in the Western Mediterranean Basin: Past and Present Limits to Gene Movement Over a Broad Geographic Scale." *Molecular Phylogenetics and Evolution* 55: 805-815. https://doi.org/10.1016/j.ympev.2010.03.001.

Hadfield, J. D. 2010. "MCMCglmm: MCMC Generalised Linear Mixed Models, Version 2.36."

Hadfield, J. D., A. J. Wilson, D. Garant, B. C. Sheldon, and L. E. B. Kruuk. 2010. "The Misuse of BLUP in Ecology and Evolution." *American Naturalist* 175: 116-125. https://doi.org/10.1086/648604.

Hereford, J. 2009. "A Quantitative Survey of Local Adaptation and Fitness Trade-Offs." *American Naturalist* 173: 579-588. https://doi.org/10.1086/597611.

Holliday, J. A., K. Ritland, and S. N. Aitken. 2010. "Widespread, Ecologically Relevant Genetic Markers Developed From Association Mapping of Climate-Related Traits in Sitka Spruce (*Picea sitchensis*)." *New Phytologist* 188: 50-1514. https://doi.org/10.1111/j.1469-8137.2010.03380.x.

Hoste, A., T. Capblancq, T. Broquet, et al. 2024. "Projection of Current and Future Distribution of Adaptive Genetic Units in an Alpine Ungulate." *Heredity* 132: 54-66. https://doi.org/10.1038/s41437-023-00661-2.

James, G., D. Witten, T. Hastie, and R. Tibshirani. 2021. "An Introduction to Statistical Learning." In *Springer Texts in Statistics*. Springer. https://doi.org/10.1007/978-1-0716-1418-1.

Joron, M., L. Frezal, R. T. Jones, et al. 2011. "Chromosomal Rearrangements Maintain a Polymorphic Supergene Controlling Butterfly Mimicry." *Nature* 477: 203-206. https://doi.org/10.1038/nature10341.

Jost, L. 2008. "G(ST) and Its Relatives Do Not Measure Differentiation." *Molecular Ecology* 17: 4015-4026. https://doi.org/10.1111/j.1365-294X.2008.03887.x.

Jump, A. S., and J. Peñuelas. 2006. "Genetic Effects of Chronic Habitat Fragmentation in a Wind-Pollinated Tree." *Proceedings of the National Academy of Sciences of the United States of America* 103: 8096-8100. https://doi.org/10.1073/pnas.0510127103.

Kass, R. E., and A. E. Raftery. 1995. "Bayes Factors." *Journal of the American Statistical Association* 90: 773-795. https://doi.org/10.1080/01621459.1995.10476572.

Kawecki, T. J., and D. Ebert. 2004. "Conceptual Issues in Local Adaptation." *Ecology Letters* 7: 1225-1241. https://doi.org/10.1111/j.1461-0248.2004.00684.x.

Keller, S. R., R. Y. Soolanayakanahally, R. D. Guy, S. N. Silim, M. S. Olson, and P. Tiffin. 2011. "Climate-Driven Local Adaptation of Ecophysiology and Phenology in Balsam Poplar, *Populus balsamifera* L. (Salicaceae)." *American Journal of Botany* 98: 99-108. https://doi.org/10.3732/ajb.1000317.

Knaus, B. J., and N. J. Grünwald. 2017. "vcfR: A Package to Manipulate and Visualize Variant Call Format Data in R." *Molecular Ecology Resources* 17: 44-53. https://doi.org/10.1111/1755-0998.12549.

Knight, J. 2022. "Scientists' Warning of the Impacts of Climate Change on Mountains." *PeerJ* 10: e14253. https://doi.org/10.7717/peerj.14253.

Kremer, A., O. Ronce, J. J. Robledo-Arnuncio, et al. 2012. "Long-Distance Gene Flow and Adaptation of Forest Trees to Rapid Climate Change." *Ecology Letters* 15: 378-392. https://doi.org/10.1111/j.1461-0248.2012.01746.

Lachmuth, S., T. Capblancq, S. R. Keller, and M. C. Fitzpatrick. 2023. "Assessing Uncertainty in Genomic Offset Forecasts From Landscape Genomic Models (And Implications for Restoration and Assisted Migration)." *Frontiers in Ecology and Evolution* 11: 1155783. https://doi.org/10.3389/fevo.2023.1155783.

Lachmuth, S., T. Capblancq, A. Prakash, S. R. Keller, and M. C. Fitzpatrick. 2024. "Novel Genomic Offset Metrics Integrate Local Adaptation Into Habitat Suitability Forecasts and Inform Assisted Migration." *Ecological Monographs* 94: e1593. https://doi.org/10.1002/ecm.1593.

Lande, R. 1993. "Risks of Population Extinction From Demographic and Environmental Stochasticity and Random Catastrophes." *American Naturalist* 142: 911-927. https://doi.org/10.1086/285580.

Láruson, Á. J., M. C. Fitzpatrick, S. R. Keller, B. C. Haller, and K. E. Lotterhos. 2022. "Seeing the Forest for the Trees: Assessing Genetic Offset Predictions From Gradient Forest." *Evolutionary Applications* 15: 403-416. https://doi.org/10.1111/eva.13354.

Lazic, D., C. Geßner, K. J. Liepe, et al. 2024. "Genomic Variation of European Beech Reveals Signals of Local Adaptation Despite High Levels of Phenotypic Plasticity." *Nature Communications* 15: 8553. https://doi.org/10.1038/s41467-024-52933-y.

Lê, S., J. Josse, and F. Husson. 2008. "FactoMineR: An R Package for Multivariate Analysis." *Journal of Statistical Software* 25: 1-18. https://doi.org/10.18637/jss.v025.i01.

Legendre, P., and L. Legendre. 2012. "Chapter 11—Canonical Analysis." In *Developments in Environmental Modelling*, edited by P. Legendre and L. Legendre, 625-710. Elsevier.

Leimu, R., and M. Fischer. 2008. "A Meta-Analysis of Local Adaptation in Plants." *PLoS One* 3: e4010. https://doi.org/10.1371/journal.pone.0004010.

Leimu, R., P. Vergeer, F. Angeloni, and N. J. Ouborg. 2010. "Habitat Fragmentation, Climate Change, and Inbreeding in Plants." *Annals of the New York Academy of Sciences* 1195: 84-98. https://doi.org/10.1111/j.1749-6632.2010.05450.x.

Linares, J. C. 2013. "Shifting Limiting Factors for Population Dynamics and Conservation Status of the Endangered English Yew (*Taxus baccata* L., Taxaceae)." *Forest Ecology and Management* 291: 119-127. https://doi.org/10.1016/j.foreco.2012.11.009.

Linck, E., and C. J. Battey. 2019. "Minor Allele Frequency Thresholds Strongly Affect Population Structure Inference With Genomic Data Sets." *Molecular Ecology Resources* 19: 639-647. https://doi.org/10.1111/1755-0998.12995.

Lind, B. M., R. Candido-Ribeiro, P. Singh, et al. 2024. "How Useful Is Genomic Data for Predicting Maladaptation to Future Climate?" *Global Change Biology* 30: e17227. https://doi.org/10.1111/gcb.17227.

Lind, B. M., and K. E. Lotterhos. 2024. "The Accuracy of Predicting Maladaptation to New Environments With Genomic Data." *Molecular Ecology Resources* 25: 14008. https://doi.org/10.1111/1755-0998.14008.
Lind, B. M., M. Menon, C. E. Bolte, T. M. Faske, and A. J. Eckert. 2018. "The Genomics of Local Adaptation in Trees: Are We out of the Woods Yet?" *Tree Genetics & Genomes* 14: 29. https://doi.org/10.1007/s11295-017-1224-y.

Liu, L., Z. Wang, Y. Su, and T. Wang. 2021. "Population Transcriptomic Sequencing Reveals Allopatric Divergence and Local Adaptation in *Pseudotaxus chienii* (Taxaceae)." *BMC Genomics* 22: 388. https://doi.org/10.1186/s12864-021-07682-3.

Liu, W., L. Zheng, and D. Qi. 2020. "Variation in Leaf Traits at Different Altitudes Reflects the Adaptive Strategy of Plants to Environmental Changes." *Ecology and Evolution* 10: 8166-8175. https://doi.org/10.1002/ece3.6519.

Longo-Minnolo, G., S. Consoli, and M. Tessitori. 2025. "Evaluation of the Impact of Climate Change on *Fagus sylvatica* Dieback-A Combined Approach With ERA5-Land Data and Landsat Imagery." *Remote Sensing* 17: 873. https://doi.org/10.3390/rs17050873.

Lotterhos, K. E. 2024a. "Interpretation Issues With "Genomic Vulnerability" Arise From Conceptual Issues in Local Adaptation and Maladaptation." *Evolution Letters* 8: 331-339. https://doi.org/10.1093/evlett/qrae004.

Lotterhos, K. E. 2024b. "Principles in Experimental Design for Evaluating Genomic Forecasts." *Methods in Ecology and Evolution* 15: 1466-1482. https://doi.org/10.1111/2041-210X.14379.

Marchi, M., G. Bucci, P. Iovieno, and D. Ray. 2024. "ClimateDT: A Global Scale-Free Dynamic Downscaling Portal for Historic and Future Climate Data." *Environments* 11: 82. https://doi.org/10.3390/envir onments11040082.

Marcora, P. I., P. A. Tecco, M. C. Ferrero, et al. 2021. "Are Populations of *Polylepis australis* Locally Adapted Along Their Elevation Gradient?" *Neotropical Biodiversity* 7: 246-256. https://doi.org/10.1080/23766808. 2021.1940049.

Maroso, F., M. Vera, J. Ferreiro, et al. 2021. "Genetic Diversity and Structure of *Taxus baccata* From the Cantabrian-Atlantic Area in Northern Spain: A Guide for Conservation and Management Actions." *Forest Ecology and Management* 482: 118844. https://doi.org/10.1016/j. foreco.2020.118844.

Mayol, M., M. Riba, S. Cavers, et al. 2020. "A Multiscale Approach to Detect Selection in Nonmodel Tree Species: Widespread Adaptation Despite Population Decline in *Taxus baccata* L." *Evolutionary Applications* 13: 143-160. https://doi.org/10.1111/eva.12838.

Mayol, M., M. Riba, S. C. González-Martínez, et al. 2015. "Adapting Through Glacial Cycles: Insights From a Long-Lived Tree (*Taxus baccata*)." *New Phytologist* 208: 973-986. https://doi.org/10.1111/nph. 13496.

Mayr, E. 1982. *The Growth of Biological Thought: Diversity, Evolution, and Inheritance*. Belknap Press of Harvard Univ. Press.

McLennan, E. A., T. G. L. Kovacs, L. W. Silver, et al. 2025. "Genomics Identifies Koala Populations at Risk Across Eastern Australia." *Ecological Applications* 35: e3062. https://doi.org/10.1002/eap.3062.

Melzack, R. N., and D. Watts. 1982. "Cold Hardiness in the Yew (*Taxus baccata* L.) in Britain." *Journal of Biogeography* 9: 231-241. https://doi.org/10.2307/2844666.

Mendoza, I., L. Gómez-Aparicio, R. Zamora, and L. Matías. 2009. "Recruitment Limitation of Forest Communities in a Degraded Mediterranean Landscape." *Journal of Vegetation Science* 20: 367-376. https://doi.org/10.1111/j.1654-1103.2009.05705.x.

Metzger, M. J. 2018. *The Environmental Stratification of Europe*. University of Edinburgh. https://doi.org/10.7488/ds/2356.

Moir, A. K. 1999. "The Dendrochronological Potential of Modern Yew (*Taxus baccata*) With Special Reference to Yew From Hampton Court

Palace, UK.” *New Phytologist* 144: 479–488. https://doi.org/10.1046/j. 1469-8137.1999.00545.x.

Mysterud, A., and E. Østbye. 2004. "Roe Deer (*Capreolus capreolus*) Browsing Pressure Affects Yew (*Taxus baccata*) Recruitment Within Nature Reserves in Norway." *Biological Conservation* 120: 545–548. https://doi.org/10.1016/j.biocon.2004.03.027.

Oksanen, J., G. L. Simpson, F. G. Blanchet, et al. 2025. "vegan: Community Ecology Package." 2:6-8.

Petkova, D., J. Novembre, and M. Stephens. 2016. "Visualizing Spatial Population Structure With Estimated Effective Migration Surfaces." *Nature Genetics* 48: 94–100. https://doi.org/10.1038/ng.3464.

Pritchard, J. K., M. Stephens, and P. Donnelly. 2000. "Inference of Population Structure Using Multilocus Genotype Data." *Genetics* 155: 945-959. https://doi.org/10.1093/genetics/155.2.945.

Quan, Q., N. He, R. Zhang, et al. 2024. "Plant Height as an Indicator for Alpine Carbon Sequestration and Ecosystem Response to Warming." *Nature Plants* 10: 890–900. https://doi.org/10.1038/s4147 7-024-01705-z.

R Core Team. 2023. *R: A Language and Environment for Statistical Computing*. R Foundation for Statistical Computing.

Rellstab, C., B. Dauphin, and M. Exposito-Alonso. 2021. "Prospects and Limitations of Genomic Offset in Conservation Management." *Evolutionary Applications* 14: 1202-1212. https://doi.org/10.1111/eva. 13205.

Rellstab, C., F. Gugerli, A. J. Eckert, A. M. Hancock, and R. Holderegger. 2015. "A Practical Guide to Environmental Association Analysis in Landscape Genomics." *Molecular Ecology* 24: 4348-4370. https://doi.org/10.1111/mec.13322.

Rellstab, C., S. Zoller, L. Walthert, et al. 2016. "Signatures of Local Adaptation in Candidate Genes of Oaks (*Quercus* spp.) With Respect to Present and Future Climatic Conditions." *Molecular Ecology* 25: 5907- 5924. https://doi.org/10.1111/mec.13889.

Rhoné, B., D. Defrance, C. Berthouly-Salazar, et al. 2020. "Pearl Millet Genomic Vulnerability to Climate Change in West Africa Highlights the Need for Regional Collaboration." *Nature Communications* 11: 5274. https://doi.org/10.1038/s41467-020-19066-4.

Ruprecht, H., A. Dhar, B. Aigner, G. Oitzinger, R. Klumpp, and H. Vacik. 2010. "Structural Diversity of English Yew (*Taxus baccata* L.) Populations." *European Journal of Forest Research* 129: 189-198. https:// doi.org/10.1007/s10342-009-0312-4.

Sanchez-Martinez, P., A. Marcer, M. Mayol, and M. Riba. 2021. "Shaping the Niche of *Taxus baccata*, a Modelling Exercise Using Biologically Meaningful Information." *Forest Ecology and Management* 501: 119688. https://doi.org/10.1016/j.foreco.2021. 119688.

Savolainen, O., M. Lascoux, and J. Merilä. 2013. "Ecological Genomics of Local Adaptation." *Nature Reviews. Genetics* 14: 807-820. https://doi.org/10.1038/nrg3522.

Sebbenn, A. 2011. "Effects of Forest Fragmentation on the Effective and Realized Gene Flow of Neotropical Tree Species: Implications for Genetic Conservation." *BMC Proceedings* 5: O6. https://doi.org/10.1186/ 1753-6561-5-S7-O6.

Silvestro, R., C. Mura, D. Alano Bonacini, et al. 2023. "Local Adaptation Shapes Functional Traits and Resource Allocation in Black Spruce." *Scientific Reports* 13: 21257. https://doi.org/10.1038/s41598-023- 48530-6.

Smith, T. B., T. L. Fuller, Y. Zhen, et al. 2021. "Genomic Vulnerability and Socio-Economic Threats Under Climate Change in an African Rainforest Bird." *Evolutionary Applications* 14: 1239-1247. https://doi.org/10.1111/eva.13193.
Theraroz, A., C. Guadaño-Peyrot, J. Archambeau, et al. 2024. "The Genetic Consequences of Population Marginality: A Case Study in Maritime Pine." *Diversity and Distributions* 30: e13910. https://doi.org/10.1111/ddi.13910.

Thomas, P. A., and A. Polwart. 2003. "*Taxus baccata* L." *Journal of Ecology* 91: 489–524. https://doi.org/10.1046/j.1365-2745.2003.00783.x.

Warnes, G., G. Gorjanc, F. Leisch, and M. Man. 2022. "genetics: Population Genetics (Version 1.3.8.1.3). R Package." https://doi.org/10.32614/CRAN/package.genetics.

Winter, D. J. 2012. "Mmod: An R Library for the Calculation of Population Differentiation Statistics." *Molecular Ecology Resources* 12: 1158–1160. https://doi.org/10.1111/j.1755-0998.2012.03174.x.

Wu, N., Q. Xiao, Z. Liao, et al. 2025. "Landscape Genomics Provides Insights Into Climate Change-Driven Vulnerability in Torrent Frogs (Ranidae: Amolops).." *Molecular Ecology* 34: e17807. https://doi.org/10.1111/mec.17807.

Yeaman, S., and S. P. Otto. 2011. "Establishment and Maintenance of Adaptive Genetic Divergence Under Migration, Selection, and Drift." *Evolution* 65: 2123–2129. https://doi.org/10.1111/j.1558-5646.2011.01277.x.

Young, A., T. Boyle, and T. Brown. 1996. "The Population Genetic Consequences of Habitat Fragmentation for Plants." *Trends in Ecology & Evolution* 11: 413–418. https://doi.org/10.1016/0169-5347(96)10045-8.

Zhu, K., C. W. Woodall, and J. S. Clark. 2012. "Failure to Migrate: Lack of Tree Range Expansion in Response to Climate Change." *Global Change Biology* 18: 1042–1052. https://doi.org/10.1111/j.1365-2486.2011.02571.x.

Zu, K., Z. Wang, X. Zhu, et al. 2021. "Upward Shift and Elevational Range Contractions of Subtropical Mountain Plants in Response to Climate Change." *Science of the Total Environment* 783: 146896. https://doi.org/10.1016/j.scitotenv.2021.146896.

## Supporting Information

Additional supporting information can be found online in the Supporting Information section. **Data S1:** eva70160-sup-0001-DataS1.docx.