# Linking ecological niche models and common garden experiments to predict phenotypic differentiation in stressful environments: Assessing the adaptive value of marginal populations in an alpine plant

Javier Morente-López, Jamie Kass, Carlos Lara-Romero, Josep Serra-Diaz, José Carmen Soto-Correa, Robert Anderson, José Iriondo

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

## To cite this version:

Javier Morente-López, Jamie Kass, Carlos Lara-Romero, Josep Serra-Diaz, José Carmen Soto-Correa, et al.. Linking ecological niche models and common garden experiments to predict phenotypic differentiation in stressful environments: Assessing the adaptive value of marginal populations in an alpine plant. Global Change Biology, 2022, 28 (13), pp.4143-4162. 10.1111/gcb.16181. hal-04066774

HAL Id: hal-04066774

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

Submitted on 12 Apr 2023

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
## RESEARCH ARTICLE

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

WILEY

# Linking ecological niche models and common garden experiments to predict phenotypic differentiation in stressful environments: Assessing the adaptive value of marginal populations in an alpine plant

Javier Morente-López<fnref n="1" />,<fnref n="2" /> ID

Jamie M. Kass<fnref n="3" />,<fnref n="4" />,<fnref n="5" /> ID

Carlos Lara-Romero<fnref n="1" /> ID

Josep M. Serra-Diaz<fnref n="6" /> ID

José Carmen Soto-Correa<fnref n="7" />

Robert P. Anderson<fnref n="3" />,<fnref n="4" />,<fnref n="8" />,* ID

José M. Iriondo<fnref n="1" />,* ID

### Correspondence

Javier Morente-López, Grupo de Ecología y Evolución en Islas, Dpto. Ciencias de la Vida y de la Tierra. Instituto de Productos Naturales y Agrobiología. Consejo Superior de Investigaciones Científicas. Avda. Astrofísico Sánchez, 3, 38206, San Cristóbal de La Laguna, Santa Cruz de Tenerife, Spain.

Email: javimorente@gmail.com

### Funding information

This work was funded by the projects AdAptA (CGL2012-33528) and EVA (CGL2016-77377-R) of the Spanish Ministry of Science and Innovation. JML was supported by a predoctoral fellowship (BES-2013-064951) and by a short research stay fellowship at the City College of New York (EEBB-I-16-11347) from the Spanish Ministry of Science and Innovation. JML was also supported by the REMEDINAL TE-CM project postdoctoral fellowship (S2018/EMT-

### Abstract

Environmental variation within a species' range can create contrasting selective pressures, leading to divergent selection and novel adaptations. The conservation value of populations inhabiting environmentally marginal areas remains in debate and is closely related to the adaptive potential in changing environments. Strong selection caused by stressful conditions may generate novel adaptations, conferring these populations distinct evolutionary potential and high conservation value under climate change. On the other hand, environmentally marginal populations may be genetically depauperate, with little potential for new adaptations to emerge. Here, we explored the use of ecological niche models (ENMs) linked with common garden experiments to predict and test for genetically determined phenotypic differentiation related to contrasting environmental conditions. To do so, we built an ENM for the alpine plant *Silene ciliata* in central Spain and conducted common garden experiments, assessing flowering phenology changes and differences in leaf cell resistance to extreme temperatures. The suitability patterns and response curves of the ENM led to the predictions that:

<footnote n="1">¹Área de Biodiversidad y Conservación, Depto. de Biología, Geología, Física y Química Inorgánica, ESCET, Universidad Rey Juan Carlos (URJC), Madrid, Móstoles, Spain</footnote>

<footnote n="2">²Island Ecology and Evolution Research Group, Institute of Natural Products and Agrobiology, Consejo Superior de Investigaciones Científicas (IPNA-CSIC), San Cristóbal de La Laguna, Tenerife, Spain</footnote>

<footnote n="3">³Department of Biology, City College of New York, City University of New York, New York, New York, USA</footnote>

<footnote n="4">⁴Ph.D. Program in Biology, Graduate Center, City University of New York, New York, New York, USA</footnote>

<footnote n="5">⁵Biodiversity and Biocomplexity Unit, Okinawa Institute of Science and Technology Graduate University, Kunigami-gun, Okinawa, Japan</footnote>

<footnote n="6">⁶Université de Lorraine, AgroParisTech, INRAE, Nancy, France</footnote>

<footnote n="7">⁷Facultad de Ciencias Naturales, Universidad Autónoma de Querétaro (FCN-UAQ), Santa Rosa Jáuregui, Querétaro, Mexico</footnote>

<footnote n="8">⁸Division of Vertebrate Zoology (Mammalogy), American Museum of Natural History, New York, New York, USA</footnote>

<footnote>*Equal contribution of senior authors.</footnote>
4338). CLR was supported by a Juan de la Cierva post-doctoral fellowship (MINECO: FJCI-2015-24712). JMK received support from the City College of New York, the Graduate Center of the City University of New York, a Japan Society for the Promotion of Science Postdoctoral Fellowship for Foreign Researchers, and the Okinawa Institute of Science and Technology Graduate University. The Consejo Superior de Investigaciones Científicas covered open access charges of the publication.

(1) the environmentally marginal populations experiencing less snowpack and higher minimum temperatures would have delayed flowering to avoid risks of late-spring frosts and (2) those with higher minimum temperatures and greater potential evapo-transpiration would show enhanced cell resistance to high temperatures to deal with physiological stress related to desiccation and heat. The common garden experiments revealed the expected genetically based phenotypic differentiation in flowering phenology. In contrast, they did not show the expected differentiation for cell resistance, but these latter experiments had high variance and hence lower statistical power. The results highlight ENMs as useful tools to identify contrasting putative selective pressures across species ranges. Linking ENMs with common garden experiments provides a theoretically justified and practical way to study adaptive processes, including insights regarding the conservation value of populations inhabiting environmentally marginal areas under ongoing climate change.

KEYWORDS

adaptation, climate change, ecological niche models, environmental gradients, marginal areas, phenotypic variation, physiology, response curves, species distribution models

# 1 | INTRODUCTION

Understanding key factors determining species distributions has long been of major interest in biology and now holds greater urgency given ongoing global environmental changes (Brown, 1984; Chen et al., 2011; Grinnell, 1917; Hutchinson, 1957; Parmesan, 2006; Root et al., 2003; Walther et al., 2002). One key concept, the center-periphery hypothesis, posits that geographically disparate populations across a species' range may have important differences, including variation in demographic, genetic, and phenotypic characteristics (Abeli et al., 2014; Gaston, 2009; Pironon et al., 2017; Sexton et al., 2009; Soule, 1973; Zakharov & Hellmann, 2008). Understanding the effects of such variation holds clear relevance under rapidly changing climate, land cover, and other environmental conditions because it may be crucial to facilitate adaptation and overcome the resulting threats (Abeli et al., 2014; Abeli & Orsenigo, 2018; Eckert et al., 2008; Gaston, 2009; Holt & Keitt, 2005; Kawecki, 2008; Pennington et al., 2021; Sexton et al., 2009; Swab et al., 2015). In this context, better predictions of future biodiversity patterns at and below the species level depend on two overlapping issues: (1) the relationship between geographic peripherality and ecological marginality, which directly influences patterns of genetic and phenotypic variation and (2) the potential adaptive value (or cost) of ecologically marginal populations under environmental change.

First, it is often assumed that environmental conditions are optimal in the center of a species' range and suboptimal at the periphery (Brown, 1984; Hoffmann & Blows, 1994; Holt & Keitt, 2005), with geographically peripheral populations also being ecologically marginal. Accordingly, it has been proposed that peripheral populations exhibit lower individual fitness, have lower abundance, and experience declining demographic trends (e.g., fertility, mortality) (Angert, 2006; Brown, 1984; Hengeveld & Haeck, 1982; Herlihy & Eckert,

2005; Villellas et al., 2013). However, geographical and ecological gradients do not necessarily coincide (Hardie & Hutchings, 2010; Pironon et al., 2015; Soley-Guardia et al., 2014; Soule, 1973), and the corresponding patterns regarding fitness, abundance, and demographic trends are not always consistent (Abeli et al., 2014; Hellmann et al., 2008; Holt, 2020; Osorio-Olvera et al., 2019; Pironon et al., 2017). Hence, environmentally marginal areas can exist not only in the geographical periphery of a species' range but also within it (Pironon et al., 2017). Defining environmental marginality is thus key for studying the ecological niches, adaptive processes, and conservation relevance of populations in such areas (Kawecki, 2008).

Second, opposing viewpoints have been presented regarding the genetic patterns and properties of geographically peripheral and/or environmentally marginal populations (Eckert et al., 2008; Kawecki, 2008). On one hand, lower survival and reproductive rates of such populations would lead to reduced effective population size and lower gene flow—and consequently, reductions in neutral and adaptive genetic diversity through genetic drift and inbreeding (Eckert et al., 2008; Kawecki, 2008). As a result, ecologically marginal populations may experience a significant genetic load due to maladaptation via low heterosis and an accumulation of deleterious alleles. On the other hand, strong selective pressures caused by environmentally stressful conditions may yield novel adaptations, conferring the population distinct evolutionary potential (Bontrager et al., 2021; Kawecki, 2008; Rolland et al., 2015). Indeed, the conservation value (or cost) of populations inhabiting environmentally marginal areas constitutes a long-standing debate now of high importance given rapid global environmental changes (Abeli et al., 2014; Abeli & Orsenigo, 2018; Channell, 2004; Hunter & Hutchinson, 1994; Lesica & Allendorf, 1995; Millar & Libby, 1991) because they may represent a source of valuable adaptations. Therefore, identifying genetically based
phenotypic variation related to contrasting environmental condi-
tions across species ranges is a key step in assessing local adap-
tive processes and their potential to influence range shifts under
changing climates (O'Neil et al., 2014; Razgour et al., 2019).

Correlative ecological niche models (ENMs) (or species dis-
tribution models; hereafter ENMs) have the potential to guide
experiments aimed at testing for genetically based phenotypic dif-
ferentiation. Correlative ENMs estimate habitat suitability based
on associations between occurrence records of the species and
predictor variables that characterize aspects of its environment
(Elith & Leathwick, 2009; Guevara et al., 2018; Phillips & Dudík,
2008). Ecological niche models model environmental suitability via
observed occurrence patterns, which may show mismatches with
related population-level phenomena such as abundance (Dallas &
Hastings, 2018), productivity (Serra-Diaz et al., 2013), and func-
tional traits (Bohner & Diez, 2020; Pagel et al., 2020; Pironon
et al., 2018). Hence, they do not necessarily represent an unbiased
characterization of the species' underlying niche, which defines
the environmental requirements of its full life cycle (Bohner &
Diez, 2020; Holt, 2009; Pagel & Schurr, 2012; Pagel et al., 2020;
Pironon et al., 2018; Soberón, 2007).

Nevertheless, under clear assumptions ENMs can be used to
study many important issues in ecology, evolution, biogeogra-
phy, and global change (Anderson, 2013). These include the cli-
matic drivers of species distributions (Guisan & Thuiller, 2005;
Peterson et al., 2011), including adaptive processes and responses
to climate change (Banta et al., 2012; Bush et al., 2016; Chardon
et al., 2020; Csergő et al., 2017; Ikeda et al., 2017; Lee-Yaw et al.,
2016; Lira-Noriega & Manthey, 2014; Mercer et al., 2016; Rolland
et al., 2015; Swab et al., 2015). Such studies typically assume no
functional genetic variation across the species' range for factors
related to its distributional ecology (Anderson, 2013; Hällfors
et al., 2016; Valladares et al., 2014; but see Capblancq et al.,
2020; Fitzpatrick & Keller, 2015). However, an adaptive response
can occur if contrasting selective pressures exist across the spe-
cies' range and relevant genes hold enough variation (Blanquart
et al., 2013; Kawecki & Ebert, 2004). In addition to indicating
low-suitability areas, ENMs can also identify the particular envi-
ronmental conditions associated with them, or putative limiting
factors. By doing so, they can aid in testing specific hypotheses
regarding local selection pressures, including predictions of the
genetically based phenotypic differentiation that would be ex-
pected (Dixon & Busch, 2017).

1.1 | Proposed approach and study system

Here, we propose an integrated framework aimed at identifying
putative environmental stressors using ENMs and testing for ex-
pected patterns of phenotypic differentiation via common garden
experiments. We implement this framework using a temperate-
latitude alpine plant. We first use ENMs to characterize envi-
ronmentally optimal and marginal areas and identify particular

putative stressors for the latter. We then undertake common gar-
den experiments to test for the genetically based phenotypic dif-
ferentiation predicted under each of the stressors for populations
inhabiting marginal areas. We consider as alpine those montane
areas above treeline (Körner, 2003). As alpine ecosystems are lo-
cally highly variable (i.e., small differences in elevation can cause
large changes in temperature, humidity, and exposure) (Körner,
2003), selective pressures can differ greatly over short distances
(Herrera & Bazaga, 2008). Accordingly, in applying this framework,
we take advantage of a well-studied Mediterranean cushion plant,
*Silene ciliata*, which is especially suited for evaluating evolution-
ary responses to ongoing climatic warming (e.g., García-Fernández
et al., 2012, 2013; Giménez-Benavides et al., 2011b, 2018). We
do so in an alpine ecosystem with a set of particularly challenging
environmental conditions that can prompt a wide variety of phe-
notypic adaptations (Körner, 2003).

