INTRODUCTION
In the near future, global climate change will surpass habitat
destruction as the leading threat to biodiversity (Leadley et al.,
2010). Given such changes, the ability of species to persist will rely
on either redistribution, evolution towards new behavioural and
physiological optima, or both. While a number of range redistributions
have already been observed (e.g., Laliberte & Ripple, 2004; Parmesan,
2006), the few species that are capable of such large-scale dispersal to
track changing environments are unlikely to redistribute to landscapes
containing their complete suite of optimal habitat conditions. For
example, where species redistribute to maintain equilibrium with one
habitat feature (e.g., prey abundance or distribution, temperature),
they might simultaneously become discordant with other critical features
(e.g., novel predators, photoperiod). Thus, it is likely that for most
species, some extent of microevolution will be necessary to allow for
their persistence in a changing climate (Visser, 2008).
With the advent of novel molecular technologies and developments in
bioinformatics, acquisition and analyses of large-scale data has made
the investigation of genomic elements of adaptation possible for both
model and non-model species (Harrison, Pavlova, Telonis-Scott, &
Sunnucks, 2014). In wild populations, adaptive genomic change will be
dependent on the level of standing genetic variation or the rate ofde novo mutations in the existing population (Carlson,
Cunningham, & Westley, 2014). While standing genetic variation is
considered most important for contemporary adaptation as it is readily
available for selection (Barrett & Schluter, 2008), genetic variation
in ecologically important traits may be absent in widely dispersed
species, making them targets for extirpation. Specifically, panmictic
species experiencing high rates of dispersal and gene flow may be
limited in the extent to which adaptive genetic mutations can accumulate
and contribute to adaptation. However, while gene flow is predicted to
preclude the local adaptation of wide-ranging species (Lenormand, 2002),
recent empirical evidence provides support for adaptive genetic
divergence in the face of ongoing gene flow (e.g., André et al., 2011;
Feder, Egan, & Nosil, 2012; Hemmer-Hansen, Nielsen, Frydenberg, &
Loeschcke, 2007; Nielsen et al., 2009). This is often attributed to
adaptations along ecological clines, including clines in salinity
(Defaveri, Johnsson, & Merilӓ, 2013), temperature (Aitken, Yeaman,
Holliday, Wang, & Curtis-McLane, 2008), or altitude, (Keller, Taverna,
& Seehausen, 2011). Characterizing this diversity across species
distributions can allow us to better understand the “adaptive
potential” of species, and contribute to conservation plans where
needed (Keller et al., 2011; Meröndun, Murray, & Shafer, 2019).
In lieu of high standing genetic variability, some species will rely
largely on de novo mutations for adaptation, a mechanism that may
not be effective over short time-frames. Specifically, rates of mutation
are likely to be too low (Fan & Chu, 2007), most new mutations will be
either neutral or deleterious rather than providing an adaptive benefit
(Eyre-Walker & Keightley, 2007), and in small or declining populations,
deleterious mutations can accumulate (Willi, Griffin, & Van Buskirk,
2012). Further, the fate of new mutations depends on generation time and
population size, where species with longer generation times and smaller
population sizes have fewer mutational opportunities (Carlson et al.,
2014). In cases where beneficial mutations do arise, they will have to
persist until they reach high frequencies, during which time the chances
of losing beneficial alleles to genetic drift or demographic
stochasticity can be high (Whitlock, 2000).
Considering the importance of rapid adaptation, coding trinucleotide
repeats (cTNRs) may be important targets of selection due to their
containment within exonic regions of the genome and elevated mutation
rates in comparison to more conventionally studied markers (ranging from
10-6 to 10-2 for microsatellites
versus 10-9 for SNPs; Fan & Chu, 2007). These repeats
are found within genes that code for biologically relevant functions,
and may facilitate the fine-tuning of gene function and expression
patterns (Gemayel, Vinces, Legendre, & Verstrepen, 2010; Gemayel, Cho,
Boeynaems, & Verstrepen, 2012; Haerty & Golding, 2010). Most
interestingly, however, is the high mutability of cTNRs, which has been
attributed to their inherent instability and may allow for higher rates
of evolution in response to environmental stress (Gemayel et al., 2010,
2012). For example, the presence of variable cTNRs in genes coding for
body morphology has implicated these repeats as the driver of the
rapidly evolved canid skeleton in an evolutionarily short time span of
fifty years (Laidlaw et al., 2007). This provides compelling evidence
for the role of cTNRs in the rapid morphological adaptation of large
mammals.
The candidate gene approach to identifying loci subjected to natural
selection consists of the selection of previously characterized genes
that are hypothesized to play a role in local adaptation, or genes that
may be biologically influenced by environmental heterogeneity (Rellstab,
Gugerli, Eckert, Hancock, & Holderegger, 2015). This type of approach
is especially appealing for non-model species with uncharacterized
genomes, as genomic information from closely related model species
(e.g., homologous genes) can be used to select genes of known function
and design primers for applications in non-model species of interest
(Primmer, Papakostas, Leder, Davis, & Ragan, 2013). Further, the
candidate gene approach is beneficial for traits that are not
phenotypically apparent but suspected to be of adaptive significance for
adaptation to climate change (Keller et al., 2011). The environmental
heterogeneity associated with latitude implies that adaptation should
involve a range of traits in species spanning latitudinal gradients. For
example, adaptation to varying temperatures and photoperiods, which are
often used to time life-history events (e.g., breeding and migration),
should be of particular importance to many species (O’Malley, Camara, &
Banks, 2007). Thus, adaptation of these traits will be critical for
species persistence in changing environments.
Although species can track shifts in temperature via redistribution,
exposure to novel photoperiods will likely result in challenges for
species relying on photoperiodic cues, as photoperiod-dependent
behaviours may not change at the same rate as the change in photoperiod
that occurs during a range shift, resulting in photoperiodic mismatches
(Milligan, Holt, & Lloyd, 2009). As a result, it is possible that the
rate of evolution in traits controlling photoperiod-dependent behaviours
will be more important than traits involved in temperature sensitivity
(Bradshaw & Holzapfel, 2006, 2008), and the ability of species to adapt
their behavioural response to new photoperiods, or eliminate
photoperiod-dependent behaviours in favour of opportunism, will become
critical (Bronson, 2009). In these cases, microevolution of genes
responding to photoperiod is likely necessary to ensure species
persistence in the face of climate change (Bellard, Bertelsmeier,
Leadley, Thullier, & Courchamp, 2012). Species inhabiting seasonal
environments often respond to photoperiod cues via circadian clocks
(Goldman, 2001), making clock genes good candidates for characterizing
the potential genetic responses of species to shifting seasonal
conditions and novel environments (Kondratova, Dubrovsky, Antoch, &
Kondratov, 2010). cTNRs have been observed in a number of clock genes,
and emerging studies have begun to demonstrate the evolutionary and
adaptive importance of clock gene cTNRs in a range of species (e.g.,
Johnsen et al., 2007; Liedvogel, Szulkin, Knowles, Wood, & Sheldon,
2009; O’Malley, Ford, & Hard, 2010; Prentice et al., 2017a).
Our objective was to explore the relationship between space, the
environment, and a putatively adaptive cTNR within the candidate clock
gene, Nuclear Receptor Subfamily 1, Group D, Member 1 (NR1D1 ), in
Canada lynx (Lynx canadensis ). The NR1D1 gene is a nuclear
receptor that links circadian rhythms to the transcriptional control of
metabolic pathways, and functions as a core component of the circadian
clock. Its influence on numerous pathways, including memory and
learning, behaviour, immune function, metabolism, and mitochondrial
content and function (i.e., exercise capacity), suggest that this gene
likely plays a role in adaptation to new environments (Everett & Lazar,
2014). Most interestingly, the NR1D1 gene has demonstrated
involvement in circadian regulation of food entrainment (Delezie et al.,
2016), and cold tolerance (Gerhart-Hines et al., 2013). The cTNR within
the NR1D1 gene has been previously shown to exhibit signatures of
positive selection in Canada lynx (Prentice et al., 2017a). Lynx are a
suitable species for investigating selection on clock genes because they
show high gene flow at neutral genetic markers across their mainland
range (Row et al., 2012), and thus, spatial patterns in NR1D1 are likely
to be due to selection rather than neutral processes such as dispersal
and gene flow.
Here, we assessed the influence of spatial and environmental patterns on
the variation and distribution of NR1D1 cTNR alleles in Canada
lynx, while controlling for the influence of background population
structure using a dataset of neutral microsatellite markers. First, we
examine inter-specific differences in the distribution and frequencies
of NR1D1 alleles between Canada lynx, and the closely related
species, bobcat (Lynx rufus ) (Johnson & O’Brien, 1997). Although
the evolutionary history of Canada lynx and bobcat is unclear, the
largely allopatric distributions of these species suggest that they are
likely to be differentially adapted to their respective environments.
The lynx is largely distributed across the northern boreal, sub-boreal
and western montane forests of North America (Quinn & Parker, 1987),
occupying a narrower niche than expected if randomly distributed (Peers,
Thornton, & Murray, 2012). The bobcat is broadly distributed throughout
southern North America in the contiguous USA, and occupies a broader
niche than expected (Peers et al., 2012), suggesting that it is a
habitat generalist (Anderson & Lovallo, 2003). We hypothesized that the
differing habitats used by the two species will result in a divergence
in selection pressures experienced by lynx and bobcat, and predicted
that we should observe a corresponding divergence between lynx and
bobcat alleles at the NR1D1 cTNR locus [consistent with the
findings of Prentice et al., (2017a)]. We also characterized this
locus in introgressed lynx and bobcat individuals (Koen, Bowman, Lalor,
& Wilson, 2014) to explore the possibility of selection in the
lynx-bobcat contact zone. We hypothesized that selection will favour
lynx alleles in bobcats inhabiting lynx habitats (i.e., individuals at
the north end of the bobcat range). We contrasted genetic
differentiation in our candidate cTNR locus and neutral microsatellite
datasets to test for differences in pairwise population differentiation
at a putatively selected locus versus neutral loci. As differing
processes are responsible for influencing neutral versus functional
genetic variation (gene flow and genetic drift versus natural selection,
respectively), we hypothesized that estimates of differentiation at our
putatively selected cTNR locus would be outside of the range reported
for neutral loci, representing population pairs that were under either
more similar (NR1D1 FST < neutral
FST) or different (NR1D1 FST> neutral FST) selection pressures than
would be expected based on gene flow alone.. Further, we assessed the
association of environmental variables with the observed genetic
structure of the NR1D1 locus, by assessing relationships betweenNR1D1 alleles and environmental and spatial variables, while
controlling for background genetic structure with our neutral
microsatellite dataset. We hypothesized that environmental heterogeneity
is influencing selection at the NR1D1 locus, and predicted that
we would detect associations between environmental variables associated
with lynx habitat selection and NR1D1 alleles, after removing the
effects of background neutral genetic structure.