The evolutionary pathways for local adaptation in mountain hares

Abstract

Understanding the evolution of local adaptations is a central aim of evolutionary biology and key for the identification of unique populations and lineages of conservation relevance. By combining RAD sequencing and whole-genome sequencing, we identify genetic signatures of local adaptation in mountain hares (Lepus timidus) from isolated and distinctive habitats of its wide distribution: Ireland, the Alps and Fennoscandia. Demographic modelling suggested that the split of these mountain hares occurred around 20 thousand years ago, providing the opportunity to study adaptive evolution over a short timescale. Using genome-wide scans, we identified signatures of extreme differentiation among hares from distinct geographic areas that overlap with area-specific selective sweeps, suggesting targets for local adaptation. Several identified candidate genes are associated with traits related to the uniqueness of the different environments inhabited by the three groups of mountain hares, including coat colour, ability to live at high altitudes and variation in body size. In Irish mountain hares, a variant of ASIP, a gene previously implicated in introgression-driven winter coat colour variation in mountain and snowshoe hares (L. americanus), may underlie brown winter coats, reinforcing the repeated nature of evolution at ASIP moulding adaptive seasonal colouration. Comparative genomic analyses across several hare species suggested that mountain hares' adaptive variants appear predominantly species-specific. However, using coalescent simulations, we also show instances where the candidate adaptive variants have been introduced via introgressive hybridization. Our study shows that standing adaptive variation, including that introgressed from other species, was a crucial component of the post-glacial dynamics of species.

Keywords: Lepus timidus; adaptive evolution; genome scans; introgression.

FIGURE 1

Mountain hare sampling and genetic relationships. (a) Distribution of the mountain hare in western Europe (black), according to Mitchell‐Jones et al. (1999), and sampled regions: Ireland, IRL; Alps, ALP; Fennoscandia, FSC. (b) Principal component analysis (PCA) based on whole‐genome data, inferred using 125,196 polymorphic sites, sampled at least 20 kb apart along the genome. (c) Neighbour‐joining tree based on individual pairwise genetic distances estimated from 128,019 SNPs derived from whole‐genome datasets. Nodes with bootstrap support values ≥0.90 are shown as black circles. (d) mtDNA phylogeny. Ages of key nodes are displayed (mean and 95% HDI) and nodes with posterior probability ≥0.90 are shown as black circles

FIGURE 2

Demographic history of the mountain hare populations. Models were inferred from RAD‐seq data (see detailed inferenced and 95% high density intervals in Table S4). (a) The simultaneous split (STM) model. (b) The sequential split (DTM) model. Parameter estimates are shown: N _e, effective population size in number of diploid individuals (I, initial; ANC, ancestral); t, time of splits in years, considering 2 years per generation (Marboutin & Peroux, 1995). The two models showed similar fit to the data with ΔLhood of 3870 and 3802 for the STM and DTM model, respectively. Ireland, IRL; Alps, ALP; Fennoscandia, FSC

FIGURE 3

Population branch statistic and selective sweep (DCMS) scans along the autosomes for each mountain hare population (A, Ireland; B, Alps; C, Fennoscandia). The peaks marked with black triangles represent regions with at least two consecutive windows with PBS p ≤ .01 (blue dashed horizontal lines indicate thresholds based on the STM model simulations, and black dashed horizontal lines indicate thresholds based on the DTM model simulation), confirmed by BayPass contrast analysis (see Figure S3). Codes indicate genes within the outlier regions (see Table S5). Dashed rectangles indicate PBS peaks that overlap with signatures of population‐specific selective sweeps, evidenced by DCMS analysis (see Figure S6 and Table S10). See Figure S5 for PBS scans along chromosome X

FIGURE 4

Maximum‐likelihood trees for the whole genome and local trees of candidate regions for local adaptation in the mountain hare (genes in the regions are noted, unless no annotation is known, in which case the coordinate is noted). (a) Global tree based on the concatenated whole genome. Local phylogenies of candidate regions in (b–d) Alpine mountain hares, (e–f) Irish mountain hares, and (g–h) Fennoscandian mountain hares. In the ASIP phylogeny, notations WW (winter‐white), WB (winter‐brown) and WG (winter‐grey), indicate winter colouration morphs of the sequenced individual

FIGURE 5

Signatures of introgression. Left panels: scan of d_XY (top; horizontal dotted line indicates 1% lowest values in interspecific d_XY simulated under a model with selection; IRL, Ireland; ALP, Alps; FSC, Fennoscandia; EUR, L. europaeus; CAS, L. castroviejoi) and RND between focal mountain hare and potential donor species (bottom; horizontal dotted line indicates 1% lowest values in the observed chromosome‐wide distribution), with vertical dashed lines representing the candidate regions. Right panels: distribution of d_XY between the focal mountain hare population and the putative donor species, including the observed distribution along the whole chromosome and the candidate region, and simulated distribution under the model with and without selection. Top horizontal bars represent 95% range for each distribution. (a–b) Irish L. timidus ‐ candidate region chr01:177.8 Mb and possible introgression from L. europaeus, EUR. (c–d) Alpine L. timidus ‐ candidate region chr08:101.7 Mb, overlapping PCCA gene, and possible introgression from L. castroviejoi, CAS