Discussion
Our results point to a discordance between the absence of genomic
divergence and consistent phenotypic differentiation in traits whose
adaptive value seems obvious (Iraeta et al. 2006, 2010, 2011, 2013). The
fact that, as previously found, in the montane population the load of
ectoparasites is related to secondary sexual characteristics of males,
and that there are no infestations in the lowland population, shows that
the presence of ticks is associated with different selective pressures
in both populations (Llanos-Garrido et al. 2017). Apart from their
negative effect on the survival of individuals (Salvador et al. 1996),
ticks appear to have other significant impacts on fitness, influencing
the expression of sexual ornaments that are important for reproduction
(Salvador et al. 1996, Llanos-Garrido et al. 2017). In addition, the
altitudinal gradient involves many other selective pressures that, while
less easy to identify than the one imposed by ticks (but still evident
given the environmental differences shown by our PCA scores), generate
an undeniable impact on fitness-related traits such as growth rate, egg
size or clutch size (Iraeta et al. 2006, 2010, 2011, 2013). In addition,
the adaptive phenotypic differences found in this and earlier studies
span the whole life cycle of lizards (e.g. incubation time, hatchling
size, juvenile growth rate, and sexual ornamentation and reproductive
investment of adults), which should involve a considerable variety of
divergent selective pressures underlying their phenotypic
differentiation (Roff 2007). Moreover, the fact that these traits
maintain their differentiation in common garden conditions (Iraeta et
al. 2013) or in reciprocal transplant experiments (Iraeta et al. 2006)
reveals a divergent genetic basis between the two populations (De Kort
et al. 2014).
However, the background genetic differentiation between these two
populations is very low (Díaz et al. 2017, Verdú-Ricoy et al. 2010),
even when analytical power is increased with genomic scans. Irrespective
of phenotypic differenciation, other populations of Psammodromus
algirus showed marked genetic isolation despite being much closer (both
geographycally and environmentally) than the ones presented here (Díaz
et al. 2017, Llanos-Garrido et al. 2019). This is true even at a
fragmented metapopulation where suitable habitat patches are as close as
hundreds of meters apart, and where subpopulations have diverged less
than 100 years ago (Santos et al. 2008, Pérez-Tris et al. 2019).
However, this latter scenario occurs at the northern edge of the
species’ range, where genetic diversity might be lower than at El Pardo
or Navacerrada, which occupy the center of the Iberian distribution
range (Holt and Keitt 2000, Takahashi et al. 2016). Such decrease in
genetic diversity and/or effective population size could have led to a
subsequent reduction of migration rate and increased genetic structure
in this marginal population (Holt and Keitt 2000, Zalewski et al. 2016,
Li et al. 2018; see Deng et al. 2020), whereas the opposite might be
true at the center of the species’ distribution range. However, as
migration rate between El Pardo and Navacerrada is high enough (or
divergence time recent enough) to blur any hint of genetic structure, it
is nearly impossible to properly estimate gene flow among undefined
populations (Pfenninger and Posada, 2002). Thus, elucidating how gene
flow affects local adaptation in this system constitutes a distinct
challenge, really difficult to circumvent even by increasing our
analytical power or by getting longer assemblages.
However, the lack of overall genetic structure should not hinder the
genetic differentiation of specific variants on which the basis of
adaptive divergence is located. There are other systems in which the
background genetic differentiation between phenotypically
well-differentiated groups is low, and yet islands of selection whose
local divergence outstands above the rest of the genome can still be
found (Aguillon et al. 2018). These genome fractions may persist
differentiated even in continuous gene flow scenarios between the study
groups (Philips et al. 2008, Moody et al. 2015, Poelstra et al. 2014,
Shaner et al. 2015). In our case, it seems that the divergence between
populations is so recent that it has not translated into a genomic
imprint to a scale large enough to be detected with the approach we have
used (Aguillon et al al. 2018). However, there is compelling evidence
that points to the existence of at least a part of the genome that has
diverged, giving rise to the observed phenotypic differences. One
possible explanation is that the only divergent regions are those
expressed as phenotypes (Safran et al. 2016). These regions would be
located in very specific areas of the genome and could be really scarce
(Campagna et al. 2015), to the extent that low coverage genotyping
approaches could miss them. Another possibility is that the SNPs
responsible for the observed adaptive divergence could be located in
regions of the genome with low genetic diversity, which would express a
peak of divergence at a local genomic scale, but which would escape the
outliers detection analysis used in this study (Campbell et al 2018) if
the amount of genetic variability is distributed heterogeneously
throughout the genome (Feulner et al 2015, Nosil et al. 2009, Poelstra
et al. 2014). However, such heterogeneity does not seem to be an issue
in our system, because background divergence is very low across all the
genome, which suggests that SNPs with almost any degree of divergence
could have been detected as outliers.
Finally, there could be no significant divergence at all on a
genome-wide scale. It is true that previous studies point to a genetic
basis, especially those that included reciprocal transplants or common
garden experiments, but it could be the case that the observed
differences were the result of a heritable epigenetic response to the
environment (Groot et al. 2018, Trerotola et al. 2015). Such epigenetic
marks could be inherited only up to a specific number of generations,
after which the signal would be lost (Trerotola et al. 2015). To discard
this hypothesis, we could confirm previous findings by extending common
garden experiments a few generations more (Groot et al 2018), or we
could perform a transcriptomic differentiation study (McGirr y Martin
2018).
In summary, we succeeded to demonstrate that genome-wide differentiation
was not associated with local adaptation in an altitudinal gradient with
marked environment differences between two lizard populations. The
phenotypic differentiation pattern of these populations was obviously
concordant with a scenario of ecologically-driven divergence, but it did
not lead to any detectable genotypic isolation between them. Moreover,
and for reasons that remain unveiled, genome-wide differenciation seems
much lower than expected by geographic distance. Thus, we uncovered
another example of locally (and divergently) adapted populations which
remain undifferentiated on a genome-wide scale, a phenomenon that could
be more frequent than previously thought due to publication bias against
this apparent lack of positive results (Hendry 2009). Highlighting
natural systems where divergent adaptation occurs with different degrees
of genomic-wide differenciation (including none) is of paramount
importance to understand how frequently isolation by adaptation occurs,
and to what extent the first incipient stages of ecological speciation
are widespread (Feder et al. 2012, Hendry 2009, Shafer and Wolf 2013,
Sexton et al. 2014, Krohn et al. 2018).