Substantial natural history and experimental information exists for *S. ciliata*, informing general putative selection pressures. Adult plants of *S. ciliata* face critical issues regarding annual phenological timing for flowering and seed production. Plants must avoid flowering too early after spring snow-thaw due to the risk of losing their flowers to late-spring frosts (Wheeler et al., 2014), yet they must flower quickly enough to produce fruits and set seeds before summer drought and early-autumn frosts. Additionally, the species is dependent on water availability to complete the various stages of its life cycle, including flowering and fruiting as well as subsequent germination and establishment (Giménez-Benavides et al., 2007a, 2018). *Silene ciliata* plants also confront particular environmental stressors during their first stages of development. They must survive the dry conditions typical of summer and early autumn as well as the occasional early autumn frosts. Furthermore, those plants go on to experience an extended cold period during the winter and late spring. The species also must deal with the extremely hot conditions and in situ heat-related damage to plants in open alpine environments in the summer (Gauslaa, 1984; Körner, 2003; Körner & Larcher, 1988; Neuner et al., 2000). Thus, the leaf tissues of young plants must deal with physiological stress caused by drought and extreme temperatures (García-Fernández et al., 2013; Giménez-Benavides et al., 2007a, 2018).

Based on existing knowledge of environmental pressures
acting in Mediterranean alpine areas, we identified an initial set
of environmental variables. These best represent the extreme
environmental conditions experienced by the species during its
different life stages, show associations with its distributional
patterns, and characterize intra-range environmental variation
(minimum temperature of coldest month, precipitation of driest
and wettest months, mean snowpack in thaw months, and annual
potential evapotranspiration, see Appendix 1 for further detail).
Taking into account the species' natural history (García-Fernández
et al., 2014; Giménez-Benavides et al., 2007a, 2007b, 2018; Lara-
Romero et al., 2014), we then hypothesized that: (1) the populations
experiencing contrasting values of these five environmental vari-
ables would exhibit different genetically determined phenological
patterns and (2) populations experiencing contrasting values of these environmental variables would have leaf tissues in juvenile plants that exhibit genetically determined differentiation in cell resistance to extreme temperatures.

Guided by these general hypotheses for S. ciliata, we couple ENMs and common garden experiments to identify contrasting putative selection pressures across the species' distribution and test for them. We posit that the ENM developed for this species will be a valuable tool to define environmental marginality across geography and identify specific stressors that may cause phenotypic differentiation of populations. Additionally, in interpreting the results in the context of other available information for the species, we make inferences regarding the adaptive potential of populations inhabiting environmentally marginal areas. Altogether, this approach could help future research aimed at disentangling species' phenotypic differentiation and adaptive processes, as well as improve the accuracy of climate change forecasts of biodiversity.

# 2 | MATERIALS AND METHODS

## 2.1 | Research design

This research design links ENMs and common garden experiments following three major phases: pre-modeling, modeling, and post-modeling (Figure 1). We first formulated general hypotheses that identify potential environmental stressors likely to cause genetically based phenotypic differentiation between populations of the species (Step 1, Introduction). Subsequently, we selected the environmental variables to be used in ecological niche modeling according to these hypotheses (Step 2, Introduction). Next, we built an optimized ENM to predict habitat suitability across the landscape and delimitated environmentally marginal areas based on these values (Steps 3 and 4, Section 2.3). We then identified the environmental variables that held the strongest responses in the ENM (and thus would be most likely to exert local selective pressures) and assessed their response curves (Step 5, Section 2.3). Based on the hypotheses formulated in Step 1 and the modeled responses of the environmental variables identified via the ENM in Step 5 (in particular for the conditions present in the ecologically marginal areas), we made more specific predictions of expected genetically based phenotypic differentiation between populations regarding flowering phenology and cell resistance to extreme temperatures (Step 6, Section 3.2). We assessed these specific predictions via common garden experiments using populations from environmentally marginal and optimal areas to detect any genetically based phenotypic differentiation associated with relevant functional traits (Step 7, Section 2.4). Finally, we conducted statistical tests to determine if the results match the patterns expected (Step 8, Section 2.4). Taken together, these analyses tested whether populations occupying environmentally marginal and optimal environmental conditions showed the expected phenotypes linked to potentially adaptive genetic

differences that could be of conservation value under climate change.

## 2.2 | Study species

Silene ciliata Pourret (Caryophyllaceae) is a dwarf cushion, long-lived perennial plant (Giménez-Benavides et al., 2011) that is one of the dominant species of Mediterranean alpine habitats (i.e., those above tree line; Figure 1, Step 1). The study was undertaken in the Sistema Central of the Iberian Peninsula, a 500 km long mountain range that comprises the southernmost limit of the species' distribution (see Kyrkou et al., 2015; Tutin et al., 1964) (Appendix 2A). More specifically, we studied populations in each of its three main mountain subranges: Béjar (BJR), Gredos (GRD), and Guadarrama (GDM) (Figure A2.1). Silene ciliata experiences specific environmental stressors of the Mediterranean alpine ecosystem over its life stages and responds to stark environmental gradients over small geographic distances. In this region, populations that grow at lower elevations (1850 m) start and end the flowering period (late July until late August) earlier than those occurring at higher elevations (2400 m) (Morente-López et al., 2019). Seeds require exposure to cold conditions to break dormancy and germinate (i.e., cold stratification) (Baskin & Baskin, 1998; Giménez-Benavides et al., 2005). Effective pollen and seed dispersal distances are low and relatively invariant across populations (Lara-Romero et al., 2014, 2016).

## 2.3 | Ecological niche modeling

### 2.3.1 | Input data

We used species occurrence data in the Sistema Central obtained from the Global Biodiversity Information Facility (www.GBIF.org accessed October 2013), complemented with 10 new records from our fieldwork (with specimens deposited in the Real Jardín Botánico de Madrid herbarium; Appendix 2B). Silene ciliata is easily distinguishable from the other sympatric species of the same genus in the alpine areas of the Sistema Central (Silene boryi Boiss., which inhabits crag environments). Occurrence records with >10 m precision in geographic coordinates were removed. After cleaning the data, 120 total occurrence records were used for model building (Appendix 2C). To account for sampling biases (Boria et al., 2014; Merow et al., 2014) and associated artifactual spatial autocorrelation and environmental biases (Radosavljevic & Anderson, 2014), we used the R package spThin 0.1.0.1 to spatially thin our occurrence records by 0.75 km, yielding 87 records (Aiello-Lammens et al., 2015; Appendix 3). A model of the species' full range and inhabited conditions might produce different variable importance values and modeled response curves. However, the populations occurring in central Spain may differ in niche from those in other parts of the species' range. Because of this, and to
<figure><img src="image_4.png" /></figure>

<figure><img src="image_5.png" /><figcaption>FIGURE 1 Workflow for linking an ecological niche model (ENM) and common garden experiments, illustrated for analyses of the alpine plant *Silene ciliata* in the Sistema Central mountain range of the Iberian Peninsula. The steps fall in three major phases: pre-modeling (blue), modeling (yellow), and post-modeling (light red). Step 1: Generate hypotheses of potential environmental stressors that may act as selective pressures for the species (such as snowpack, precipitation regimes, or temperatures), based on previous knowledge. Step 2: Select fine-scale environmental variables that characterize the hypothesized stressors, to use as predictors in an ENM. Step 3: Build an optimized ENM aimed at characterizing the environmental gradients within the species' distribution in this region, as well as modeled responses to the variables. Step 4: Define environmental marginality (optimal and marginal areas) based on the habitat suitability values estimated by the ENM. The map shows the habitat suitability values from the ENM in the Gredos subrange and indicates the location of the particular populations (in optimal and marginal environments) used in later steps of the study. Step 5: Identify the environmental variables with the strongest effects on the ENM, assessing their response curves to determine the values associated with marginal areas. The plot shows the response curve for snowpack of thaw months, indicating the position of optimal and marginal populations. Step 6: Formulate specific predictions of genetically based phenotypic differentiation between optimal and marginal populations that would be expected under selection pressures associated with the identified variables. The image shows the predicted phenotypic differentiation in flowering phenology related with differences in snow cover. Step 7: Undertake common garden experiments to assess specific predictions of genetically based phenotypic differentiation. The photo shows the common garden experiment assessing flowering phenology. Step 8: Undertake statistical tests to determine if the results match the expected patterns. The plot shows differences in flowering rate across time between optimal and marginal populations. Photos in Steps 1 and 7 by Javier Morente-López</figcaption></figure>

focus on identifying and characterizing fine-grain environmental gradients within this area (i.e., those most likely to drive local adaptation), we selected one mountain range (Sistema Central) and characterized it as accurately as possible.

Our aim for modeling was explanatory rather than predictive, as we sought to identify putative driving variables and environmentally differentiated portions of the range rather than

just characterize suitability (Araújo et al., 2019; Araújo & Guisan, 2006; Elith & Leathwick, 2009; Merow et al., 2014). Thus, we selected as model predictors the five environmental variables that we hypothesized could best explain this species' distributional patterns and characterize intra-range variation and extremes for environmental conditions (see Section 1 and Appendix 1). All variables had a 200 m spatial resolution and were obtained or derived
<figure><img src="image_6.png" /></figure>

from the climatic digital atlas of the Iberian Peninsula (Ninyerola et al., 2005, accessed October 2016). Mean snowpack in thaw months was calculated following the methodology proposed by López-Moreno et al. (2007), reflecting the presence and potential magnitude of snow accumulation in a variable that ranges from 0 to 1. Environmental variables were constructed and managed using the R packages raster 2.8-19 (Hijmans & van Etten, 2012) and envirem 1.4 (Title & Bemmels, 2018). The pairwise correlation values of all variables were below 0.7, indicating relatively low collinearity (see Appendix 1 for detailed information on variable creation and selection).

### 2.3.2 | Ecological niche model parameterization, evaluation, and selection

We used the algorithm Maxent 3.3.3k (Phillips et al., 2006), a presence-background machine learning technique, to estimate the ecological niche of S. ciliata in the study region. We explored a range of settings affecting model complexity (combinations of feature classes and regularization multipliers) with the R package ENMeval 0.3.0 (Muscarella et al., 2014). To optimize model complexity, we first selected the candidate models with the lowest value of the Akaike Information Criterion corrected for small sample sizes ($AIC_c$), specifically those within 2 units of the best model (Warren & Seifert, 2011). We then evaluated these models via cross-validation performance metrics using a spatial checkerboard partition method ("checkerboard2" from the ENMeval R package 0.3.0; aggregation factor 2) (Muscarella et al., 2014). Threshold-dependent (omission rate calculated on withheld validation data after applying the 10 percent training omission threshold, or OR10) and threshold-independent metrics (AUC calculated on withheld validation data, or $AUC_{TEST}$ in ENMeval 0.3.0) were used to evaluate the models in terms of overfitting and discrimination, respectively (see Anderson & Gonzalez, 2011). Details of model development and evaluation including background selection, model settings used (such as combinations of feature classes and regularization multipliers), data partitioning, and performance metrics can be found in Appendix 4 (Anderson & Gonzalez, 2011; Anderson & Raza, 2010; Araújo & Guisan, 2006; Galante et al., 2017; Guevara et al., 2018; Merow et al., 2014; Phillips & Dudik, 2008; Saupe et al., 2012).

### 2.3.3 | Environmental variable importance: Defining environmental marginality, identifying putative selective pressures, and developing specific predictions

We visualized the optimal model in the full study area to indicate habitat suitability (Maxent's logistic transformation, bounded by 0 and 1; Phillips et al., 2006; Figure 1 and Figure A5.2), which we used to distinguish the degree of environmental marginality of

different areas. We classified portions of the study area as "optimal" and "marginal" based on the predicted suitability values for the spatially thinned localities. The optimal category was assigned to areas with habitat suitability values at or above the highest 33rd percentile of the suitability value distribution of ranked pixels, whereas the marginal category was assigned to those at or below the lowest 33rd percentile. As we focused on contrasting environmental conditions, we did not consider intermediate suitability values further. We then plotted the occurrence records in environmental space defined by the variables with non-zero coefficients in the model (see Appendix 5) to visualize environmental differences among the populations selected for addressing phenotypic differentiation (see below). We also identified which mountain subrange each of these occurrence records were associated with to determine if there was any environmental differentiation among subranges within the Sistema Central.

We then used the ENMs to develop specific predictions of the phenotypic patterns that would be expected under genetic differentiation of environmentally marginal populations. Based on the variable importance (Maxent's "permutation importance") and parameter weights (i.e., lambda values) from the ENM, we identified and focused on the environmental variables that contributed most to explaining the species' distribution (Galante et al., 2017; Phillips, 2005; Searcy & Shaffer, 2016); we did not consider further the variables that showed little importance for the model. We then plotted marginal response curves (model predictions for one variable when all other variables are set to their means) to determine the direction of the species' modeled responses to the environmental variables, which potentially describe local selective pressures. This allowed us to pinpoint the particular environmental conditions predicted to be stressful for populations inhabiting environmentally marginal localities (in contrast to optimal ones). Harnessing this information and existing knowledge of the crucial life stages of the species, we developed specific predictions regarding genetically based phenotypic differentiation and outlined expected patterns for phenological and physiological traits, as described in the next section.

### 2.4 | Genetically based phenotypic differentiation analysis

To assess associations of phenotype with habitat suitability values and environmental variables identified by the ENM, we selected functional traits with the potential to indicate phenological and physiological adaptations. First, we examined flowering phenological traits, which have been shown to be strongly affected by adaptation to changing environmental conditions (Bradshaw & Holzapfel, 2008; Chuine, 2010; Hoffmann & Sgró, 2011). Second, we characterized physiological traits associated with leaf cell resistance to cold and hot temperatures, representing the species' capacity to withstand extreme environmental pressures (Nobel et al., 2002).
<figure><img src="image_7.png" /></figure>

## 2.4.1 | Common garden experiments and functional trait measurements

We then intersected the ENM results with a set of *S. ciliata* populations considered in an ongoing project addressing adaptive patterns associated with genetically based phenotypic differentiation (Morente-López et al., 2020). Genetically based phenotypic differentiation can be detected through the variation in phenotypes found between different groups of individuals grown under the same controlled environmental conditions (since the differences that emerge are related to genetic composition of the groups). Specifically, we characterized the examined populations of *S. ciliata* as either “optimal” or “marginal” according to the habitat suitability values predicted by the ENM. In each of the three mountain subranges, we considered two populations: one located in an area predicted as “optimal” (high habitat suitability) and one predicted as “marginal” (low habitat suitability but still sufficient to support reproductive populations, see Section 2.3 for thresholding implementation) (Table 1; Figure 1 – Steps 3 and 4, Figure A5.2).

### Flowering phenology

In the late summer of 2013, we collected a minimum of 30 randomly selected adult plants from each of the six populations. From these, we obtained rooted cuttings to generate a sufficient number of genetically identical individuals to overcome potential losses during rooting and growth. The entire collection was allowed to grow in the Universidad Rey Juan Carlos CULTIVE laboratory (690 m, https://urjc-cultive.webnode.es) under common garden conditions for 7 months to minimize carry-over effects from the original environment (Appendix 6A). The diameter of all selected individuals of each population (Table 1) was measured on the date when the first plant started flowering (15 April 2014). From then until 25 September 2014, when the last new flower appeared, we counted the total number of new flowers of each plant weekly (n = 23). From these data, we calculated four individual-based flowering traits as representations of the temporal patterns of reproductive investment: flowering onset, peak, end, and duration (details in Appendix 6A).

### Leaf cell resistance to extreme temperatures

The purpose of this study was not to obtain accurate measures of absolute temperature tolerances of individual plants in their natural habitats, but rather to assess whether relative differences existed in assays of temperature tolerance performed under shared artificial conditions, between plants originating from populations in optimal versus marginal environments. With this goal, and recognizing that laboratory conditions will not match the complex environmental sequences experienced in the field, we designed a simple experiment to evaluate possible patterns of genetically based differentiation with all plants exposed to one common state. To do so, we collected seeds from 30 randomly selected plants from each of the six populations. Seeds were germinated and seedlings grown for 8 months (juvenile plants for this species) under

<table><caption>TABLE 1 Six *Silene ciliata* populations were chosen for phenotypic differentiation analyses under common garden conditions (Population and ID), corresponding to one environmentally marginal and one optimal population for each mountain subrange (Environment; Figure 1). Environmental classification in optimal and marginal areas was made based on habitat suitability (HS) values from the optimal Maxent ecological niche model (logistic transformation, see Methods). Geographic coordinates and elevation are provided for each population. Phenology (flowering onset, peak, end and duration) and leaf cell resistance to extreme temperatures (lethal temperature, LT50) were measured in each of the populations for given numbers of individuals and genets under common garden conditions (see Methods and Appendix 6)</caption><thead><tr><th rowspan="2">Population</th><th rowspan="2">ID</th><th rowspan="2">Mountain subrange</th><th rowspan="2">Environment</th><th rowspan="2">HS values</th><th rowspan="2">Latitude</th><th rowspan="2">Longitude</th><th colspan="3">Phenology trait analysis</th><th colspan="2">LT50 traits analysis</th></tr><tr><th>Nº</th><th>Individuals</th><th>Nº Genets</th><th>Nº Individuals</th><th>Nº Genets</th></tr></thead><tbody><tr><td>Canchal Negro</td><td>NEG</td><td>Béjar</td><td>Optimal</td><td>0.80</td><td>40°20'19.97"N</td><td>5°41'22.27"W</td><td>2360</td><td>54</td><td>31</td><td>6</td><td>6</td></tr><tr><td>Pico El Aguila</td><td>AGI</td><td>Béjar</td><td>Marginal</td><td>0.51</td><td>40°21'12.36"N</td><td>5°41'46.52"W</td><td>1950</td><td>79</td><td>46</td><td>6</td><td>6</td></tr><tr><td>Altos del Morezón</td><td>ZON</td><td>Gredos</td><td>Optimal</td><td>0.85</td><td>40°14'57.5"N</td><td>5°16'8.3"W</td><td>2380</td><td>57</td><td>29</td><td>6</td><td>6</td></tr><tr><td>El Sestil</td><td>SES</td><td>Gredos</td><td>Marginal</td><td>0.43</td><td>40°16'24.45"N</td><td>5°14'54.93"W</td><td>1900</td><td>103</td><td>33</td><td>6</td><td>6</td></tr><tr><td>Pico de Peñalara</td><td>PEÑ</td><td>Guadarrama</td><td>Optimal</td><td>0.87</td><td>40°51'2.11"N</td><td>3°57'24.02"W</td><td>2400</td><td>65</td><td>35</td><td>6</td><td>6</td></tr><tr><td>Najarra baja</td><td>NAJ</td><td>Guadarrama</td><td>Marginal</td><td>0.29</td><td>40°49'23.46"N</td><td>3°49'52.53"W</td><td>1850</td><td>74</td><td>43</td><td>6</td><td>6</td></tr></tbody></table>
common garden conditions in the CULTIVE laboratory. Although leaves may change physiologically over the seasons, the simplified experiment conducted here used two treatments that were subjected to all plants at a single growth stage. Tolerance of juvenile plant leaf cells to extreme temperatures was determined in May–June 2014 by measuring electrolyte leakage, an indicator of cell membrane integrity (Didden-Zopfy & Nobel, 1982; Drennan, 2009; Soto-Correa et al., 2013). The lethal temperature, defined as the temperature at which half of the maximum electrolyte leakage was reached, was estimated for each plant separately for both low (hereafter LT50<sub>LOW</sub>) and high (hereafter LT50<sub>HIGH</sub>) incubation temperature sequences (see Nobel et al., 2002 detailed description). LT50 values have been found to be good predictors of extreme temperature tolerance in the field (Drennan, 2009; Nobel et al., 2002). Each low (5 to -20°C) and high (30 to 55°C) sequence was composed of six temperatures (5°C intervals), for which the initial temperature acts as a baseline of the electrolytic conductivity for undamaged tissue. The temperatures selected can be reached in this ecosystem under natural conditions (see Section 1 and Appendix 6). Six plants per population were randomly selected for the analysis, and six 2-mm leaf discs were taken from each plant for testing at multiple temperatures (see Appendix 6 for detailed information).

### 2.4.2 | Data analysis

#### Categorical approach

We used two different statistical approaches to analyze whether the optimal and marginal environmental categories identified through ENMs were associated with differences among populations in the functional traits. For flowering phenology, we used mixed effects Cox proportional hazards regression models implemented in the *coxme* function from the R package *coxme* 2.2-7 (Therneau, 2018). We modeled rate of flowering over time with plant size as a covariate, environmental class as a fixed factor with two levels (optimal/marginal), and population of origin as a random factor. In the case of leaf cell resistance to extreme temperatures, we tested the association between physiological traits (LT50<sub>HIGH</sub> and LT50<sub>LOW</sub>) and environmental classes using linear mixed models implemented in the *lmer* function from the R package *lme4* 1.1-17 (Bates et al., 2014) with the same model design explained above for phenological traits (Appendix 7).

#### Continuous approach

We used linear models (R Core Team, 2013) to assess whether habitat suitability values predict genetically based phenotypic differentiation. We assessed the association of: (1) flowering phenology (onset, peak, duration, and end of flowering) and (2) cell resistance to extreme temperatures (LT50<sub>HIGH</sub> and LT50<sub>LOW</sub>) with continuous habitat suitability values and plant size (Appendix 7).

## 3 | RESULTS

### 3.1 | Ecological niche model selection, evaluation, and spatial projection

Four ENMs were identified as co-optimal. They all had delta AIC<sub>c</sub> values <2, similar values for other metrics of performance (averages of test omission rates and AUC, each calculated on data withheld in cross-validations), and the same feature class combination (Table A8.1). For further analysis, we selected the model with the lowest AIC<sub>c</sub> value, which was also one of the simplest overall (see model 1 in Table A8.1). When we projected this ENM to the full background spatial extent, the highest habitat suitability values were concentrated along mountain crests around 2400 m, and lower elevations corresponded to decreasing habitat suitability values that were <0.1 below 1800 m (Figure 1 Step 3; Figure A5.2).

### 3.2 | Environmental variables and specific predictions

The modeled responses for the predictor variables differed in weights and complexity. Snowpack constituted the most important variable, followed by minimum temperature of the coldest month, potential evapotranspiration, and precipitation of the driest month (Table 2). Marginal response curves showed the shapes of the relationships between each of the environmental variables and habitat suitability (logistic transformation, Figure 2; raw outputs in Appendix 8B, Figure

<table><caption>TABLE 2 The five input environmental variables used to develop the final Maxent ecological niche model (ENM) for *Silene ciliata* in the Sistema Central mountain range of the Iberian Peninsula. These variables reflect factors considered likely to affect habitat suitability based on prior ecological knowledge of the species. These include snowpack (mean snow accumulation ratio in the thaw months); minimum temperature (annual minimum temperature of the coldest month); annual potential evapotranspiration (PET); precipitation driest (precipitation of the driest month); and precipitation wettest (precipitation of the wettest month). These environmental variables were then ranked based on their permutation importance obtained from the ENM. Feature classes used and lambda values (coefficient weights) from the ENM are provided for each variable (see Section 2 and Appendix 4)</caption><thead><tr><th>Variable</th><th>Permutation importance</th><th>Feature class</th><th>Lambda value coefficient</th></tr></thead><tbody><tr><td rowspan="2">Snowpack</td><td rowspan="2">76.89</td><td>Linear</td><td>0.88</td></tr><tr><td>Hinge</td><td>-2.48</td></tr><tr><td>Minimum temperature</td><td>10.32</td><td>Linear</td><td>-4.52</td></tr><tr><td>PET</td><td>7.19</td><td>Quadratic</td><td>-3.27</td></tr><tr><td>Precipitation driest</td><td>5.60</td><td>Quadratic</td><td>-2.84</td></tr><tr><td>Precipitation wettest</td><td>0</td><td>Linear</td><td>0</td></tr></tbody></table>
<figure><img src="image_8.png" /></figure>

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

FIGURE 2 Response curves for four environmental variables from the Maxent ecological niche model (ENM) for *Silene ciliata* in the Sistema Central mountain range of the Iberian Peninsula: (a) minimum temperature of coldest month, (b) snowpack in thaw months, (c) annual potential evapotranspiration, and (d) precipitation of driest month were the four variables identified based on the permutation importance of the final ENM (Table 2). These marginal response curves show the shape of the relationship between the environmental variable and habitat suitability values ranging from 0.0 (lowest) to 1.0 (highest; Maxent's logistic transformation, see Section 2 and Appendix 8B) when all other variables are set to their means. This allows identification of the environmental conditions predicted to be stressful for populations inhabiting environmentally marginal localities. Based on this and existing knowledge regarding crucial life stages of the species, we generated specific predictions of expected genetically based phenotypic differentiation of phenological and physiological traits (see Section 2 and Figure 3). Vertical bars represent values for each of the six populations selected for the common garden experiments (Figure 1 and Table 1; for environmental and habitat suitability values associated with each population, see Appendix 11, Table A11.1). Note that: As minimum temperature of the coldest month increases (a), habitat suitability values decrease. Habitat suitability values are close to zero until snowpack in thaw months is positive (b; presence of snow cover). Precipitation of the driest month (c) and potential evapotranspiration (d) showed similar responses, with gradually decreasing habitat suitability values as each variable increased. Minimum temperature of the coldest month (a) and mean snowpack in the thaw months (b) showed the greatest contrast in suitability between optimal and marginal populations. Thus, these two variables were crucial to generate specific expectations regarding phenotypic differentiation

A8.1). Minimum temperature of the coldest month and snowpack showed the strongest contrast between optimal and marginal populations, making these variables the most likely to distinguish populations experiencing contrasting environmental conditions (Figure 2). As minimum temperature of the coldest month increased, habitat suitability decreased, dropping rapidly around -8°C and approaching zero at around 3°C (Figure 2a). Habitat suitability was close to zero until snowpack values were positive (i.e., presence of snow cover), where habitat suitability spiked suddenly and reached a maximum at snowpack values of 1 (Figure 2b). Potential evapotranspiration and precipitation of the driest month showed responses similar to each other (Figure 2c and d), with gradually decreasing habitat

suitability as each variable increased. However, whereas potential evapotranspiration showed a consistent trend regarding the habitat suitability of optimal vs. marginal populations, precipitation of the driest month showed no clear contrast between them (Figure 2d). Tellingly, populations selected as “optimal” and “marginal” fell in different subsets of the environmental space defined by the two most important variables, snowpack and minimum temperature of coldest month (Appendix 5, Figure A5.1a). No major environmental differentiation was found among mountain subranges (Appendix 5, Figure A5.1b).

The response curves led to specific predictions for functional traits related to flowering phenology and leaf cell resistance to extreme
temperatures, refinements of the original general hypotheses (see Section 2.2). Considering the environmental variables with high importance for the model, response curves indicated that lack of snowpack, higher minimum temperature of the coldest month, and higher potential evapotranspiration led to lower suitability. Translation of these modeled patterns into specific predictions for genetically based phenotypic differentiation produced the following expected differences between populations occurring in contrasting environments:

## Flowering phenology

Environmentally marginal areas defined by the model have higher temperatures of the coldest month (which prompt a more active metabolism) and lack of snowpack (which would have provided effective frost protection). This combination can result in higher mortality risk for flowers when late-spring frosts occur (Hennon et al., 2012; Kreyling, 2010; Schaberg et al., 2008). Therefore, we predicted that populations in such areas experience selective pressure to delay flowering. If populations have evolved according to this pressure, we would expect them to have a later flowering onset, peak, and end, as well as a longer flowering duration.

## Leaf cell resistance to extreme temperatures

Environmentally marginal and optimal areas identified by the model differed substantially in temperature and potential evapotranspiration (which is highly correlated with maximum temperature of the warmest month, see Appendix 1). Marginal areas correspond to warmer conditions and higher potential evapotranspiration values (which prompt physiological cell stress to high temperatures and desiccation) (Zandalinas et al., 2018). In contrast, optimal ones are consistently colder and with lower values of PET (with modeled suitability that remains high even for the coldest temperatures within the study region). High temperatures and higher values of PET, especially in late summer and early autumn when water availability is very low, can lead to high physiological stress related to desiccation and greater mortality risk for plants (García-Fernández et al., 2013). Therefore, we predicted that populations in environmentally marginal areas would experience selective pressure for leaf cells that are resistant to higher temperatures—but that populations in optimal areas would not experience analogous pressure for greater resistance to extreme cold. If populations have evolved according to this pressure, we would expect marginal populations to have greater values of LT50HIGH but that the optimal populations would not have lower values of LT50LOW (see Methods 2.5).

## 3.3 | Functional traits

### Flowering phenology

Results of both the categorical approach (environmental classification into optimal and marginal areas) and the continuous one

(directly using habitat suitability values) matched the predictions for flowering phenology (Section 3.2). The effect of environmental classification based on the ENM was significant in the expected direction for all phenological variables (Figure 3; Appendix 9, Table A9.1). Samples from populations located in environmentally marginal areas displayed a later flowering onset, peak, and end than those inhabiting optimal areas (Figure 3a-c; Table A9.1). Additionally, flowering season duration was longer in samples from populations in marginal areas (Figure 3d; Table A9.1). Similarly, a significant negative relationship existed between habitat suitability values and each of the four flowering phenology variables (Table A9.1): as suitability values increased, the onset, peak, and end of flowering were all earlier, and the duration of flowering decreased (Figure 4). Parameter estimates of all developed models are in Appendix 10 (Tables A10.1 and A10.2).

## Leaf cell resistance to extreme temperatures

Results for both the environmental classification into optimal and marginal areas (categorical approach) and the habitat suitability values (continuous approach) matched predictions for LT50LOW but not for LT50HIGH (Section 3.2). As expected, no significant associations were found between environmental classification or habitat suitability values and cell resistance to extreme low temperatures (LT50LOW); however, contrary to expectations for cell resistance to extreme high temperatures, no significant associations existed for either characterization of suitability (LT50HIGH, Figure 5; Appendix 9, Table A9.2). Parameter estimates of all developed models are in Appendix 10 (Tables A10.3 and A10.4).

# 4 | DISCUSSION

## 4.1 | Overall findings and utility of this framework

The ENM characterized specific putative selective pressures in environmentally marginal areas, and common garden experiments confirmed most of the expected differences in phenotypic traits. Specifically, the weights and responses of the model allowed refining the original hypotheses and developing specific predictions of expected phenotypic differentiation between marginal and optimal populations in flowering phenology and leaf cell resistance to extreme temperatures. The common garden experiments confirmed the expectation that marginal populations delay and elongate the flowering period, as predicted in response to selective pressure to avoid the negative effects of late-spring frost. Experiments regarding leaf cell tolerance to extreme cold and hot temperatures matched the expectation for cell resistance for optimal populations but not for the marginal ones, based on predicted selective pressures to avoid the negative effects of high temperatures and desiccation. In relevant sections below, we discuss caveats to these interpretations (including variation and statistical power), many of which point toward fruitful directions for future research.
<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /><figcaption>FIGURE 3 Differences in flowering phenology in common garden experiments between populations of *Silene ciliata* from the Sistema Central of the Iberian Peninsula from environmentally marginal and optimal conditions (categorical approach). Environmental classification into optimal and marginal areas was based on the habitat suitability values of the ecological niche model (Table 1; Figure A5.2). Samples from populations located in marginal areas displayed significantly later flowering onset (a), peak (b), and end (c) than those from populations inhabiting optimal areas; furthermore, flowering season duration (d) was significantly longer for samples from populations in marginal areas (*p* < .05, **p** < .001). These results matched predictions made for all four phenological traits measured, based on the modeled response curves and existing knowledge of the species (Figure 2).</figcaption></figure>

Taken together, the results support the proposal that at least in some cases ENMs can be useful for identifying areas within a species' range that are environmentally marginal, as well as for developing refined predictions regarding expected functional genetic

differentiation in such areas. Hence, ENMs can guide experiments with the potential for identifying adaptive population divergence driven by divergent selective pressures in marginal areas (Hällfors et al., 2016). Extrapolating the results found in this study to other
<figure><img src="image_12.png" /></figure>

<figure><img src="image_13.png" /><figcaption>FIGURE 4 Phenological traits related to flowering based on common garden experiments for *Silene ciliata* in the Sistema Central of the Iberian Peninsula, in relation to habitat suitability values from the optimal ecological niche model (continuous approach). Habitat suitability had a significant negative relationship with each of the four phenology variables (***p < .001). As habitat suitability decreased, the onset (a), peak (b), and end (c) of flowering were all later, and the duration (d) of flowering increased. These results matched predictions made for all phenological traits measured, based on the modeled response curves and existing knowledge of the species (Figure 2)</figcaption></figure>

areas or species should be undertaken with caution, due to the
species- and context-dependency of functional differentiation
(Bohner & Diez, 2020; Pagel et al., 2020; Pironon et al., 2018;
Thuiller et al., 2010; Wittmann et al., 2016). Indeed, this framework
aims to guide the identification of mismatches between suitability
and underlying genetics and/or functional performance. Hence,
it adds to a growing body of research addressing interesting—and

important—mismatches between suitability and factors such as
survival, reproductive success, population abundance and growth
rates, and both adaptive and neutral genetic patterns (Bohner &
Diez, 2020; Dallas & Hastings, 2018; Pagel et al., 2020; Pironon
et al., 2018; Serra-Diaz et al., 2013). The methodological framework
proposed here could be applicable across species and study systems
that vary widely in the amount of previous knowledge available.
<figure><img src="image_14.png" /></figure>

<figure><img src="image_15.png" /><figcaption>FIGURE 5 Traits related to leaf cell resistance to extreme high and low temperatures (LT50HIGH and LT50LOW, see Section 2.5) based on common garden experiments for *Silene ciliata* in the Sistema Central of the Iberian Peninsula. Results are shown in relation to habitat suitability values from the final ecological niche model (a and c; continuous approach), and from environmentally marginal and optimal conditions (b and d; categorical approach). Environmental classification into optimal and marginal areas was based on the habitat suitability values of the ecological niche model (Table 1, Figure A5.2). Neither habitat suitability values nor environmental classification had any significant relationship with LT50HIGH or LT50LOW (*p* > .005)</figcaption></figure>

When little is known, ENMs can be useful for identifying novel hy-
potheses for exploration. In contrast, if previous knowledge exists
and general hypotheses can be created *a priori* as in this study, the
models can be useful for selecting which hypothesized stressors are
most likely and for refining their corresponding predictions.

## 4.2 | Interpretation of responses to putative selective pressures

In this study, the key factor leading to low modeled habitat suitabil-
ity (i.e., the strongest putative driver of selective pressure) was low
snowpack in thaw months. Snowpack in thaw months is driven by

spring temperature and precipitation (López-Moreno et al., 2007)
and is related to the timing and success of crucial life stages such
as plant flowering in S. ciliata (Giménez-Benavides et al., 2007a,
2011) and other temperate alpine plant species (Dunne et al., 2003;
Giménez-Benavides et al., 2018; Inouye, 2008; Inouye et al., 2002;
Lara-Romero et al., 2014). Moreover, strong evidence has been
found for the role of snowpack in local adaptative processes, particu-
larly in alpine plants (Anderson & Wadgymar, 2020; Inouye, 2020).
Snowpack thaw directly affects vulnerability to spring freezing of
various plant tissues and the growth patterns of alpine plant species
(Jonas et al., 2008; Wheeler et al., 2014). In this context, the late-
flowering strategy (seed-risk sensu Molau, 1993) followed in general
by S. ciliata would be especially prevalent in populations inhabiting
areas with shorter snow cover duration and, consequently, higher risk of late-spring frosts (Wheeler et al., 2014). This result highlights the ability of ENMs to identify one or more particular environmental factors strongly associated with great intra-range differences in habitat suitability, even for species inhabiting landscapes with high local environmental heterogeneity. Notably, snowpack is greatly affected by microtopographic variation, especially in alpine areas where wind scouring occurs on exposed ridges (DeBeer & Pomeroy, 2009); this could be a source of uncertainty in estimating habitat suitability values, especially at fine grains. In any case, marked ecological niche differences related with fine-grain environmental variables also have been found in other alpine species between populations inhabiting central and marginal areas (Giménez-Benavides et al., 2007b; Papuga et al., 2018).

Indeed, genetic differentiation was found in *S. ciliata* for phenological traits related to flowering, as expected for adaptation to avoid late-spring frosts in areas with low snowpack and in concordance with earlier studies. When grown in common garden settings, samples from populations in marginal areas exhibited later flowering start and peak and flowered for a longer time than those from populations inhabiting optimal ones. It is important to consider the possibility that the phenological variation found could be due at least in part to effects of the environments of origin because plants used in the common garden experiment derive from cuttings obtained from individuals that grew in natural conditions. However, the patterns found are congruent with other analyses of selection, based on the effects of phenological trait differences in fitness for the same species at the same sites, which showed that flowering onset and duration are under selection in populations at different elevations (Giménez-Benavides et al., 2011b, see also Morente-López et al., 2019). The patterns of phenological differentiation found are also congruent with results obtained in subsequent experiments implemented with plants grown from seeds derived from the cuttings used in the present study (unpublished results). Moreover, a parallel transcriptomic study with the same populations identified SNPs (single nucleotide polymorphisms) with unusually high allele-frequency differentiation between marginal and optimal environments, for candidate genes in the flower development pathway (Sacristán-Bajo et al., 2019). In addition, the phenotypic differentiation found in this study adds to recent evidence of the adaptive potential of marginal populations at early life stages (i.e., seedlings). In a reciprocal field seed-sowing experiment involving populations used in this study (Morente-López et al., 2020), seeds from marginal populations showed higher seedling survival rate when sown in both optimal and marginal areas. Additionally, in that study gene flow between marginal populations (by artificial pollination) improved seedling survival, but gene flow from optimal to marginal populations did not.

Although common garden experiments such as this one only can demonstrate the existence of genetically based differentiation (not adaptation *per se*), the information from other studies mentioned above provides additional evidence toward the conclusion that the phenological differentiation found here likely is the result of divergent selection (Giménez-Benavides et al., 2011; Morente-López

et al., 2019; Sacristán-Bajo et al., 2019). Alternatively, such differentiation could result from neutral processes unrelated to selection, such as genetic drift (Garcia-Ramos & Kirkpatrick, 1997). Therefore, future research to corroborate the adaptive nature of the differentiation found here could involve the use of genomic approaches to quantify the magnitude of genetic differentiation and the nature of the genes involved, as well as reciprocal in situ common garden experiments with fitness measurements.

Similar patterns also have been observed in co-occurring species in the same ecosystems (Lara-Romero et al., 2014), suggesting that phenological reproductive traits are likely to be under strong selection. In fact, they have been shown to be critical for plant fitness (Fox, 2003; Levin, 2006), capable of evolving rapidly under climatic variation (Franks et al., 2007), and influencing species distributions and potential responses to climate change (Banta et al., 2012; Chuine, 2010). This also agrees with studies showing that these traits generally show high evolvability (Kawakami et al., 2011; Méndez-Vigo et al., 2013; Volis, 2011) and with a meta-analysis indicating consistent selection patterns (Munguía-Rosas et al., 2011).

Overall, the environmentally marginal populations of *S. ciliata* that show genetically determined phenotypes expected to perform better with low snowpack may be advantageous under climate change. In the Iberian Peninsula, temperature is expected to increase, and precipitation decrease, causing a reduction in snow cover in the coming decades. Nevertheless, as discussed above, extrapolating the patterns found in this study to other parts of the distribution of *S. ciliata* or to other species should be undertaken only with caution, because functional differentiation of plant populations related with abiotic factors is commonly species- and context-dependent (Leuschner et al., 2009; Papuga et al., 2018; Wagner et al., 2011).

In contrast to the differentiation found for flowering phenology, the lack of expected association between low habitat suitability values and functional differentiation in the cell resistance of leaves may reflect the absence of strong selective pressures concerning high temperatures *per se*, but other non-exclusive possibilities also exist. Notably, the importance of temperature and PET to the model (which led to the expectations regarding cell stress to extreme temperatures and desiccation) was much lower than that of snowpack (which more directly led to the expectations regarding flowering phenology). Previous research in the same ecosystem demonstrated no differences in leaf functional traits linked to water-use efficiency (e.g., carbon and nitrogen isotope ratios of organic material) or economic spectrum (e.g., specific leaf area and leaf dry matter content), along an elevational gradient for *S. ciliata* and co-occurring species (Pescador et al., 2015). Other studies in both plants and animals have shown that phenological traits are more sensitive to divergent selection than physiological ones (Bradshaw & Holzapfel, 2008; Hoffmann & Sgró, 2011). For instance, evolved differences in thermal resistance in organisms are usually weak or absent among populations occupying even a broad range of climatic conditions, which may be related to a lack of genetic variation for these traits at the population level (Bradshaw & Holzapfel, 2008; Hoffmann & Sgró, 2011, but see
Sinclair et al., 2012). A previous study performing a drought experiment using the same species, although limited to one population per habitat, observed that plants derived from marginal populations had greater survival than those from optimal ones (García-Fernández et al., 2013). This fact highlights that decoupling the effects of drought and high temperature on plants is very difficult, particularly for seedlings and juvenile plants. Additionally, the statistical power of the analyses developed in the present study regarding temperature tolerances was hampered by the limited sample size (number of plants studied per population) and the unexpectedly high variance of the results. To reach more conclusive results, additional experiments would be required with a larger sample size. In addition, it would be relevant to test resistance to extreme temperatures at different developmental stages of the plant (e.g., seedlings, mature plants). Indeed, for traits more related to drought resistance physiological differences may be more likely to arise in earlier crucial life stages (such as germination or seedling establishment), (García-Fernández et al., 2013, 2014; Giménez-Benavides et al., 2007a, 2008), whereas other traits may differ at later developmental stages.

## 5 | CONCLUSIONS AND FUTURE DIRECTIONS

Using ecological niche modeling in conjunction with laboratory and field experiments represents a fruitful way to study adaptive processes, including the evolutionary capabilities of species confronted with climate change (Dixon & Busch, 2017; Swab et al., 2015). Our findings demonstrated experimentally that ENMs have the potential not only to define environmental marginality (Lira-Noriega & Manthey, 2014; Madeira de Madeiros et al., 2018), but also to predict genetically determined variation in a functional trait among populations (Manel et al., 2012; Thuiller et al., 2010; Wittmann et al., 2016). Moreover, this study highlights the value of conserving marginal populations during a time of rapidly changing global climate (Channell, 2004) based on their adaptive potential, which can derive from ecologically relevant genetic uniqueness (Papuga et al., 2018). Although the present results correspond to a single species, the implications and applicability of the proposed methodological framework for evolutionary biology—especially in the context of climate change—should be generalizable for other species and spatial scales.

In addition to correlative models, mechanistic ones may improve predictions of species' responses to climate change. Such models can consider factors such as survival and reproductive rates, dispersal capability, competition-facilitation relationships, and measured physiological responses of species. However, they require much more data and the establishment of larger-scale and more complex common garden experiments implemented across a variety of environmental conditions (Benito Garzón et al., 2011; Benito-Garzón et al., 2019; Kearney & Porter, 2009; Valladares et al., 2014). Although the approach we propose does not include the data

necessary to parameterize a mechanistic model, it can provide fundamental information and insights that can lead to further experiments that better facilitate mechanistic models. Either approach—or a combination—likely will prove feasible for taxa conducive to common garden treatments (e.g., plants or arthropods, Pelini et al., 2009; Williams et al., 2012). Paired with landscape genetics and genomic approaches, studies such as this one that combine modeling and experiments should lead to a greater unified understanding of adaptive processes under changing climates and other anthropogenic alterations of the environment (Exposito-Alonso et al., 2018, 2019; Fournier-Level et al., 2011).

## ACKNOWLEDGMENTS

We thank Alfredo García Fernández, María Luisa Rubio Teso, Raquel Ruiz Checa, and Pablo Tabarés Siville for their help with fieldwork; José Margalet and Carlos Diaz for support and technical recommendations in the CULTIVE facility laboratory greenhouse; and Greta Carrete Vega and Pilar Hurtado Aragües for assistance in making figures. We thank Lázaro Guevara López for support and useful discussions during the conception of this study. We thank the staff of the Parque Nacional de Guadarrama (Dirección General de Medio Ambiente, Comunidad de Madrid) and the Delegación Territorial de Salamanca y Ávila (Servicio Territorial de Medio Ambiente, Junta Castilla y León) for permission to conduct fieldwork in the area. For funding, we thank the Spanish Ministry of Science and Innovation (AdAptA and EVA projects) and the Madrid Autonomous Region Government (REMEDINAL TE-CM project). Insightful discussions with Brent Sinclair and Jessica Hellmann regarding modeled environmental stressors and expected physiological adaptation inspired our integrative research approach.

## CONFLICT OF INTEREST

The authors have no conflict of interest to declare.

## AUTHOR CONTRIBUTIONS

JML, JMI, RPA, and JMK designed the framework integrating niche models and phenotypic assays. JML, JMI and CLR designed the experimental study, collected the samples, and conducted the field work. JSC undertook the leaf cell resistance experiment laboratory work. JML developed the niche models in collaboration with JMK, RPA, and JMSD. JML analyzed phenotypic data with the help of CLR. JML wrote the paper with the help of RPA, JMI, JMK, CLR, and JMSD. All authors reviewed the paper and approved the final manuscript.

## DATA AVAILABILITY STATEMENT

The three databases that support the findings of this study are openly available in public repositories. Presence data is stored in figshare at https://doi.org/10.6084/m9.fishshare.13562618.v1. Cell resistance to extreme temperatures data is stored in figshare at https://doi.org/10.6084/m9.fishshare.19368872.v1. Flowering phenology data is stored in URJC e-cienciaDatos public repository at https://doi.org/10.21950/XLYLQK. [Correction added on 14-May-2022, after first online publication: DATA AVAILABILITY STATEMENT has been added.]
<figure><img src="image_16.png" /></figure>

## ORCID

Javier Morente-López https://orcid.org/0000-0001-9141-8581

Jamie M. Kass https://orcid.org/0000-0002-9432-895X

Carlos Lara-Romero https://orcid.org/0000-0003-0962-0567

Josep M. Serra-Diaz https://orcid.org/0000-0003-1988-1154

Robert P. Anderson https://orcid.org/0000-0002-7706-4649

José M. Iriondo https://orcid.org/0000-0003-2710-3889

## REFERENCES

Abeli, T., Gentili, R., Mondoni, A., Orsenigo, S., & Rossi, G. (2014). Effects of marginality on plant population performance. *Journal of Biogeography*, *41*(2), 239–249. https://doi.org/10.1111/jbi.12215

Abeli, T., Vamosi, J. C., & Orsenigo, S. (2018). The importance of marginal population hotspots of cold-adapted species for research on climate change and conservation. *Journal of Biogeography*, *45*(5), 977–985. https://doi.org/10.1111/jbi.13196

Aiello-Lammens, M. E., Boria, R. A., Radosavljevic, A., Vilela, B., & Anderson, R. P. (2015). spThin: An R package for spatial thinning of species occurrence records for use in ecological niche models. *Ecography*, *38*(5), 541–545. https://doi.org/10.1111/ecog.01132

Anderson, J. T., & Wadgymar, S. M. (2020). Climate change disrupts local adaptation and favours upslope migration. *Ecology Letters*, *23*(1), 181–192. https://doi.org/10.1111/ele.13427

Anderson, R. P. (2013). A framework for using niche models to estimate impacts of climate change on species distributions. *Annals of the New York Academy of Sciences*, *1297*(1), 8–28. https://doi.org/10.1111/nyas.12264

Anderson, R. P., & Gonzalez, I. (2011). Species-specific tuning increases robustness to sampling bias in models of species distributions: An implementation with Maxent. *Ecological Modelling*, *222*(15), 2796–2811. https://doi.org/10.1016/j.ecolmodel.2011.04.011

Anderson, R. P., & Raza, A. (2010). The effect of the extent of the study region on GIS models of species geographic distributions and estimates of niche evolution: Preliminary tests with montane rodents (*genus Nephelomys*) in Venezuela. *Journal of Biogeography*, *37*(7), 1378–1393. https://doi.org/10.1111/j.1365-2699.2010.02290.x

Angert, A. L. (2006). Demography of central and marginal populations of monkeyflowers (*Mimulus cardinalis* and *M. lewisii*). *Ecology*, *87*(8), 2014–2025.

Araújo, M. B., Anderson, R. P., Barbosa, A. M., Beale, C. M., Dormann, C. F., Early, R., Garcia, R. A., Guisan, A., Maiorano, L., Naimi, B., O'Hara, R. B., Zimmermann, N. E., & Rahbek, C. (2019). Standards for distribution models in biodiversity assessments. *Science Advances*, *5*(1), 1–10. https://doi.org/10.1126/sciadv.aat4858

Araújo, M. B., & Guisan, A. (2006). Five (or so) challenges for species distribution modelling. *Journal of Biogeography*, *33*(10), 1677–1688. https://doi.org/10.1111/j.1365-2699.2006.01584.x

Banta, J. A., Ehrenreich, I. M., Gerard, S., Chou, L., Wilczek, A., Schmitt, J., & Purugganan, M. D. (2012). Climate envelope modelling reveals intraspecific relationships among flowing phenology, niche breadth and potential range size in *Arabidopsis thaliana*. *Ecology Letters*, *15*(8), 769–777. https://doi.org/10.1111/j.1461-0248.2012.01796.x

Baskin, C. C., & Baskin, J. M. (1998). *Seeds: Ecology, biogeography, and evolution of dormancy and germination*. Elsevier.

Bates, D., Mächler, M., Bolker, B., & Walker, S. (2014). Fitting linear mixed-effects models using lme4. *ArXiv Preprint*, ArXiv:1406.5823.

Benito Garzón, M., Alía, R., Robson, T. M., & Zavala, M. A. (2011). Intraspecific variability and plasticity influence potential tree species distributions under climate change. *Global Ecology and Biogeography*, *20*(5), 766–778. https://doi.org/10.1111/j.1466-8238.2010.00646.x

Benito-Garzón, M., Robson, T. M., & Hampe, A. (2019). ΔTraitSDMs: Species distribution models that account for local adaptation and

phenotypic plasticity. *New Phytologist*, *222*(4), 1757–1765. https://doi.org/10.1111/nph.15716

Blanquart, F., Kaltz, O., Nuismer, S. L., & Gandon, S. (2013). A practical guide to measuring local adaptation. *Ecology Letters*, *16*(9), 1195–1205. https://doi.org/10.1111/ele.12150

Bohner, T., & Diez, J. (2020). Extensive mismatches between species distributions and performance and their relationship to functional traits. *Ecology Letters*, *23*(1), 33–44. https://doi.org/10.1111/ele.13396

Bontrager, M., Usui, T., Lee-Yaw, J. A., Anstett, D. N., Branch, H. A., Hargreaves, A. L., Muir, C. D., & Angert, A. L. (2021). Adaptation across geographic ranges is consistent with strong selection in marginal climates and legacies of range expansion. *Evolution*, *75*(6), 1316–1333. https://doi.org/10.1111/evo.14231

Boria, R. A., Olson, L. E., Goodman, S. M., & Anderson, R. P. (2014). Spatial filtering to reduce sampling bias can improve the performance of ecological niche models. *Ecological Modelling*, *275*, 73–77. https://doi.org/10.1016/j.ecolmodel.2013.12.012

Bradshaw, W. E., & Holzapfel, C. M. (2008). Genetic response to rapid climate change: It's seasonal timing that matters. *Molecular Ecology*, *17*(1), 157–166. https://doi.org/10.1111/j.1365-294X.2007.03509.x

Brown, J. H. (1984). On the relationship between abundance and distribution of species. *The American Naturalist*, *124*(2), 255–279. https://doi.org/10.1086/284267

Bush, A., Mokany, K., Catullo, R., Hoffmann, A., Kellermann, V., Sgrò, C., McEvey, S., & Ferrier, S. (2016). Incorporating evolutionary adaptation in species distribution modelling reduces projected vulnerability to climate change. *Ecology Letters*, *19*(12), 1468–1478. https://doi.org/10.1111/ele.12696

Capblancq, T., Fitzpatrick, M. C., Bay, R. A., Exposito-Alonso, M., & Keller, S. R. (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/annurev-econsys-020720-042553

Channell, R. (2004). The conservation value of peripheral populations: the supporting science. In *Proceedings of the species at risk 2004 pathways to recovery conference* (pp. 1–17). Species at Risk 2004 Pathways to Recovery Conference Organizing Committee.

Chardon, N. I., Pironon, S., Peterson, M. L., & Doak, D. F. (2020). Incorporating intraspecific variation into species distribution models improves distribution predictions, but cannot predict species traits for a wide-spread plant species. *Ecography*, *43*(1), 60–74. https://doi.org/10.1111/ecog.04630

Chen, I. C., Hill, J. K., Ohlemüller, R., Roy, D. B., & Thomas, C. D. (2011). Rapid range shifts of species associated with high levels of climate warming. *Science*, *333*(6045), 1024–1026. https://doi.org/10.1126/science.1206432

Chuine, I. (2010). Why does phenology drive species distribution? *Philosophical Transactions of the Royal Society B: Biological Sciences*, *365*(1555), 3149–3160. https://doi.org/10.1098/rstb.2010.0142

Csergő, A. M., Salguero-Gómez, R., Broennimann, O., Coutts, S. R., Guisan, A., Angert, A. L., Welk, E., Stott, I., Enquist, B. J., McGill, B., Svenning, J.-C., Violle, C., & Buckley, Y. M. (2017). Less favourable climates constrain demographic strategies in plants. *Ecology Letters*, *20*(8), 969–980. https://doi.org/10.1111/ele.12794

Dallas, T. A., & Hastings, A. (2018). Habitat suitability estimated by niche models is largely unrelated to species abundance. *Global Ecology and Biogeography*, *27*(12), 1448–1456. https://doi.org/10.1111/geb.12820

DeBeer, C. M., & Pomeroy, J. W. (2009). Modelling snow melt and snowcover depletion in a small alpine cirque, Canadian Rocky Mountains. *Hydrological Processes*, *227*(23), 2584–2599. https://doi.org/10.1002/hyp

Didden-Zopfy, B., & Nobel, P. S. (1982). High temperature tolerance and heat acclimation of *Opuntia bigelovii*. *Oecologia*, *52*(2), 176–180. https://doi.org/10.1007/BF00363833
Dixon, A. L., & Busch, W. (2017). Common garden test of range limits as predicted by a species distribution model in the annual plant Mimulus bicolor. *American Journal of Botany*, 104(6), 817–827.
https://doi.org/10.3732/ajb.1600414

Drennan, P. M. (2009). Temperature influences on plant species of arid a semi-arid regions with emphasis on CAM succulents. In E. Barrera, & W. K. Smith (Eds.), *Perspectives in biophysical plant ecophysiology* (pp. 57–94). UNAM.

Dunne, J. A., Harte, J., & Taylor, K. J. (2003). Subalpine meadow flowing phenology responses to climate change: Integrating experimental and gradient methods Author (s): Jennifer A. Dunne, John Harte and Kevin J. Taylor Published by: Wiley on behalf of the Ecological Society of America Stable. *Ecological Monographs*, 73(1), 69–86.

Eckert, C. G., Samis, K. E., & Lougheed, S. C. (2008). Genetic variation across species' geographical ranges: The central-marginal hypothesis and beyond. *Molecular Ecology*, 17(5), 1170–1188. https://doi.org/10.1111/j.1365-294X.2007.03659.x

Elith, J., & Leathwick, J. R. (2009). Species distribution models: Ecological explanation and prediction across space and time. *Annual Review of Ecology, Evolution, and Systematics*, 40, 677–697. https://doi.org/10.1146/annurev.ecolsys.110308.120159

Exposito-Alonso, M., Exposito-Alonso, M., Gómez Rodríguez, R., Barragán, C., Capovilla, G., Chae, E., & Weigel, D. (2019). Natural selection on the *Arabidopsis thaliana* genome in present and future climates. *Nature*, 573(7772), 126–129. https://doi.org/10.1038/s41586-019-1520-9

Exposito-Alonso, M., Vasseur, F., Ding, W., Wang, G., Burbano, H. A., & Weigel, D. (2018). Genomic basis and evolutionary potential for extreme drought adaptation in *Arabidopsis thaliana*. *Nature Ecology & Evolution*, 2(2), 352–358. https://doi.org/10.1038/s41559-017-0423-0

Fitzpatrick, M. C., & Keller, S. R. (2015). Ecological genomics meets community-level modelling of biodiversity: Mapping the genomic landscape of current and future environmental adaptation. *Ecology Letters*, 18(1), 1–16. https://doi.org/10.1111/ele.12376

Fournier-Level, A., Korte, A., Cooper, M. D., Nordborg, M., Schmitt, J., & Wilczek, A. M. (2011). A map of local adaptation in *Arabidopsis thaliana*. *Science*, 334(6052), 86–89. https://doi.org/10.1126/science.1209271

Fox, G. A. (2003). Assortative mating and plant phenology: Evolutionary and practical consequences. *Evolutionary Ecology Research*, 5(1), 1–18.

Franks, S. J., Sim, S., & Weis, A. E. (2007). Rapid evolution of flowering time by an annual plant in response to a climate fluctuation. *PNAS*, 104(4), 1278–1282. https://doi.org/10.1073/pnas.0608379104

Galante, P. J., Alade, B., Muscarella, R., Jansa, S. A., Goodman, S. M., & Anderson, R. P. (2017). The challenge of modeling niches and distributions for data-poor species: A comprehensive approach to model complexity. *Ecography*, 40(5), 1–10. https://doi.org/10.1111/ecog.02909

García-Fernández, A., Escudero, A., Lara-Romero, C., & Iriondo, J. M. (2014). Effects of the duration of cold stratification on early life stages of the Mediterranean alpine plant *Silene ciliata*. *Plant Biology*, 17(2), 344–350. https://doi.org/10.1111/plb.12226

García-Fernández, A., Iriondo, J. M., Bartels, D., & García, A. (2013). Response to artificial drying until drought-induced death in different elevation populations of a high-mountain plant. *Plant Biology*, 15, 93–100. https://doi.org/10.1111/j.1438-8677.2012.00638.x

García-Fernández, A., Segarra-Moragues, J. G., Widmer, A., Escudero, A., & Iriondo, J. M. (2012). Unravelling genetics at the top: Mountain islands or isolated belts? *Annals of Botany*, 110(6), 1221–1232. https://doi.org/10.1093/aob/mcs195

García-Ramos, G., & Kirkpatrick, M. (1997). Genetic models of adaptation and gene flow in peripheral populations. *Evolution*, 51(1), 21. https://doi.org/10.2307/2410956

Gaston, K. J. (2009). Geographic range limits: Achieving synthesis. *Proceedings of the Royal Society B: Biological Sciences*, 276(1661), 1395–1406. https://doi.org/10.1098/rspb.2008.1480

Gauslaa, Y. (1984). Heat resistance and energy budget in different Scandinavian plants. *Ecography*, 7(1), 5–6. https://doi.org/10.1111/j.1600-0587.1984.tb01098.x

Giménez-Benavides, L., Albert, M. J., Iriondo, J. M., & Escudero, A. (2011). Demographic processes of upward range contraction in a long-lived Mediterranean high mountain plant. *Ecography*, 34(1), 85–93. https://doi.org/10.1111/j.1600-0587.2010.06250.x

Giménez-Benavides, L., Escudero, A., García-Camacho, R., García-Fernández, A., Iriondo, J. M., Lara-Romero, C., & Morente-López, J. (2018). How does climate change affect regeneration of Mediterranean high-mountain plants? An integration and synthesis of current knowledge. *Plant Biology*, 20, 50–62. https://doi.org/10.1111/plb.12643

Giménez-Benavides, L., Escudero, A., & Iriondo, J. M. (2007a). Local adaptation enhances seedling recruitment along an altitudinal gradient in a high mountain mediterranean plant. *Annals of Botany*, 99(4), 723–734. https://doi.org/10.1093/aob/mcm007

Giménez-Benavides, L., Escudero, A., & Iriondo, J. M. (2007b). Reproductive limits of a late-flowering high-mountain Mediterranean plant along an elevational climate gradient. *New Phytologist*, 173(2), 367–382. https://doi.org/10.1111/j.1469-8137.2006.01932.x

Giménez-Benavides, L., Escudero, A., & Iriondo, J. M. (2008). What shapes the altitudinal range of a high mountain Mediterranean plant? Recruitment probabilities from ovule to seedling stage. *Ecography*, 31(6), 731–740. https://doi.org/10.1111/j.0906-7590.2008.05509.x

Giménez-Benavides, L., Escudero, A., & Pérez-García, F. (2005). Seed germination of high mountain Mediterranean species: Altitudinal, interpopulation and interannual variability. *Ecological Research*, 20(4), 433–444. https://doi.org/10.1007/s11284-005-0059-4

Giménez-Benavides, L., García-Camacho, R., Iriondo, J. M., & Escudero, A. (2011). Selection on flowering time in Mediterranean high-mountain plants under global warming. *Evolutionary Ecology*, 25(4), 777–794. https://doi.org/10.1007/s10682-010-9440-z

Grinnell, J. (1917). The niche-relationships of the California Thrasher. *The Auk*, 34(4), 427–433. https://doi.org/10.2307/4072271

Guevara, L., Gerstner, B. E., Kass, J. M., & Anderson, R. P. (2018). Toward ecologically realistic predictions of species distributions: A cross-time example from tropical montane cloud forests. *Global Change Biology*, 24(4), 1511–1522. https://doi.org/10.1111/gcb.13992

Guisan, A., & Thuiller, W. (2005). Predicting species distribution: Offering more than simple habitat models. *Ecology Letters*, 8(9), 993–1009. https://doi.org/10.1111/j.1461-0248.2005.00792.x

Hällfors, M. H., Liao, J., Dzurisin, J., Grundel, R., Hyvärinen, M., Towle, K., Wu, G. C., & Hellmann, J. J. (2016). Addressing potential local adaptation in species distribution models: Implications for conservation under climate change. *Ecological Applications*, 26(4), 1154–1169. https://doi.org/10.1890/15-0926

Hardie, D. C., & Hutchings, J. A. (2010). Evolutionary ecology at the extremes of species' ranges. *Environmental Reviews*, 18, 1–20. https://doi.org/10.1139/A09-014

Hellmann, J. J., Pelini, S. L., Prior, K. M., & Dzurisin, J. D. K. (2008). The response of two butterfly species to climatic variation at the edge of their range and the implications for poleward range shifts. *Oecologia*, 157(4), 583–592. https://doi.org/10.1007/s00442-008-1112-0

Hengeveld, R., & Haeck, J. (1982). The distribution of abundance. *I. Measurements. Journal of Biogeography*, 9, 303–316. https://doi.org/10.2307/2844717

Hennon, P. E., D'Amore, D. V., Schaberg, P. G., Wittwer, D. T., & Shanley, C. S. (2012). Shifting climate, altered niche, and a dynamic conservation strategy for yellow-cedar in the north pacific coastal rainforest. *BioScience*, 62(2), 147–158. https://doi.org/10.1525/bio.2012.62.2.8
<figure><img src="image_17.png" /></figure>

Herlihy, C. R., & Eckert, C. G. (2005). Evolution of self-fertilization at geographic range margins? A comparison of demographic, floral, and mating system variables in central vs. peripheral populations of *Aquilegia canadensis* (Ranunculaceae). *American Journal of Botany*, 92(4), 744–751.

Herrera, C. M., & Bazaga, P. (2008). Population-genomic approach reveals adaptive floral divergence in discrete populations of a hawk moth-pollinated violet. *Molecular Ecology*, 17(24), 5378–5390.
https://doi.org/10.1111/j.1365-294X.2008.04004.x

Hijmans, R. J., & van Etten, J. (2012). raster: Geographic analysis and modeling with raster data. R Package Version, 2.

Hoffmann, A. A., & Blows, M. W. (1994). Species borders: ecological and evolutionary perspectives. *Trends in Ecology and Evolution*, 9(6), 223–227.
https://doi.org/10.1016/0169-5347(94)90248-8

Hoffmann, A. A., & Sgró, C. M. (2011). Climate change and evolutionary adaptation. *Nature*, 470(7335), 479–485.
https://doi.org/10.1038/nature09670

Holt, R. D. (2009). Bringing the Hutchinsonian niche into the 21st century: Ecological and evolutionary perspectives. *Proceedings of the National Academy of Sciences of the United States of America*, 106(Supplement 2), 19659–19665.
https://doi.org/10.1073/pnas.0905137106

Holt, R. D. (2020). Reflections on niches and numbers. *Ecography*, 43(3), 387–390.
https://doi.org/10.1111/ecog.04828

Holt, R. D., & Keitt, T. H. (2005). Species' borders: A unifying theme in ecology. *Oikos*, 108(1), 3–6.

Hunter, M. L., & Hutchinson, A. (1994). The virtues and shortcomings of parochialism: Conserving species that are locally rare, but globally common. *Conservation Biology*, 8(4), 1163–1165.
https://doi.org/10.1046/j.1523-1739.1994.08041163.x

Hutchinson, G. E. (1957). Concluding remarks cold spring harbor symposia on quantitative biology. *GS SEARCH*, 22, 415–427.

Ikeda, D. H., Max, T. L., Allan, G. J., Lau, M. K., Shuster, S. M., & Whitham, T. G. (2017). Genetically informed ecological niche models improve climate change predictions. *Global Change Biology*, 23(1), 164–176.
https://doi.org/10.1111/gcb.13470

Inouye, D. W. (2008). Effects of climate change on phenology, frost damage and floral abundance of montane wild flowers. *Ecology*, 89(2), 353–362.
https://doi.org/10.1890/06-2128.1

Inouye, D. W. (2020). Effects of climate change on alpine plants and their pollinators. *Annals of the New York Academy of Sciences*, 1469(1), 26–37.
https://doi.org/10.1111/nyas.14104

Inouye, D. W., Morales, M. A., & Dodge, G. J. (2002). Variation in timing and abundance of flowering by *Delphinium barbeyi* Huth (Ranunculaceae): the roles of snowpack, frost, and La Niña, in the context of climate change. *Oecologia*, 130(4), 543–550.
https://doi.org/10.1007/s00442-001-0835-y

Jonas, T., Rixen, C., Sturm, M., & Stoeckli, V. (2008). How alpine plant growth is linked to snow cover and climate variability. *Journal of Geophysical Research: Biogeosciences*, 113(3), 1–10.
https://doi.org/10.1029/2007JG000680

Kawakami, T., Morgan, T. J., Nippert, J. B., Ocheltree, T. W., Keith, R., Dhakal, P., & Ungerer, M. C. (2011). Natural selection drives clinal life history patterns in the perennial sunflower species, *Helianthus maximiliani*. *Molecular Ecology*, 20(11), 2318–2328.
https://doi.org/10.1111/j.1365-294X.2011.05105.x

Kawecki, T. J. (2008). Adaptation to marginal habitats. *Annual Review of Ecology, Evolution, and Systematics*, 39(1), 321–342.
https://doi.org/10.1146/annurev.ecolsys.38.091206.095622

Kawecki, T. J., & Ebert, D. (2004). Conceptual issues in local adaptation. *Ecology Letters*, 7(12), 1225–1241.
https://doi.org/10.1111/j.1461-0248.2004.00684.x

Kearney, M., & Porter, W. (2009). Mechanistic niche modelling: Combining physiological and spatial data to predict species' ranges. *Ecology Letters*, 12(4), 334–350.
https://doi.org/10.1111/j.1461-0248.2008.01277.x

Körner, C. (2003). *Alpine plant life. Functional plant ecology of high mountain ecosystems*. Springer-Verlag. https://doi.org/10.1007/978-3-642-18970-8

Körner, C., & Larcher, W. (1988). Plant life in cold climates. In *Symposia of the society for experimental biology* (Vol. 42, pp. 25–57).

Kreyling, J. (2010). Winter climate change: A critical factor for temperature vegetation performance. *Ecology*, 91(7), 1939–1948.
https://doi.org/10.1890/09-1160.1

Kyrkou, I., Iriondo, J. M., & García-Fernández, A. (2015). A glacial survivor of the alpine Mediterranean region: phylogenetic and phylogeographic insights into *Silene ciliata* Pourr. (Caryophyllaceae). *PeerJ*, 3, 1–19.
https://doi.org/10.7717/peerj.1193

Lara-Romero, C., García-Camacho, R., Escudero, A., & Iriondo, J. M. (2014). Genetic variation in flowering phenology and reproductive performance in a Mediterranean high-mountain specialist, *Armeria caespitosa* (Plumbaginaceae). *Botanical Journal of the Linnean Society*, 176(3), 384–395.
https://doi.org/10.1111/boj.12208

Lara-Romero, C., García-Fernández, A., Robledo-Arnuncio, J. J., Roumet, M., Morente-López, J., López-Gil, A., & Iriondo, J. M. (2016). Individual spatial aggregation correlates with between-population variation in fine-scale genetic structure of *Silene ciliata* (Caryophyllaceae). *Heredity*, 116(5), 417–423.
https://doi.org/10.1038/hdy.2015.102

Lara-Romero, C., Robledo-Arnuncio, J. J., García-Fernández, A., & Iriondo, J. M. (2014). Assessing intraspecific variation in effective dispersal along an altitudinal gradient: A test in two Mediterranean high-mountain plants. *PLoS One*, 9(1), e87189.
https://doi.org/10.1371/journal.pone.0087189

Lee-Yaw, J. A., Kharouba, H. M., Bontrager, M., Mahony, C., Csergő, A. M., Noreen, A. M. E., Li, Q., Schuster, R., & Angert, A. L. (2016). A synthesis of transplant experiments and ecological niche models suggests that range limits are often niche limits. *Ecology Letters*, 19(6), 710–722.
https://doi.org/10.1111/ele.12604

Lesica, P., & Allendorf, F. W. (1995). When are peripheral populations valuable for conservation? *Conservation Biology*, 9(4), 753–760.
https://doi.org/10.1046/j.1523-1739.1995.09040753.x

Leuschner, C., Köchemann, B., & Buschmann, H. (2009). Forest ecology and management abundance, niche breadth, and niche occupation of Central European tree species in the centre and at the margin of their distribution range. *Forest Ecology*, 258(7), 1248–1259.
https://doi.org/10.1016/j.foreco.2009.06.020

Levin, D. A. (2006). Flowering phenology in relation to adaptive radiation. *Systematic Botany*, 31(2), 239–246.
https://doi.org/10.1600/03636406777585928

Lira-Noriega, A., & Manthey, J. D. (2014). Relationship of genetic diversity and niche centrality: A survey and analysis. *Evolution*, 68(4), 1082–1093.
https://doi.org/10.1111/evo.12343

López-Moreno, J. I., Vicente-Serrano, S. M., & Lanjeri, S. (2007). Mapping snowpack distribution over large areas using GIS and interpolation techniques. *Climate Research*, 33(3), 257.
https://doi.org/10.3354/cr033257

Madeira de Madeiros, C., Hernández-Lambraño, R. E., Felix Ribeiro, K. A., & Sánchez Agudo, J. Á. (2018). Living on the edge: Do central and marginal populations of plants differ in habitat suitability? *Plant Ecology*, 219(9), 1029–1043.

Manel, S., Gugerli, F., Thuiller, W., Alvarez, N., Legendre, P., Holderegger, R., Gielly, L., & Taberlet, P. (2012). Broad-scale adaptive genetic variation in alpine plants is driven by temperature and precipitation. *Molecular Ecology*, 21(15), 3729–3738.
https://doi.org/10.1111/j.1365-294X.2012.05656.x

Marcer, A., Méndez-Vigo, B., Alonso-Blanco, C., & Picó, F. X. (2016). Tackling intraspecific genetic structure in distribution models better reflects species geographical range. *Ecology and Evolution*, 6(7), 2084–2097.
https://doi.org/10.1002/ece3.2010

Méndez-Vigo, B., Gomaa, N. H., Alonso-Blanco, C., & Picó, F. X. (2013). Among- and within-population variation in flowering time of
Iberian *Arabidopsis thaliana* estimated in field and glasshouse conditions. *New Phytologist*, 197(4), 1332-1343. https://doi.org/10.1111/nph.12082

Merow, C., Smith, M. J., Edwards, T. C., Guisan, A., McMahon, S. M., Normand, S., Thuiller, W., Wüest, R. O., Zimmermann, N. E., & Elith, J. (2014). What do we gain from simplicity versus complexity in species distribution models? *Ecography*, 37(12), 1267-1281. https://doi.org/10.1111/ecog.00845

Millar, C. I., & Libby, W. J. (1991). Strategies for conserving clinal, ecotypic, and disjunct population diversity in widespread species. *Genetics and Conservation of Rare Plants*, 149, 170.

Molau, U. (1993). Relationships between flowering phenology and life history strategies in tundra plants. *Arctic and Alpine Research*, 25(4), 391. https://doi.org/10.2307/1551922

Morente-López, J., Lara-Romero, C., García-Fernandez, A., Rubio Teso, M. L., Prieto-Benítez, S., & Iriondo, J. M. (2020). Gene flow effects on populations inhabiting marginal areas: Origin matters. *Journal of Ecology*, 109(1), 139-153. https://doi.org/10.1111/1365-2745.13455

Morente-López, J., Scheepens, J. F., Lara-Romero, C., Ruiz-Checa, R., Tabarés, P., & Iriondo, J. M. (2019). Past selection shaped phenological differentiation among populations at contrasting elevations in a Mediterranean alpine plant. *Environmental and Experimental Botany*, 170, 103894. https://doi.org/10.1016/j.envexpbot.2019.103894

Munguía-Rosas, M. A., Ollerton, J., Parra-Tabla, V., & De-Nova, J. A. (2011). Meta-analysis of phenotypic selection on flowering phenology suggests that early flowering plants are favoured. *Ecology Letters*, 14(5), 511-521. https://doi.org/10.1111/j.1461-0248.2011.01601.x

Muscarella, R., Galante, P. J., Soley-Guardia, M., Boria, R. A., Kass, J. M., Uriarte, M., & Anderson, R. P. (2014). ENMeval: An R package for conducting spatially independent evaluations and estimating optimal model complexity for Maxent ecological niche models. *Methods in Ecology and Evolution*, 5(11), 1198-1205. https://doi.org/10.1111/2041-210X.12261

Neuner, G., Buchner, O., & Braun, V. (2000). Short-term changes in heat tolerance in the alpine cushion plant *Silene acaulis* ssp. excapa [All.] J. Braun at different altitudes. *Plant Biology*, 2(06), 677-683.

Ninyerola, M., Pons, X., & Roure, J. M. (2005). *Atlas climático digital de la Península Ibérica: metodología y aplicaciones en bioclimatología y geobotánica*. Universitat Autónoma de Barcelona, Departament de Biologia Animal, Biologia Vegetal i Ecologia (Unitat de Botánica).

Nobel, P. S., Barrera, E. D., Beilman, D. W., Doherty, J. H., Zutta, R., Nobel, P. S., & Zutta, B. R. (2002). Temperature limitations for cultivation of edible cacti in California. *California Botanical Society*, 49(4), 228-236.

O'Neil, S. T., Dzurisin, J. D. K., Williams, C. M., Lobo, N. F., Higgins, J. K., Deines, J. M., Carmichael, R. D., Zeng, E., Tan, J. C., Wu, G. C., Emrich, S. J., & Hellmann, J. J. (2014). Gene expression in closely related species mirrors local adaptation: Consequences for responses to a warming world. *Molecular Ecology*, 23(11), 2686-2698. https://doi.org/10.1111/mec.12773

Osorio-Olvera, L., Soberón, J., & Falconi, M. (2019). On population abundance and niche structure. *Ecography*, 42(8), 1415-1425. https://doi.org/10.1111/ecog.04442

Pagel, J., & Schurr, F. M. (2012). Forecasting species ranges by statistical estimation of ecological niches and spatial population dynamics. *Global Ecology and Biogeography*, 21(2), 293-304. https://doi.org/10.1111/j.1466-8238.2011.00663.x

Pagel, J., Treurnicht, M., Bond, W. J., Kraaij, T., Nottebrock, H., Schutte-Vlok, A. L., Tonnabel, J., Esler, K. J., & Schurr, F. M. (2020). Mismatches between demographic niches and geographic distributions are strongest in poorly dispersed and highly persistent plant species. *Proceedings of the National Academy of Sciences of the United States of America*, 117(7), 3663-3669. https://doi.org/10.1073/pnas.1908684117

Papuga, G., Gauthier, P., Pons, V., Farris, E., & Thompson, J. D. (2018). Ecological niche differentiation in peripheral populations: A

comparative analysis of eleven Mediterranean plant species. *Ecography*, 41(10), 1650-1664. https://doi.org/10.1111/ecog.03331

Parmesan, C. (2006). Ecological and evolutionary responses to recent climate change. *Annual Review of Ecology, Evolution, and Systematics*, 37, 637-669. https://doi.org/10.2307/annurev.ecols ys.37.091305.30000024

Pelini, S. L., Dzurisin, J. D. K., Prior, K. M., Williams, C. M., Marsico, T. D., Sinclair, B. J., & Hellmann, J. J. (2009). Translocation experiments with butterflies reveal limits to enhancement of poleward populations under climate change. *Proceedings of the National Academy of Sciences of the United States of America*, 106(27), 11160-11165. https://doi.org/10.1073/pnas.0900284106

Pennington, L. K., Slatyer, R. A., Ruiz-Ramos, D. V., Veloz, S. D., & Sexton, J. P. (2021). How is adaptive potential distributed within species ranges? *Evolution*, 75(9), 2152-2166. https://doi.org/10.1111/ evo.14292

Pescador, D. S., De Bello, F., Valladares, F., & Escudero, A. (2015). Plant trait variation along an altitudinal gradient in mediterranean high mountain grasslands: Controlling the species turnover effect. *PLoS One*, 10(3), 1-16. https://doi.org/10.1371/journal.pone.0118876

Peterson, A. T., Soberón, J., Pearson, R. G., Anderson, R. P., Martínez-Meyer, E., Nakamura, M., & Araújo, M. B. (2011). *Ecological niches and geographic distributions*. Princeton University Press. (MPB-49)

Phillips, S. J. (2005). A brief tutorial on Maxent. *AT&T Research*, 190(4), 231-259.

Phillips, S. J., Anderson, R. P., & Schapire, R. E. (2006). Maximum entropy modeling of species geographic distributions. *Ecological Modelling*, 190(3-4), 231-259. https://doi.org/10.1016/j.ecolmod.2005.03.026

Phillips, S. J., & Dudík, M. (2008). Modeling of species distributions with Maxent: New extensions and a comprehensive evaluation. *Ecography*, 31(2), 161-175. https://doi.org/10.1111/j.2007.0906-7590.05203.x

Pironon, S., Papuga, G., Villellas, J., Angert, A. L., García, M. B., & Thompson, J. D. (2017). Geographic variation in genetic and demographic performance: New insights from an old biogeographical paradigm. *Biological Reviews*, 92(4), 1877-1909. https://doi.org/10.1111/brv.12313

Pironon, S., Villellas, J., Morris, W. F., Doak, D. F., & García, M. B. (2015). Do geographic, climatic or historical ranges differentiate the performance of central versus peripheral populations? *Global Ecology and Biogeography*, 24(6), 611-620. https://doi.org/10.1111/geb.12263

Pironon, S., Villellas, J., Thuiller, W., Eckhart, V. M., Geber, M. A., Moeller, D. A., & García, M. B. (2018). The 'Hutchinsonian niche' as an assemblage of demographic niches: Implications for species geographic ranges. *Ecography*, 41(7), 1103-1113. https://doi.org/10.1111/ecog.03414

R Core Team. (2013). *R: A language and environment for statistical computing*. R Foundation for Statistical Computing.

Radosavljevic, A., & Anderson, R. P. (2014). Making better Maxent models of species distributions: Complexity, overfitting and evaluation. *Journal of Biogeography*, 41(4), 629-643. https://doi.org/10.1111/jbi.12227

Razgour, O., Forester, B., Taggart, J. B., Bekaert, M., Juste, J., Ibáñez, C., Puechmaille, S. J., Novella-Fernandez, R., Alberdi, A., & Manel, S. (2019). Considering adaptive genetic variation in climate change vulnerability assessment reduces species range loss projections. *Proceedings of the National Academy of Sciences of the United States of America*, 116(21), 10418-10423. https://doi.org/10.1073/pnas.18206633116

Rolland, J., Lavergne, S., & Manel, S. (2015). Combining niche modelling and landscape genetics to study local adaptation: A novel approach illustrated using alpine plants. *Perspectives in Plant Ecology, Evolution and Systematics*, 17(6), 491-499. https://doi.org/10.1016/j. ppees.2015.07.005

Root, T. L., Price, J. T., Hall, K. R., Schneider, S. H., Rosenzweig, C., & Pounds, J. A. (2003). Fingerprints of global warming on wild animals
<figure><img src="image_18.png" /></figure>

and plants. *Nature*, 421(6918), 57–60. https://doi.org/10.1038/nature01333

Sacristán-Bajo, S., García-Fernández, A., Iriondo, J. M., & Lara-Romero, C. (2019). Transcriptome assembly and polymorphism detection in *Silene ciliata* (Caryophyllaceae). *Plant Genetic Resources*, 17(5), 452–455. https://doi.org/10.1017/S1479262119000157

Saupe, E. E., Barve, V., Myers, C. E., Soberón, J., Barve, N., Hensz, C. M., & Owens, H. L. (2012). Variation in niche and distribution model performance: The need for a priori assessment of key causal factors. *Ecological Modelling*, 237, 11–22. https://doi.org/10.1016/j.ecolmodel.2012.04.001

Schaberg, P. G., Hennon, P. E., D'amore, D. V., & Hawley, G. J. (2008). Influence of simulated snow cover on the cold tolerance and freezing injury of yellow-cedar seedlings. *Global Change Biology*, 14(6), 1282–1293. https://doi.org/10.1111/j.1365-2486.2008.01577.x

Searcy, C. A., & Shaffer, H. B. (2016). Do ecological niche models accurately identify climatic determinants of species ranges? *The American Naturalist*, 187(4), 423–435. https://doi.org/10.1086/685387

Serra-Diaz, J. M., Keenan, T. F., Ninyerola, M., Sabaté, S., Gracia, C., & Lloret, F. (2013). Geographical patterns of congruence and incongruence between correlative species distribution models and a process-based ecophysiological growth model. *Journal of Biogeography*, 40(10), 1928–1938. https://doi.org/10.1111/jbi.12142

Sexton, J. P., McIntyre, P. J., Angert, A. L., & Rice, K. J. (2009). Evolution and ecology of species range limits. *Annual Review of Ecology, Evolution, and Systematics*, 40(1), 415–436. https://doi.org/10.1146/annurev.ecolsys.110308.120317

Sinclair, B. J., Williams, C. M., & Terblanche, J. S. (2012). Variation in thermal performance among insect populations. *Physiological and Biochemical Zoology*, 85(6), 594–606. https://doi.org/10.1086/665388

Soberón, J. (2007). Grinnellian niches and geographic distributions of species. *Ecology Letters*, 10(12), 1115–1123. https://doi.org/10.1111/j.1461-0248.2007.01107.x

Soley-Guardia, M., Radosavljevic, A., Rivera, J. L., & Anderson, R. P. (2014). The effect of spatially marginal localities in modelling species niches and distributions. *Journal of Biogeography*, 41(7), 1390–1401. https://doi.org/10.1111/jbi.12297

Soto-Correa, J. C., Saenz-Romero, C., Lindig-Cisneros, R., & de la Barrera, E. (2013). The neotropical shrub *Lupinus elegans*, from temperate forests, may not adapt to climate change. *Plant Biology*, 15(3), 607–610. https://doi.org/10.1111/j.1438-8677.2012.00716.x

Soule, M. (1973). The epistasis cycle: A theory of marginal populations. *Annual Review of Ecology and Systematics*, 4(1), 165–187. https://doi.org/10.1146/annurev.es.04.110173.001121

Swab, R. M., Regan, H. M., Matthies, D., Becker, U., & Bruun, H. H. (2015). The role of demography, intra-species variation, and species distribution models in species' projections under climate change. *Ecography*, 38(3), 221–230. https://doi.org/10.1111/ecog.00585

Therneau, T. M. (2018). coxe: Mixed effects cox models. R Package Version, 2, 2–7.

Thuiller, W., Albert, C. H., Dubuis, A., Randin, C., & Guisan, A. (2010). Variation in habitat suitability does not always relate to variation in species' plant functional traits. *Biology Letters*, 6(1), 120–123. https://doi.org/10.1098/rsbl.2009.0669

Title, P. O., & Bemmels, J. B. (2018). ENVIREM: An expanded set of bioclimatic and topographic variables increases flexibility and improves performance of ecological niche modeling. *Ecography*, 41(2), 291–307. https://doi.org/10.1111/ecog.02880

Tutin, T. G., Heywood, V. H., Burges, N. A., Valentine, D. H., Walters, S. M., & Webb, D. A. (1964). *Flora europaea*. Cambridge University Press.

Valladares, F., Matesanz, S., Guilhaumon, F., Araújo, M. B., Balaguer, L., Benito-Garzón, M., Cornwell, W., Gianoli, E., Kleunen, M., Naya, D. E., Nicotra, A. B., Poorter, H., & Zavala, M. A. (2014). The effects of phenotypic plasticity and local adaptation on forecasts of species

range shifts under climate change. *Ecology Letters*, 17(11), 1351–1364. https://doi.org/10.1111/ele.12348

Villellas, J., Ehren, J., Olesen, J. M., Braza, R., & Garc, B. (2013). Plant performance in central and northern peripheral populations of the widespread *Plantago coronopus*. *Ecography*, 36(2), 136–145. https://doi.org/10.1111/j.1600-0587.2012.07425.x

Volis, S. (2011). Adaptive genetic differentiation in a predominantly self-pollinating species analyzed by transplanting into natural environment, crossbreeding and $Q_{ST}-F_{ST}$ test. *New Phytologist*, 192(1), 237–248. https://doi.org/10.1111/j.1469-8137.2011.03799.x

Wagner, V., von Wehrden, H., Wesche, K., Fedulin, A., Sidorova, T., & Hensen, I. (2011). Similar performance in central and range-edge populations of a Eurasian steppe grass under different climate and soil pH regimes. *Ecography*, 34, 498–506. https://doi.org/10.1111/j.1600-0587.2010.06658.x

Walther, G.-R., Post, E., Convey, P., Menzel, A., Parmesan, C., Beebee, T. J. C., Fromentin, J.-M., Hoegh-Guldberg, O., & Bairlein, F. (2002). Ecological responses to recent climate change. *Nature*, 416(6879), 389–395. https://doi.org/10.1038/416389a

Warren, D. L., & Seifert, S. N. (2011). Ecological niche modeling in Maxent: The importance of model complexity and the performance of model selection criteria. *Ecological Applications*, 21(2), 335–342. https://doi.org/10.1890/10-1171.1

Wheeler, J. A., Hoch, G., Cortés, A. J., Sedlacek, J., Wipf, S., & Rixen, C. (2014). Increased spring freezing vulnerability for alpine shrubs under early snowmelt. *Oecologia*, 175(1), 219–229. https://doi.org/10.1007/s00442-013-2872-8

Williams, C. M., Marshall, K. E., MacMillan, H. A., Dzurisin, J. D. K., Hellmann, J. J., & Sinclair, B. J. (2012). Thermal variability increases the impact of autumnal warming and drives metabolic depression in an overwintering butterfly. *PLoS One*, 7(3), https://doi.org/10.1371/journal.pone.0034470

Wittmann, M. E., Barnes, M. A., Jerde, C. L., Jones, L. A., & Lodge, D. M. (2016). Confronting species distribution model predictions with species functional traits. *Ecology and Evolution*, 6(4), 873–879. https://doi.org/10.1002/ece3.1898

Zakharov, E. V., & Hellmann, J. J. (2008). Genetic differentiation across a latitudinal gradient in two co-occurring butterfly species: Revealing population differences in a context of climate change. *Molecular Ecology*, 17(1), 189–208. https://doi.org/10.1111/j.1365-294X.2007.03488.x

Zandalinas, S. I., Mittler, R., Balfagón, D., Arbona, V., & Gómez-Cadenas, A. (2018). Plant adaptations to the combination of drought and high temperatures. *Physiologia Plantarum*, 162(1), 2–12. https://doi.org/10.1111/ppl.12540

## SUPPORTING INFORMATION

Additional supporting information may be found in the online version of the article at the publisher's website.

**How to cite this article:** Morente-López, J., Kass, J. M., Lara-Romero, C., Serra-Diaz, J. M., Soto-Correa, J. C., Anderson, R. P., & Iriondo, J. M. (2022). Linking ecological niche models and common garden experiments to predict phenotypic differentiation in stressful environments: Assessing the adaptive value of marginal populations in an alpine plant. *Global Change Biology*, 28, 4143–4162. https://doi.org/10.1111/gcb.16181