# Genome-wide SNP assessment of contemporary European red deer genetic structure highlights the distinction of peripheral populations and the main admixture zones in Europe

Juan Carranza, Javier Pérez-González, Gabriel Anaya, Menno de Jong, Camilla Broggini, Frank E Zachos, Allan D Mcdevitt, Magdalena Niedzialkowska, Maciej Sykut, Sándor Csányi, et al.

► To cite this version:

Juan Carranza, Javier Pérez-González, Gabriel Anaya, Menno de Jong, Camilla Broggini, et al. Genome-wide SNP assessment of contemporary European red deer genetic structure highlights the distinction of peripheral populations and the main admixture zones in Europe. Molecular Ecology, 2024, 33 (18), 10.1111/mec.17508. http://hal-05111892

HAL Id: hal-05111892

https://ofb.hal.science/hal-05111892v1

Submitted on 13 Jun 2025

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_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Genome-wide SNP assessment of contemporary European red deer genetic structure highlights the distinction of peripheral populations and the main admixture zones in Europe

Juan Carranza¹ ID

Javier Pérez-González²

Gabriel Anaya¹

Menno de Jong³

Camilla Broggini¹

Frank E. Zachos⁴,⁵,⁶,⁷

Allan D. McDevitt⁸

Magdalena Niedziałkowska⁹ ID

Maciej Sykut⁹

Sándor Csányi¹⁰

Norber Bleier¹¹

Lázló Csirke¹¹

Knut Røed¹²

Christine Saint-Andrieux¹³

Aurelie Barboiron¹³

Araceli Gort-Esteve¹⁴

Jordi Ruiz-Olmo¹⁵

Jose Manuel Seoane¹

Jose Antonio Godoy¹⁶ ID

Paweł Mackiewicz¹⁷

Eva de la Peña¹,¹⁸

Giovanni Vedel¹

S. Eryn McFarlane¹⁹,²⁰ ID

Josephine Pemberton²⁰

Alberto Membrillo¹

**Correspondence**

Juan Carranza, Wildlife Research Unit
(UIRCP), University of Cordoba, Cordoba
14071, Spain.

Email: jcarranza@uco.es

**Funding information**

Andalusian Government and ERDF of
European Union, Grant/Award Number:
UCO-1262830 and P18-RT-4642

Handling Editor: David Coltman

**Abstract**

Genome-wide technologies open up new possibilities to clarify questions on genetic structure and phylogeographic history of taxa previously studied with microsatellite loci and mitochondrial sequences. Here, we used 736 individual red deer (*Cervus elaphus*) samples genotyped at 35,701 single nucleotide polymorphism loci (SNPs) to assess the population structure of the species throughout Europe. The results identified 28 populations, with higher degrees of genetic distinction in peripheral compared to mainland populations. Iberian red deer show high genetic differentiation, with lineages in Western and Central Iberia maintaining their distinctiveness, which supports separate refugial ranges within Iberia along with little recent connection between Iberian and the remaining Western European populations. The Norwegian population exhibited the lowest variability and the largest allele frequency differences from mainland European populations, compatible with a history of bottlenecks and drift during post-glacial colonization from southern refugia. Scottish populations showed high genetic distance from the mainland but high levels of diversity. Hybrid zones were found between Eastern and Western European lineages in Central Europe as well as in the Pyrenees, where red deer from France are in close contact with Iberian red deer. Anthropogenic restocking has promoted the Pyrenean contact zone, admixture events in populations on the Isle of Rum and in the Netherlands, and at least partly the admixture of the two main lineages in central-eastern Europe. Our analysis

Juan Carranza and Javier Pérez-González are first authors, equal contribution.

For affiliations refer to page 16.

This is an open access article under the terms of the [Creative Commons Attribution License](https://creativecommons.org/licenses/by/4.0/), which permits use, distribution and reproduction in any medium, provided the original work is properly cited.

© 2024 The Author(s). *Molecular Ecology* published by John Wiley & Sons Ltd.
enabled detailed resolution of population structure of a large mammal widely distrib-
uted throughout Europe and contributes to resolving the evolutionary history, which
can also inform conservation and management policies.

KEYWORDS

Cervus elaphus, glacial lineages, phylogeography, population genomics, ungulates, wildlife management

**1 | INTRODUCTION**

Our understanding of the factors contributing to the composition of
fauna across extensive regions has significantly advanced through
the analysis of the genetic structure of species with wide geographic
distribution ranges (De Kort et al., 2021). The field of population ge-
nomics is assuming a progressively pivotal role in elucidating both
natural and anthropogenic influences on the evolutionary trajecto-
ries of such widely distributed species. Additionally, it serves as a
crucial tool for formulating effective conservation and management
strategies (Hohenlohe et al., 2021).

The red deer (*Cervus elaphus*) is one of the most widespread
wild large mammals in the Palearctic, still present in most of its nat-
ural range from the Pleistocene (Mattioli et al., 2022). As such, the
species has experienced the influence of glacial pulses that have
defined its current phylogeographic patterns (Doan et al., 2017,
2022; Krojerová-Prokešová et al., 2015; Mahmut et al., 2002; Meiri
et al., 2013; Niedziałkowska et al., 2011; Skog et al., 2009; Zachos &
Hartl, 2011). It has been used as a model species to address questions
related to population genetic structure shaped by anthropogenic fac-
tors (Zachos & Hartl, 2011) and to investigate the location of glacial
refugia throughout its wide geographical distribution area and their
influence on post-glacial expansions (Doan et al., 2022; Krojerová-
Prokešová et al., 2015; Meiri et al., 2013; Niedziałkowska, Doan,
et al., 2021; Queirós et al., 2019; Skog et al., 2009). Due to its overall
importance and wide geographic distribution, it can help to answer
global and general questions, including testing phylogeographic pat-
terns of temperate species in Europe (e.g. the expansion-contrac-
tion model; Hewitt, 1996, 1999; Taberlet et al., 1998), refugial areas
(Schmitt, 2007) and hybridization events (Allendorf et al., 2001).

Studies addressing the phylogeography, genetic structure, con-
servation and evolutionary history of European red deer have pre-
viously used mitochondrial DNA sequences (e.g. Doan et al., 2022;
Ludt et al., 2004; Mackiewicz et al., 2022; Meiri et al., 2013,
2018; Niedziałkowska et al., 2011; Rey-Iglesia et al., 2017; Skog
et al., 2009), neutral microsatellite markers (e.g. Hoffmann
et al., 2016; Pérez-González et al., 2012, 2023; Zachos et al., 2016) or
both (e.g. Carranza et al., 2016; Niedziałkowska et al., 2012; Queirós
et al., 2019). The building of large SNP datasets, not limited to neu-
tral loci and distributed throughout the whole genome, is expected
to increase the precision and power in the estimation of population
genetic parameters and to more fully elucidate phylogeographic and

evolutionary patterns (Allendorf et al., 2010; Brumfield et al., 2003;
Väli et al., 2008). In terms of studies on phylogeographic patterns
in European mammals, recent large-scale investigations on widely
distributed species such as the red fox (Vulpes vulpes) and bank vole
(Myodes glareolus) using thousands of SNPs have clarified outstand-
ing questions on post-glacial re-colonization patterns (Marková
et al., 2020; McDevitt et al., 2022; Niedziałkowska et al., 2023). All
these studies show that, at a broad spatial scale, there can be high
concordance between genome-wide approaches and traditional ge-
netic markers (Vaux et al., 2023). However, the use of a high num-
ber of SNPs may grant higher precision in the study of population
genetic substructure (Marková et al., 2020; McDevitt et al., 2022;
Niedziałkowska et al., 2023). Therefore, new assessments with
genome-wide tools in previously studied systems can help to im-
prove the knowledge of genetic structure and the underlying histor-
ical processes along with enhancing genetically based decisions on
conservation and management of wildlife.

The development of a medium-density genotyping chip for
red deer (Brauning et al., 2015), even with reference to the bovine
genome sequence, has produced new insights into the genetic ar-
chitecture of red deer (e.g. Johnston et al., 2017). Then, with the
publication of the red deer reference genome and the update of
its recent alignment version (Pemberton et al., 2021), some papers
using genome-wide SNP markers have studied variability in pheno-
typic traits such as antlers (Elblinger et al., 2022; Peters et al., 2022),
birth weight (Gauzere et al., 2023), the underlying processes leading
to runs of homozygosity (ROH) (Hewett et al., 2023) and the genetic
structure of red deer populations including following anthropogenic
disturbances (de Jong et al., 2020).

However, no study so far has addressed the current genetic
structure of red deer at the European level by using genome-wide
SNP markers. Some questions posed by traditional genetic tools on
the genetic structure of red deer in Europe still deserve deeper as-
sessment. Among them are (a) the role of Iberia and other potential
refugial areas in support of alternative phylogeographic models; (b)
current differentiation between populations that may have shared
glacial refugia and experienced significant changes in their distribu-
tion ranges and contact zones throughout the European continent,
in particular, why some populations such as Norwegian or Scottish
ones differ more than others from the central European ones; (c)
current genetic diversity and inter-population admixture, in particu-
lar, where have main hybridization events have occurred in Europe,
and whether they have been natural or human mediated; (d) the maintenance of genetic distances between populations belonging to different refugial lineages that have subsequently come into contact in eastern Europe and finally (e) whether there are discrepancies between genome-wide SNPs and other genetic markers in some of these topics, in other words, elucidate what SNP genomics can contribute to the story.

In this study, we used a genome-wide approach with more than 35,000 SNPs to address these above-mentioned main points on the contemporary genetic structure of European red deer and to infer potential underlying processes.

2 | MATERIALS AND METHODS

2.1 | Study areas, sampling collection and DNA
purification

The samples collected in the present study (N=498) were from
current red deer ranges in Spain (Sierra de San Pedro: SP1 and SP2
populations, Monfragüe National Park: MO, Sierra Morena: SM1
and SM2 and the Pyrenees: PYS, PYN, PYE), Scotland (SC1 and SC2),
France (le Petite Pierre: PP), Norway (NO), Italy (Tarviso: TA), Austria
(Seefeld in Tyrol: AU), Hungary (North: HN, West: HW and South/
Gemenc: HS) and Poland (Poland North: PN, Poland South: PS), after

hunting activities or culling during ongoing management efforts. In
addition to these collected samples, we incorporated available geno-
types from the Netherlands (OVP, PWB, VZM, HGV, KRO), Scotland
(RU, KY1, KY2) and Czechia (CZ) (N=246; see below). In total, we
used information from 744 individuals from 10 European countries
(Figure 1, Table 1).

In Spain, SP1 red deer population corresponds to SP1 population in Pérez-González et al. (2012); SP2 in this study corresponds to SP2 in Pérez-González et al. (2012) and PNB in Queirós et al. (2019); SM1 and SM2 populations in this study correspond to SMR and MTR in Queirós et al. (2019; see below) respectively; MO in this study corresponds to PNM in Queirós et al. (2019).

We focused our study and the sampling effort on some of the
main areas of the European red deer distribution, and thereafter,
due to information provided by previous studies and our own re-
sults (see Figure 2, below in Results), we identified populations
within these initial areas resulting in 28 populations being ana-
lysed throughout Europe (Figure 1, Table 1). Samples from the
Pyrenees were divided into three populations based on a previous
study with microsatellite markers (Pérez-González et al., 2023):
South Pyrenees (PYS; individuals translocated in the 1970s from
Central Spain: SM2 population in this study, see below), North
Pyrenees (PYN; migrant individuals that appeared in the 1980s
and came from populations located in southern France that, in
turn, were introduced from other French areas further north) and

<figure><img src="image_2.png" /><figcaption>**FIGURE 1** Populations throughout European countries resulting from the sampled areas and the genetic structure depicted in **Figure 2**. Location codes as in **Table 1**.</figcaption></figure>
<table><caption>TABLE 1 Population names, codes, sample sizes (number of individuals) and sex ratio in our samples (females/males) for the geographical populations (see Figure 1) as resulting from the sampled areas and genetic structure depicted in Figure 2.</caption><thead><tr><th>Population name</th><th>Code</th><th>Sample size</th><th>Sex ratio (f/m)</th></tr></thead><tbody><tr><td>Spain Sierra de San Pedro 1<sup>a</sup></td><td>SP1</td><td>61</td><td>1.346</td></tr><tr><td>Spain Sierra de San Pedro 2<sup>a</sup></td><td>SP2</td><td>7</td><td>0.400</td></tr><tr><td>Spain Sierra Morena 1<sup>a</sup></td><td>SM1</td><td>75</td><td>0.442</td></tr><tr><td>Spain Sierra Morena 2<sup>a</sup></td><td>SM2</td><td>24</td><td>0.412</td></tr><tr><td>Spain Monfragüe<sup>a</sup></td><td>MO</td><td>19</td><td>0.357</td></tr><tr><td>Spain Pyrenees South<sup>a</sup></td><td>PYS</td><td>27</td><td>1.250</td></tr><tr><td>Spain Pyrenees East<sup>a</sup></td><td>PYE</td><td>59</td><td>1.565</td></tr><tr><td>Spain Pyrenees North<sup>a</sup></td><td>PYN</td><td>38</td><td>1.375</td></tr><tr><td>Scotland Rum</td><td>RU</td><td>100</td><td>0.492</td></tr><tr><td>Scotland Kintyre 1</td><td>KY1</td><td>40</td><td>0.739</td></tr><tr><td>Scotland Kintyre 2</td><td>KY2</td><td>43</td><td>0.636</td></tr><tr><td>Scotland 1<sup>a</sup></td><td>SC1</td><td>10</td><td>0.000</td></tr><tr><td>Scotland 2<sup>a</sup></td><td>SC2</td><td>10</td><td>0.111</td></tr><tr><td>Norway<sup>a</sup></td><td>NO</td><td>34</td><td>0.478</td></tr><tr><td>Netherlands North OVP</td><td>OVP</td><td>15</td><td>0.875</td></tr><tr><td>Netherlands South PWB</td><td>PWB</td><td>9</td><td>8.000</td></tr><tr><td>Netherlands South VZM</td><td>VZM</td><td>12</td><td>1.000</td></tr><tr><td>Netherlands South HGV</td><td>HGV</td><td>4</td><td>0.000</td></tr><tr><td>Netherlands Centre KRO</td><td>KRO</td><td>11</td><td>1.750</td></tr><tr><td>Czechia West</td><td>CZ</td><td>12</td><td>1.400</td></tr><tr><td>France PP<sup>a</sup></td><td>PP</td><td>25</td><td>1.273</td></tr><tr><td>Austria<sup>a</sup></td><td>AU</td><td>9</td><td>0.500</td></tr><tr><td>Poland North<sup>a</sup></td><td>PN</td><td>24</td><td>1.000</td></tr><tr><td>Poland South<sup>a</sup></td><td>PS</td><td>14</td><td>6.500</td></tr><tr><td>Hungary North<sup>a</sup></td><td>HN</td><td>13</td><td>1.600</td></tr><tr><td>Hungary West<sup>a</sup></td><td>HW</td><td>12</td><td>1.400</td></tr><tr><td>Hungary South<sup>a</sup></td><td>HS</td><td>18</td><td>0.500</td></tr><tr><td>Italy TA<sup>a</sup></td><td>TA</td><td>19</td><td>1.375</td></tr></tbody></table>

Note: Total sample size: 744 individuals from 28 populations.

East Pyrenees (PYE, a population whose individuals from North and South Pyrenees are in contact and interbreed).

For most populations (SP1, SP2, MO, SM1, SM2, PYN, PYS, PYE, SC1, SC2, NO, PP, TA, AU, HN, HW, HS, PN and PS), we collected tissue samples (e.g. fragments of ear, internal organs or muscles) from culled individuals, immediately kept at 4°C and then stored at -20°C until further use. Genomic DNA was isolated from pieces of tissue using the commercial DNA purification kit, DNeasy Blood and Tissue Kit (Qiagen, Germantown, MD, USA) according to the manufacturer's protocol.

Extracted DNA was visualized for integrity in a 0.8% agarose gel: a 3-µL aliquot of each sample was loaded with 2 µL of loading dye. Integrity of DNA samples was verified by comparing them against a

DNA Ladder Sibenzyme molecular weight marker (1 kb). After that, 2 µL of genomic DNA was measured by NanoDropTM One (Thermo Fisher Scientific Inc, Waltham, MA, USA). If the purified DNA did not have the required degree of purity (Absorbance ratios of A260/A280 and A260/230 of 1.8-2), a second extraction was carried out.

## 2.2 | SNP genotyping, dataset generation and quality control

The individuals were genotyped using the cervine 50K Illumina Infinium iSelect HD Custom BeadChip (Brauning et al., 2015; Kasarda et al., 2017) that includes 50,841 SNP markers. The chromosome location and SNP position were referenced according to the assembly and annotation of the CerEla1.1 genome version (Pemberton et al., 2021).

In addition to the genotypes obtained from samples collected in this study, we incorporated the genotypes of red deer from the Netherlands (OVP, PWB, VZM, HGV, KRO), Scotland (RU, KY2) and Czechia (CZ) published by de Jong et al. (2020) (Dryad repository: https://doi.org/10.5061/dryad.hqbzkh1dp; N=206; see Figure 1, Table 1). The genotypes from the Netherlands came from five populations (de Jong et al., 2020), renamed here: Netherlands North OVP (relocated individuals in the 1990s from Scotland and Czechia into the Dutch Nature Reserve Oostvaardersplassen), Netherlands South PWB and Netherlands South VZM (pure native populations from Planken Wambuis and Veluwezoom respectively) and Netherlands South HGV and Netherlands South KRO (private fenced estates in Hoge Veluwe and Kroondomeinen, respectively, where at the start of the 20th century the populations were supplemented with red deer from eastern Europe). Additionally, we incorporated into the dataset the genotypes (N=40) from another Scottish population located in the Kintyre Peninsula (KY1). These genotypes were published by McFarlane et al. (2020) (figshare repository: https://doi.org/10.6084/m9.figshare.9197849.v1).

Linkage disequilibrium (variance inflation factor of 1.25) was assessed with Plink (Purcell et al., 2007) and 8362 linked SNPs were removed. The remaining quality control analyses were conducted by the 'filterdata' option with the parameters indmiss=0.015, snpmiss=0.1 and min mac=1 in the R package SambaR (de Jong et al., 2021). We also compared the locus-specific minor allele frequency and locus-specific heterozygosity to identify paralogous loci (Figure S1). All loci passed the expected heterozygosity filter. However, X-chromosomal SNPs, which are all haploids for males (lower curve in Figure S1), were removed from the analyses, so we conducted all the remaining analyses with the autosomal SNPs only. Thus, after the whole selection process, 35,701 SNPs were retained. The SNPs were evenly distributed across chromosomes since there was high correlation between chromosome length and number of SNPs (quality control output of SambaR; Figure S2). Table S1 shows the mean number of SNPs per individual, grouped per population. We removed those samples with high proportion (>0.2) of missing data. Relatedness between sample pairs was estimated with the kinship coefficient (Waples

<footnote><sup>a</sup>Populations from which individuals were genotyped for this study.</footnote>
et al., 2019) and the KING-robust score (Manichaikul et al., 2010). Sample pairs with scores greater than 0.354 were considered duplicates (Manichaikul et al., 2010), so one individual for each of these pairs was excluded avoiding the possibility of analysing the same individual twice. After these quality control analyses, we retained 736 individuals that were used in the subsequent analyses.

## 2.3 | Analyses of genomic structure at the individual level

Genetic structure was assessed by using the genetic information at individual and population levels. At the individual level, we constructed an unrooted dendrogram depicting Euclidean genetic distances with the R package stats-4.2.1 (R Core Team, 2020). The genetic distances were generated with ordinary least squares (OLS) clustering applied to the Euclidean distances between individuals. In this dendrogram, the genetic distances between individuals were represented by the total path length between them. Moreover, to study the main genetic portioning among samples, principal coordinate analyses (PCoA) were performed using the function 'pcoa' of the R package ape-5.6.2 (Paradis & Schliep, 2019).

To depict possible admixed provenance of individuals, the pathlengths obtained with the OLS clustering were compared with the actual genetic distances between individuals (see de Jong et al., 2023). Elevated residual errors of the hierarchical clustering model evidence admixture origin of individuals. Admixture coefficients were calculated with the functions 'snmf' and 'Q' in the R package LEA-3.8.0 (Frichot & François, 2015). Alpha was set to 10, tolerance to 0.00001, and the number of iterations to 200.

We also estimated individual ancestry coefficients calculated with the software Admixture-1.3 (Alexander et al., 2009) and used the 'plotstructure'-function of SambaR for graphs. The optimal number of clusters (K) was determined using the elbow method on cross-entropy scores generated by the 'snmf' function, with the assumption that optimal K occurs when higher K values do not produce significant reductions of cross-entropy values. Genetic distances between individuals estimated as allele sharing distance were also calculated with SambaR.

We used the function 'calcdiversity' of the R package SambaR to quantify the genetic diversity of individuals according to individual heterozygosity. Inbreeding was estimated using the inbreeding coefficient (F), which is the probability that the two alleles at any locus of a diploid individual are identical by descent (IBD; Kardos et al., 2015).

Runs of homozygosity (ROH) were calculated with the function 'findROH' of the SambaR R package. ROHs were defined as regions with a stretch of at least 15 adjacent homozygous sites and a minimum length of 200kb. Only chromosomes with at least 500 sites were analysed. Given that the average heterozygosity per SNP for this dataset is high, at 0.275, the probability of a stretch of 15 adjacent homozygous sites equals $(1-0.275)^{15} \approx 0.008$. This means that around 99% of ROHs of 15 sites are likely to be caused by inbreeding. ROH were quantified as the number (# ROHs) and the

total proportion of the genome in runs of homozygosity per individual (F-roh). F-roh scores were compared with the number of ROHs and with the variation of F-roh scores across chromosomes to infer historical evolutionary processes of individuals. For instance, high F-roh values yet the low number of ROHs are indicative of recent inbreeding resulting in long runs of homozygosity. Recent inbreeding is also expected to leave a signature of high variation of F-roh scores across chromosomes.

## 2.4 | Analyses of genomic structure at the population level

Genetic information at the level of populations was used to complement information on genetic structure and to obtain further insight into the history of populations. Nei's genetic distances (Nei, 1972) between populations were used to generate unrooted and rooted dendrograms with the neighbour-joining (NJ) method calculated with the function 'stamppNeisD' of the R package StAMPP-1.6.3 (Pembbleton et al., 2013). As above, the total path length between populations represents the genetic distance between populations. Sambar was also used to obtain $F_{st}$ estimates between populations. Furthermore, the existence of admixture between populations was inferred by using the f3-scores obtained with the Treemix software (Pickrell & Pritchard, 2012). Low f3-scores, calculated as $(y-x) \cdot (y-z)$, indicate that allele frequencies y in population Y are intermediate between the allele frequencies x and z in populations X and Z, which indicates admixture. We represented graphically those populations for which the f3-scores were lower than 0.0015. Note that this threshold was arbitrary, but in our study system almost all f3-scores for all populations were positive (see below) and we aimed to represent those populations for which admixture was more probable, as well as for highlighting the donor populations of admixture. To determine current genetic similarity relationships between populations and possible evolutionary relationships, maximum likelihood phylogeny (based on allele frequency differences) and gene flow edges were generated with Treemix. We used the 'pure' Scottish sika deer (Cervus nippon) as an outgroup population (McFarlane et al., 2020). The optimal tree was selected with OptM (Fitak, 2021) that estimates the optimal reduced number of migration edges.

Initial geographical populations were regrouped based on the results at the individual level before conducting the analyses at the population level. In particular, our results at the individual level showed that two geographically defined populations (PYE and OVP) included individuals that belonged to two different genetic origins, evidenced by their position in the tree of Euclidean genetic distances (Figure 2). This splitting allowed us to better study the relationships between original populations and to uncover further events of admixture beyond those already evidenced. We coded these groups of individuals within populations (or subpopulations) with subscripts, that is, $PYE_1$, $PYE_2$, $OVP_1$ and $OVP_2$.

Demographic history was inferred by using GONE (Santiago et al., 2020) that estimates the effective population size $N_e$ for
<figure><img src="image_3.png" /><figcaption>FIGURE 2 Genetic structure at the individual level. (a) Dendrogram depicting Euclidean genetic distances between individuals. The diagram was generated with ordinary least squares (OLS) clustering applied to an input matrix of Euclidean distances between individuals. The genetic distances between individuals are represented by the total path length between individuals. (b) Principal coordinate analyses (PCoA). Scatterplot showing the first two ordination axes according to the PCoA applied to an input matrix of Euclidean distances between individuals. The first PC-axis has been mirrored to mimic geographical locations. (c) Residual errors of the OLS clustering model. The figure shows the heatmap on the right side depicting the discrepancy between the pathlengths in the dendrogram and the actual genetic distances. Red means attraction: The actual distances between individuals are smaller than suggested by the dendrogram. Blue means repulsion: The actual distances between individuals are greater than suggested by the dendrogram. Population codes as in Table 1, with subscripts in panel (a) for those subpopulations that split at different rooting positions in the dendrogram.</figcaption></figure>

generations back from the present on the basis of the relationship
between SNP pairs in linkage disequilibrium (LD) and distance in the
chromosome. GONE does not work properly when there is admixture
(Novo et al., 2023), so we applied it only to the least admixed popula-
tions, that is, Norway and the two main Iberian lineages (see below).
Input parameters were genetic distances between SNPs according to
the reference genome (Pemberton et al., 2021) and recombination rate
set to 1 centimorgan per megabase (cM/Mb; Johnston et al., 2017).

3 | RESULTS

## 3.1 | Genomic structure: Information at the individual level

The dendrogram depicting Euclidean genetic distances between in-
dividuals suggests the differentiation of groups of samples consist-
ent with the 28 populations (Figure 2a). There were considerable
genetic similarities among individuals within most populations, with the exceptions of eastern Pyrenees (PYE) and the Dutch population OVP that showed individuals clustering at different root positions in the dendrogram. These groups within the geographical populations are coded hereafter with subscripts (Figure 2a). Red deer sampled in Norway and in Spain are the most distant from the central node (Figure 2a) and from mainland central European populations (allele sharing distance between individuals mostly shows the genetic differentiation of the Norwegian population; Figure S3). Samples collected in Spain grouped into two distinct clusters: one containing the native populations of Sierra Morena (SM1 and SM2), San Pedro (SP1 and SP2) and Monfragüe (MO) along with the southernmost Pyrenean populations (PYS, restocked from central Spain), and the other one containing populations from northern and eastern Pyrenees (PYE). Within the first cluster, Spanish red deer from San Pedro grouped into two further clusters (SP1: eastern San Pedro; SP2: western San Pedro), while Spanish individuals from Sierra Morena fell into two groups (SM1: western Sierra Morena; SM2: eastern Sierra Morena). These subdivisions of the main Spanish populations were treated as different populations in all subsequent analyses. Regarding the remaining populations, the Norwegian individuals (NO) had short terminal branch lengths, indicating high genetic similarity. Individuals from Netherland North OVP (OVP<sub>1</sub>), a recently reintroduced population mainly stocked from Scotland, clustered with Scottish deer (more particularly with SC2). The Scottish SC2 population (and therefore also OVP<sub>1</sub>) form clusters situated between Scotland Rum (RU) and the other Scottish populations (KY1, KY2, SC1). A separate group appeared with three clusters close to the central node. One with the French (PP) population, another one containing samples from Czechia (CZ) and the Netherlands (OVP, PWB, VZM, HGV, KRO) and finally a third one with representatives from Italy (TA), Austria (AU), Hungary (HN, HS, HW) and Poland (PN, PS).

The genetic distances between individuals and the grouping of populations can be visualized in the principal coordinates analyses (PCoA) applied to an input matrix of Euclidean distances. Iberian, Norwegian and Scottish red deer were at the extreme edges/coordinates in the plot, while individuals from mainland Europe were closely grouped (Figure 2b). The first PC-axis (explaining 5.3% of the total variance) contrasts the native Iberian red deer with the other populations, in particular Norway (Figure S4). Populations from the Spanish areas of the northern and the eastern Pyrenees grouped with individuals from mainland Europe and some individuals (PE<sub>2</sub>) from eastern Pyrenees had intermediate coordinates. The second PC-axis (explaining the 3.3% of the total variance) contrasts Norwegian red deer with those from the Isle of Rum, Scotland. Because the explained variance depends on sample size, we repeated the PCoA analyses with normalized sample sizes (Figure S5). After downsampling the number of individuals to 35 for every population, Spanish and Norwegian red deer retained the same genetic differentiation, but Scottish red deer shared similar coordinates as individuals from mainland Europe in the PCoA plot (Figure S5). Figures S6 and S7 show the PCoA plots for mainland Europe (without Spanish, Norwegian and Scottish populations) and for native populations (without the Pyrenean, non-native Dutch and Rum populations) respectively. For the red deer from mainland Europe, the first PC-axis (explaining 6.5% of the variance) contrasts West-European populations with East-European populations. The second PC-axis (explaining the 4.8% of the variance) captures variation in West Europe, and contrasts Dutch populations with the French population. The third PC-axis captures variation in East Europe, and contrasts Hungary South with Poland North. On both the first and the second axis, the Dutch populations HGV and KRO clustered halfway between Czechia (CZ) and the other two Dutch populations, PWB and VZM. The PCoA plot for native populations shows three differentiated clusters: Spain, Norway and the remaining populations

<figure><img src="image_4.png" /><figcaption>FIGURE 3 Ancestry coefficients barplot. Figure shows coefficients for *K* values ranging between 2 and 7 ancestral population groupings, and for *K* = 10 (see Figure S8 for the cross-entropy scores associated with *K* values). Population codes as in Table 1.</figcaption></figure>
(Scotland and mainland Europe; **Figure S7**). Figures S5 and S7 support that, by using our SNP dataset, the Spanish and Norwegian red deer presented the highest levels of genetic differentiation among the studied populations.

The heatmap of the residual errors of OLS (Figure 2c) shows discrepancies with respect to the pathlengths depicted in the dendrogram, mostly for North and East Pyrenees, Rum, the Netherlands and eastern Europe, which may be indicative of admixed origins. Individuals from North Pyrenees showed similarities with other Spanish populations (mainly SM1). Some individuals from East Pyrenees also showed similarities with red deer from other Spanish populations, but higher genetic differentiation from Norway and eastern Europe. There is one sample in Spanish SM2 that was similar to Spanish SP1, which could be suggestive of recent gene flow. Close to the diagonal, some

populations showed mixed attraction and repulsion of some individuals to others within their own populations suggesting heterogeneous composition. This is especially the case for PYE, but also SP1, SP2, RU, OVP, NO and some of the Dutch populations except KRO, which was more attracted to eastern populations.

Figure 3 shows the ancestry coefficients for each individual using the Admixture-1.3 software. The cross-validation score (Figure S8) indicated that the optimal number of ancestral populations was *K* = 14 (i.e. there was no significant difference between *K* = 14 and higher *K* values), but also that from *K* = 10 onwards, there was little reduction of cross-entropy values. To visualize the main clusters of populations, we show here (Figure 3) the individual ancestry coefficients for *K*-values from 2 to 7 and for *K* = 10. Individuals from the native Spanish populations presented

<figure><img src="image_5.png" /><figcaption>FIGURE 4 Genomic diversity. (a) Individual genomic diversity measures. Left: Multi-locus heterozygosity; right: Inbreeding coefficient. (b) Runs of homozygosity per individual. Upper panel: Number of ROHs; lower panel: Proportion of the genome in ROHs, F-roh. (c) Comparison between measures related to ROH. Left: F-roh versus number of ROH; right: Variation (standard deviation) of F-roh across chromosomes versus F-roh. Population codes as in Table 1.</figcaption></figure>
high ancestry scores to an exclusive genetic cluster for low K values ($K=2$ to $K=4$). At $K=5$ and $K=7$, native Spanish red deer divided into two exclusive clusters (San Pedro[SP1, SP2] vs. Sierra Morena[SM1, SM2]/Monfragüe[MO]/South Pyrenees[PYS]). At $K=10$, individuals assigned to the previous Sierra Morena[SM1, SM2]/Monfragüe[MO]/South Pyrenees[PYS] cluster divided into two clusters (Sierra Morena[SM1] vs. Sierra Morena[SM2]/Monfragüe[MO]/South Pyrenees[PYS]). At $K=6$ and $K=7$, only Sierra Morena[SM2]/Monfragüe[MO]/South Pyrenees[PYS] appeared to share genetic variation with other European populations. At all K values, the figure shows the non-Spanish ancestry of north Pyrenean (PYN) individuals and that the eastern Pyrenees (PYE) comprised individuals with Spanish and non-Spanish ancestries. Moreover, Norwegian red deer (NO) formed an exclusive cluster at $K=4$ and higher K values. Scottish red deer from Rum (RU) appeared genetically different from the remaining Scottish red deer (KY1, KY2, SC1, SC2). Individuals from the Netherlands showed diverse ancestries. Finally, red deer from Poland South and Hungary shared genetic ancestry.

## 3.2 | Genomic diversity

Figure 4 shows the individual genomic diversity as multi-locus heterozygosity (MLH, Figure 4a, left) and inbreeding coefficient (Figure 4a, right). ROH analyses are shown in Figure 4b. Red deer from Norway presented the lowest estimates of genomic diversity. For Spanish samples, most individuals from San Pedro (SP1, SP2) and some individuals from the Pyrenees (mainly PYS and PYE) had relatively low genetic diversity (low MLH, high inbreeding coefficient and number of ROH), while Sierra Morena (SM1, SM2)

showed relatively high diversity for all the samples. Individuals from Scottish populations (mainly RU, as well as the Dutch population restocked from Scotland, OVP), Czechia (CZ), Austria (AU), Northern Poland (PN) and Italy (TA) exhibited the highest levels of genetic diversity. Some native populations show a relatively high number of ROH such as San Pedro (SP1 and SP2), Scotland SC2, Dutch populations PWB, VZM and HGV, France (PP) and Hungary South (HS) (Figure 4b).

The scatterplot depicting each individual F-roh score versus the total number of ROHs (Figure 4c, left) shows around 20 individuals from Spain (SP1 and the Pyrenees), Netherland South VZM and Italy TA with high F-roh values yet low number of ROHs. This pattern is indicative of recent inbreeding resulting in long runs of homozygosity. The scatterplot depicting each individual F-roh score versus the variation (standard deviation) across chromosomes (Figure 4c, right) shows low variation for Norwegian individuals which suggests that their low genetic diversity might be due to ancient population dynamics rather than recent anthropogenic events.

## 3.3 | Genomic structure: Information at population level

The dendrogram depicting Nei genetic distances between populations (Figure 5, see also Figures S9 and S10) shows the maximum distance between Iberia and eastern European populations, and the highest differentiation from the central node for Norwegian and Iberian populations. The tree is structured into four main branches, with native Spanish populations forming a clear one. The admixed population from Eastern Pyrenees, separated into the two genetic groups of individuals, appears to join with the two original clusters

<figure><img src="image_6.png" /><figcaption>**FIGURE 5** Dendrogram depicting Nei genetic distances between populations. The dendrogram was generated with the neighbour-joining (NJ) method applied to an input matrix of Nei genetic distances between populations. Population codes as in Table 1.</figcaption></figure>
of red deer from which they originated, that is PYE₁ with Iberian red deer and PYE₂ with PYN between French PP and Norway. The Scottish populations formed an exclusive branch, rooted between Iberian red deer and the other European populations, which also included the OVP₁ group of individuals restocked from Scotland, while the other group of individuals in the same Dutch population (OVP₂) is in the branch that includes Czechia (CZ) and the remaining Dutch populations. Finally, the easternmost populations from Hungary (HN, HW and HS) and Southern Poland (PS) appeared in the fourth branch at different distances from the populations from Italy (TA), Northern Poland (PN), Austria (AU), France (PP) and the French-originated Pyrenean populations (PYN and PYE₂), with Norway (NO) rooting at the base of this Central-European assembly, but on a long separate branch that denotes its high differentiation in allele frequencies.

## 3.4 | Admixture analyses

The f3-analysis (Figure 6a and Figure S11) identified two Dutch populations as the most probable to exhibit hybrid ancestry, KRO and the two groups of individuals that make up OVP (OVP₁ and OVP₂). The Netherland Centre population KRO showed a complex situation of admixture with the French portion of Eastern

Pyrenean (PYE₂), other Dutch populations, Scotland (KY2, SC1), Czechia (CZ), France (PP), Austria (AU), Italy (TA) and all eastern European (Poland and Hungary) populations. The Netherland North OVP population was also identified as a hybrid population, with OVP₁ sharing ancestry with Scottish SC2 and Hungarian HS, and with OVP₂ connecting with OVP₁ and other Dutch and eastern European ancestral populations.

The maximum likelihood admixture tree with Sika deer as an outgroup shows the inferred relationships between populations based on allele frequency differences (Figure 6b). The OptM procedure on the optimal reduced number of migration edges indicated that trees with 3 and 4 migration edges were those being better supported based on likelihoods and variance explained (Figure S12). Here, we show the tree with four migration edges, because it is more informative about admixture among red deer populations, while the tree with three edges prioritizes possible admixture relationships with the ancestral branch closer to Sika deer and is shown in the Figure S13. The tree with four gene flow edges (Figure 6b) showed a relationship between the branch of Sika deer and eastern European populations. Among red deer populations, the tree highlighted the admixture between both OVP subpopulations, between eastern European populations and the node of CZ and OVP₂ and between the French Pyrenean red deer (PYN and PYE2) and the base of the branch leading to the node of all non-Iberian red deer populations.

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

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

FIGURE 6 Admixture analyses. (a) Matrix highlighting population triplets (X, Y, Z) with near-negative autosomal f3-scores. A negative f3-score, calculated as $(y-x)-(y-z)$, indicates that allele frequencies y in population Y are intermediate between the allele frequencies x and z in populations X and Z, which is an indication of admixture. The colours in the row and column bars represent donor populations X and Z, while field colours denote admixed population Y. The obtained f3-scores (assuming an arbitrary threshold of 0.0015 instead of 0) indicate that populations 'HollandNorthOVP' and 'HollandCentreKRO' are the most likely admixed. Population codes as in Table 1. Split populations are represented by the same colour, so we have indicated the population codes for clarity. (b) Treemix admixture graph. Maximum likelihood tree of genetic distance based on allele frequency differences, generated with the software Treemix for the study populations with sika deer as an outgroup, showing four main flow edges of admixture, after prioritized by OptM for optimal number of migration edges.
## 3.5 | Estimating demographic events

For the most divergent and admixture-free populations (Norway and the two main Iberian ones), we were able to estimate the demographic history for the last c. 250 generations (Figure S14). The resulting effective population size ($N_e$) estimate at 250 generations back was very low (around 100) for Norway and for Spain West (SP1&2), while it was very high (above 10,000) for the Spain Central lineage (SM1&2), thus suggesting that demographic constraints are unlikely to be the cause of current genetic divergence for the Spain Central populations, while they can have a main role for the other two populations. In support of these demographic reconstructions, information in the caption of Figure S14 discusses points of agreement between the $N_e$ fluctuations and historical data from other sources.

## 4 | DISCUSSION

Our study provides information on the current genomic structure of red deer throughout Europe, for 736 individuals from 28 populations, genotyped with 35,701 single nucleotide polymorphism markers (SNPs). To our knowledge, this is the first study conducted on red deer genomics at the level of a major part of its distribution range on the European continent. It can contribute to evaluating previously proposed patterns regarding the genetic structure and differentiation based on neutral genetic markers, which is relevant to our knowledge on the evolutionary history of this ecologically and economically very important species, along with supporting management and conservation decisions.

Genome-wide structure highlights the large genetic difference of Norwegian red deer on the basis of allele frequencies compared to other European populations and also brings strong support for the genetic differentiation of native Iberian red deer populations, as evidenced, for example, by Euclidean distances represented in the dendrogram and by the PCoA. These two areas of European red deer distribution stand out genetically in our study, although it should be mentioned that we did not include the very divergent red deer from the Tyrrhenian islands (Sardinia and Corsica) and the Mesola region in Italy (see, e.g. Doan et al., 2022).

### 4.1 | Norwegian peculiarities

The Norwegian population stands out in almost all measured parameters. Specifically, it is the most divergent according to Euclidean distances and Nei genetic distances; formed an exclusive ancestry cluster at K=4 and higher K values; presented the lowest heterozygosity and highest mean inbreeding and showed the highest number and proportion of the genome in homozygosity runs. Additionally, the low variation among individuals in all parameters, and especially in the distribution of homozygosity runs across chromosomes, points to ancient population dynamics rather than recent anthropogenic events as the most important

contributor to its contemporary situation. Pleistocene fossil remains of red deer from Scandinavia have been dated from 4 to 6.5 ky BP (Doan et al., 2022; Rosvold, 2006), in agreement with the post-glacial recolonization dynamics from the southern refugia in south-western Europe (Hewitt, 2000, 2004; Meiri et al., 2013; Sommer et al., 2008; Taberlet et al., 1998). Our results support that the range expansion northward, to reach the most extreme locations in Europe during several thousands of years in the Holocene, may have produced selection and drift events in Norwegian red deer. Strong bottleneck signals have been found previously in several Norwegian localities (Haanes et al., 2011; Haanes, Røed, Flagstad, et al., 2010) and a microsatellite study of 27 populations across Europe found Norwegian red deer to be among the least diverse as measured by allelic richness and effective population size (Zachos et al., 2016). In addition, the isolation of the Scandinavian peninsula from mainland Europe since the Mid-Holocene (Andersen & Borns, 1994; Jonsson, 1995) may have prevented later admixture of red deer populations, provided that there were no anthropogenic restocking events (but see review in Haanes, Røed, Mysterud, et al. (2010) for possible anthropogenic influences in Swedish populations compared to Norwegian ones). Haanes et al. (2011) found shared haplotypes between British and Norwegian red deer populations and suggested that possible human introductions in recent times, even in the Viking Age, could not be ruled out as a cause of this haplotype sharing. There are reports of very recent (last decades) introductions of non-native individuals in particular populations in Norway (Haanes et al., 2013), but this appears not to have affected our analysis. Our genomic assessment for the sampled population in Norway shows genomic features (namely, the high frequency and number of small ROH widely distributed across chromosomes, as well as no signs of recent admixture) that appear more consistent with earlier events during the post-glacial recolonization than with more recent, likely anthropogenic impacts. These findings largely mirror the genetic signatures of the post-glacial colonization pattern of Scandinavia recently seen in another widespread mammal in Europe, the red fox (McDevitt et al., 2022).

### 4.2 | Iberian red deer evolutionary distinctiveness and internal structure

In the southernmost part of their European distribution, red deer in Iberia show, after Norway, the highest Euclidean genetic distances to the individuals of the remaining populations in Europe within our sample (Figure 2a). Admixture/ancestry analyses at low K (K=2) resulted in genetic clusters with high inclusion values only for the individuals from Iberian populations, even before Norway differentiated and formed an exclusive cluster (K=4). Similar results were obtained with microsatellite markers (Carranza et al., 2016; Queirós et al., 2019, 2020), although present results (e.g. PCoA) show higher differentiation than the analyses performed with microsatellite markers (e.g. factorial correspondence analyses, FCA) by Carranza
et al. (2016) and Queirós et al. (2019), in line with the higher resolution power expected from genomic versus genetic datasets.

Within the Iberian Peninsula, the populations from Sierra de San Pedro (SP1 and SP2) appear genetically different from the other native populations in their original ranges (SM1, SM2 and MO). Individuals of SP1 and SP2 maintain high Euclidean distances from the other Iberian individuals, occupy differentiated positions in the PCoA and show high assignment values to an exclusive cluster for K=5 in the admixture/ancestry analyses. Sierra de San Pedro populations do not show signs of admixture, and also appear well differentiated from the rest of Iberia in the dendrogram based on Nei genetic distances between populations. The differentiation of the Sierra de San Pedro populations was also detected with microsatellite markers by Carranza et al. (2016). Furthermore, red deer from Sierra de San Pedro were shown to belong to a genetically differentiated lineage (Carranza et al., 2016; Carranza & Martínez, 2014; Fernández-García et al., 2014) or sub-lineage (Queirós et al., 2019) based on mtDNA (D-loop) haplotypes. Thus, our results support that, within the Iberian Peninsula, at least two red deer (meta)populations occur with particular genetic composition due probably to divergent ranges during glacial periods (as evidenced by mtDNA divergence; Carranza et al., 2016) and independent evolutionary histories: on the one hand those in the Sierra de San Pedro area in Extremadura and on the other hand, those in the Central and Southern Iberian ranges, from Sierra Morena to Monfragüe and including intermediate areas such as Montes de Toledo (not sampled here, but with the same genetic composition as our SM2 and PYS populations on the basis of microsatellite analyses; Pérez-González et al., 2023). This result is in contrast to those from Queirós et al. (2019) who also found, with microsatellite markers, different genetic clusters in Spain, but higher differentiation between SMR (SM1 in our study) and MTR (SM2 in our study) than between Western and Central Iberia: here SP vs. SM). However, mitochondrial results in Queirós et al. (2019) are in agreement with our findings, since they obtained higher genetic differentiation of PNB (SP2 in our study). The presence of two or more mtDNA lineages in Iberia is a common result, also reported for other mammal species such as the rabbit *Oryctolagus cuniculus* (Branco et al., 2002) and hares *Lepus* spp. (Pérez-Suárez et al., 1994) and was argued to be due to the existence of Pleistocene local, separated refugia within the Iberian glacial refugium during the Pleistocene (Gómez & Lunt, 2007). Interestingly, the only other native deer species, the roe deer (*Capreolus capreolus*) also shows signatures of two genetically distinct (meta)populations, likely as a consequence of different glacial refugia, with their geographical distribution ranges in north-western and central-southern Iberia (Royo et al., 2007).

The Iberian red deer population from Monfragüe (MO) also shows some unique features when results with SNPs are compared with those obtained with mitochondrial sequences. Based on SNP markers, MO had a closer genetic relationship with the populations of Central Iberia (Sierra Morena; SM1 and SM2) than with the Western Iberian populations (SP1 and SP2). In the PCoA plot, MO and SM1/SM2 shared similar scores for the two axes with higher explained variance and tended to have high assignment values for

the same clusters under lower cluster values (K=2 to K=7) in the admixture/ancestry analyses. Moreover, Queirós et al. (2019) found a pattern of genetic admixture between SM1 (SMR in Queirós et al. (2019)) and SM2 (MTR in Queirós et al. (2019)) in red deer from Monfragüe. In previous analyses for the MO population, Carranza et al. (2016) found mtDNA haplotypes from the Western Iberian lineage and Queirós et al. (2019) found haplotypes that are in accordance with a genetically differentiated sub-lineage. Therefore, adding our genomic data to previous results, we can state that the MO population shows mtDNA signs of an evolutionary history of isolation from Central Iberian lineages, encompassing both private and shared haplotypes with the Western Iberian lineage, but with more recent gene flow from Central Iberian populations. This was likely favoured by its geographical location (Monfragüe National Park) located in the Extremadura region (Western Spain) but closely connected to Montes de Toledo (Central Iberian lineage) through natural corridors (see, e.g. Mingarro & Lobo, 2021), which finally shaped their genomic position as a branch between Western (SP1 and SP2) and Central (SM1 and SM2) Iberian lineages (based on Euclidean and Nei genetic distances) but with higher similarities to Central than Western Iberian populations despite its own peculiarities (based on admixture/ancestry analyses). To add further peculiarities to this population, admixture analyses, in particular the Treemix with 3 edges (Figure S13), indicate a possible relationship with ancestors at the basal node of all European red deer network. As an explanation for this, our results do not support a recent admixture, since there is no sign of connection with current European populations besides Iberian neighbours. Therefore, the potential maintenance of ancestral features due to ancient expansion waves from central areas of the continent and later isolation may be plausible and deserving of further research. Despite the Monfragüe (MO) population being analysed in previous studies using nuclear and mitochondrial markers (Carranza et al., 2016; Queirós et al., 2019), the peculiarities found here on either the discordance with respect to mtDNA results or on the suggestion of an ancestral connection with the root of European red deer had not previously been detected, evidencing that the higher resolution of genome-wide markers can reveal subtle, hitherto unnoticed patterns (Biello et al., 2022; Zimmerman et al., 2020).

At the higher *K* values in the admixture/ancestry analyses (e.g. *K*=10), four genetic clusters for Iberian populations appear. Besides the SP1+SP2 and MO populations, which were assigned to two single independent clusters, red deer from Sierra Morena at *K*=10 exhibited assignment values to two clusters that highlight the differences between Western (SM1) and Eastern (SM2) Sierra Morena at both sides of the Valle de los Pedroches area in Cordoba province. This result supports the previously found genetic substructure in southern Spain (Carranza & Martínez, 2014; Fernández-García et al., 2014; Galarza et al., 2015). From these previous papers and our own results, within the Central-South Iberian lineage (Carranza et al., 2016), SM2 shows more genetic similarity to the Central rather than to the Southern Iberian populations, and hence, in our ancestry analysis it shares adscription at *K*=10 with MO and PYS, since MO is connected to Central areas (Montes de Toledo) as indicated above, and PYS was
restocked through translocations of deer from Montes de Toledo to areas in the south of the Pyrenees as mentioned in Methods.

For the Pyrenean populations sampled in Spain, our SNP ge-
nomic assessment agrees with published results based on micro-
satellites (Pérez-González et al., 2023; Queirós et al., 2019) but
also provides further insights. PYS is mostly composed of Iberian
red deer, while PYN clusters with other non-Iberian populations,
and PYE represents a recent contact zone containing individuals
from both PYS and PYN as well as their hybrids. The French side of
the Pyrenees experienced several restocking events with individu-
als that were reported to come mainly from the Chambord and La
Petite Pierre populations (Klein, 1990; Saint-Andrieux et al., 2004).
Thus, we would expect some similarities between our PYN (also
PYE₂) and PP populations. However, our genomic analysis shows
no connection between the native French population (PP) and any
populations sampled in the Pyrenees, PP being closely related to
central European ones. Likewise, our results for the native French
population (PP) and the French population in the north Pyrenees
(PYN) do not show any signs of relationship with populations in
Spain. During the Last Glacial Maximum (LGM), red deer per-
sisted in Spain and adjacent areas in southern France (Queirós
et al., 2019; Sommer et al., 2008), so some closer genetic distance
between Spanish and French populations might have been ex-
pected, even considering that our French samples did not belong
to the original native populations from the northern Pyrenees but
were restocked from other French populations (Klein, 1990; Saint-
Andrieux et al., 2004).

West European red deer populations (currently belonging to mtDNA haplogroup A) were recolonized from southern-western glacial refugia where they survived during the LGM (Doan et al., 2022; Niedziałkowska, Doan, et al., 2021), but the origin and dynamics of these recolonizations are still unclear (Carranza et al., 2016; Doan et al., 2022; Meiri et al., 2013; Queirós et al., 2019; Rey-Iglesia et al., 2017; Skog et al., 2009). The genetic distance found here between French and Spanish populations, along with the clustering of French populations with other continental European ones, suggests a history of isolation between Iberian and other European populations in the post-LGM period. At the same time, Iberian red deer encompass mtDNA lineages with estimated divergence times older than the divergence between some Iberian and central European red deer. Estimating divergence times is difficult and highly dependent on assumptions. However, just considering the relative figures, the mtDNA divergence between the most distant Iberian lineages (Sw vs. Sc, in Carranza et al., 2016; SP1&SP2 vs. SM1&SM2 in our data) clearly pre-date the divergence time between some (mtDNA) Iberian lineages and central European ones. This means that the last connection between Iberia and central European populations took place more recently than the split of mitochondrial lineages within Iberia, and that ancestral variation in the mitochondrial genome has been retained in Iberia. But also, that non-Iberian European populations, from France to the most distant ones of haplogroup A, diverged more recently and in isolation from Iberia.

Understanding the role of the Iberian Peninsula in the red deer colonization of Europe is crucial to reconcile the possible expansion-contraction models with successive waves of refugial and recolonization ranges (Doan et al., 2022; Hewitt, 1996, 1999; Meiri et al., 2013; Rey-Iglesia et al., 2017; Sommer et al., 2008; Taberlet et al., 1998). Our SNP analysis shows that in the contemporary situation of genetic structure, the Iberian Peninsula appears well differentiated within the European continent, clearly distinct from French and other European populations. On the other hand, it shows (i) signs of ancient differentiation in some populations (SP1 and SP2) and connection with early branches of European red deer ancestors for MO in admixture trees, which evokes ancient waves of Pleistocene migration from eastern Europe or Asian-Palearctic areas. Secondly (ii), populations in Central Spain (notably SM1 and SM2) show differentiation from central European ones, compatible with an earlier (c. pre-LGM or early Holocene) connection between them and all West European lineages (haplogroup A), followed by (iii) very low gene flow in more recent periods, likely due to the Pyrenean mountains acting as a barrier along with the particular ecological conditions of Iberia during the late Holocene. Several lines of evidence indicate that Iberian lineages participated unequally in this LGM connection with the other Western European red deer. Carranza et al. (2016), on the basis of mtDNA analyses, suggested that it was Central but not Western Iberian lineages that supported this connection. A similar case may be that of the extinct mtDNA haplogroup W2, found in northern Spain and dated to just after the LGM that might have come to Iberia in some of the previous migration waves from Asia but did not contribute to the post-glacial recolonization of Europe (Doan et al., 2018, 2022; Meiri et al., 2013).

Our results from Iberia and Norway depict opposite histories.
Iberian refugial areas may have received waves of migrants, more
easily overlapping in geographical range if closer to the source of
migration (e.g. central Iberian populations), while more extreme
refugial areas (i.e. western Iberia) retain ancestral genotypes from
previous migration waves. Our results also suggest that recolonizing
expansions to northern areas such as Norway may be less common,
and may experience more bottlenecks and stronger drift events,
compared to less remote geographical ranges.

The isolated evolution of Iberian red deer may be further sup-
ported by peculiar phenotypic features described previously (e.g.
lack of neck mane, different vocal repertoire, tongue protrusion
while roaring; Frey et al., 2012; Geist, 1998; Passilongo et al., 2013;
Volodin et al., 2013), although we still need further studies to un-
derstand the reasons for these potential adaptations and whether
there may be a potential gradient of phenotypic variation through-
out European populations that could match genetic variation.

4.3 | European mainland and UK red deer

The distinctiveness of Norway and Iberia stands out more if we look
only at the native populations (see Figure S7), with the remaining
populations in Europe being much more similar to each other and
lying along a continuous gradient of variation from Scotland to Hungary in the PCoA (Figure S7).

Among the red deer sampled in Britain, the Rum population is the most differentiated, consistent with its mixed restocking origins since 1845 with individuals from different UK populations that may, in turn, have contained introgressed genes from other Scottish, mainland European and Sardinian/North-African populations (Nussey et al., 2006). Despite their admixed origin, genome-wide assessment reveals similarity within Rum samples, indicating homogenization after a high enough number of generations without further introductions. The Rum population exhibits high genetic diversity (multi-locus heterozygosity, relatively few and short ROHs). It appears closer to Scotland SC2, pointing to a source of origin, but its heterozygosity is higher than in SC2, evidencing admixture also from elsewhere. The population of Rum is perhaps the most studied red deer population in the world and has greatly contributed to our knowledge of this species (Pemberton et al., 2022). As it is an artificially founded population, genomic results may stand as a necessary basis for interpreting other related findings and drawing generalized conclusions.

Red deer populations in the Netherlands encompass different histories of anthropogenic effects. The OVP population was restocked from Scotland, KRO exhibits admixture with several European populations as revealed by the admixture analysis at the population level (Figure 6a) and individuals appear especially associated with eastern European populations (Figure 2c). The other Dutch populations (PWB, VZM, HGV) are genetically diverse but with signs of inbreeding. Thus, all admixture signals outlined here appear to be the result of recent restocking events and seem to completely hide the signals from post-glacial recolonization (see de Jong et al., 2020 for further details on Dutch populations).

## 4.4 | Peripheral versus mainland areas and hypothetical expansion-contraction routes

The present genomic analysis reveals that Iberian red deer exhibit lower genetic distances (see, e.g. Figure 5) to British red deer than to other European populations, despite its shorter geographical distance to mainland western Europe, and with no signs of recent admixture between both. This observation suggests a shared genetic legacy between Iberian and UK populations from a common migration wave, with the latter potentially preserving this genetic signature more prominently than their mainland counterparts during the recolonization events in the Holocene, which may also lend support to putative UK cryptic glacial refugia (Queirós et al., 2019). This connection between the Iberian Peninsula and the UK likely dates back to the LGM-Younger Dryas period, and affected only central Spain populations as denoted by the mtDNA similarities between them (Carranza et al., 2016), while western Spain populations appear to have preserved genetic (mtDNA) signatures indicative of an earlier migration wave. Figure 7 shows a hypothetical scenario that may be compatible with our findings. Several sources of evidence, both from our data and from previous mtDNA studies

as discussed above, suggest that West Iberian populations, as well as some areas in the extreme west of Europe (see Doan et al., 2022; Niedziałkowska, Doan, et al., 2021), may have retained ancestral signatures of ancient waves in the Pleistocene (labelled 1 in Figure 7). Contractions to refugial areas during the Last Glacial Period may have introduced populations from the northern mainland to the peripheral refugial areas of Iberia and areas situated south of Britain (lines 2 in Figure 7). Holocene recolonizations may have occurred from southern areas of Europe and the British Isles to the northern areas including Scandinavia, but probably not from Iberia (lines 3 in Figure 7), as denoted by the high genetic distance maintained between Iberia and the rest of the continent. In the Late Holocene, peripheral areas such as the British Isles, Scandinavia and Iberia remained isolated from the continental mainland due to sea channels or the Pyrenees, respectively, which may have caused the maintenance of common features between Scottish and Iberian red deer, as well as the gene flow between populations on the European mainland (line 4 in Figure 7), while Norwegian populations diverged from their origin from mainland Europe due to demographic events related to small population size.

## 4.5 | Main European lineages and putative subspecies

Our sampling areas lay within the current distribution ranges of the two main red deer lineages in Europe that belong to the refugial populations in Iberia and southern France (mtDNA haplogroup A) and the Balkans (mtDNA haplogroup C) (Doan et al., 2022; Ludt et al., 2004; Meiri et al., 2013; Niedziałkowska et al., 2011; Skog et al., 2009). Suture zones of contact between both recolonizing lineages (Taberlet et al., 1998) have been described in Czechia (Krojerová-Prokešová et al., 2015), Poland (Niedziałkowska et al., 2011) and Hungary, where also haplogroup B, typical of the Tyrrhenian islands (Corsica, Sardinia) and North Africa, was found in the south (Markov et al., 2015). Despite the overlapping range of mtDNA haplogroups, studies using microsatellite markers have evidenced little admixture between individuals belonging to these two lineages (e.g. Krojerová-Prokešová et al., 2015; Zachos et al., 2016). The occurrence of admixture at contact zones is a natural phenomenon in many species (Taylor & Larson, 2019) but in the case of red deer, because of the economic interest of the species mainly for hunting purposes, human-mediated translocations may have been common with individuals moved from one lineage into the area of another (Apollonio et al., 2014). In eastern Europe, red deer translocations may have contributed to the overlapping of the original distribution ranges, mainly because of the introduction of individuals of the western lineage into areas in eastern Europe where the species became extinct in historical times (Niedziałkowska et al., 2011, 2012; Szederjei & Szederjei, 1971). While the overall signal of western and eastern mitochondrial lineages (A and C respectively) has been confirmed for microsatellites in a continent-wide study (Zachos et al., 2016), our results show
<figure><img src="image_9.png" /><figcaption>FIGURE 7 Hypothetical migration waves and admixture zones. The arrows are not intended to be geographically precise but only to be schematic visualization of the ideas discussed in the text. Areas with recent, mostly human-mediated admixture are indicated with asterisks. Base maps from Generic Mapping Tools: https://wwwwatson mapping tools org/.</figcaption></figure>

that, at a genome-wide European level with higher resolution, the contact zone in eastern Europe has greatly blurred the genomic distinction that should have existed between red deer belonging to the Iberian and Balkan glacial lineages (see also Frank et al., 2017, 2022), which appear very similar to each other when compared for instance with the differences that Norwegian or Iberian red deer maintain from other central European populations that belong to the same glacial lineage.

The observed differences in the genetic structure between mtDNA (Doan et al., 2022; Skog et al., 2009) and autosomal SNPs (this study) in central Europe might be attributed to stronger admixture processes at autosomal loci, likely influenced by the prevailing male-biased dispersal pattern commonly observed in red deer populations (Carranza et al., 2016; Clutton-Brock et al., 1982). Besides anthropogenic admixture, natural admixture processes are expected to be more likely in the central areas of the distribution range, while genetic distinctiveness may be maintained at the edges of the distribution range, either because of the long distance and differentiation from refugial areas (e.g. Norwegian populations) or because not all refugial populations contributed to post-glacial re-colonization so that some may maintain isolation in refugial areas

(e.g. Iberian populations). Similar patterns of genetic distinctness of southern and northern European populations located at the edges of the species ranges have been detected in moose (*Alces alces*; Niedziałkowska et al., 2016), wild boar (*Sus scrofa*; e.g. Khederzadeh et al., 2019; Niedziałkowska, Tarnowska, et al., 2021) and roe deer (*Capreolus capreolus*; Plis et al., 2022).

Our sampling areas include the putative ranges of several described subspecies: *Cervus elaphus hispanicus* (Iberia), *C.e. hippelaphus* (western mainland Europe), *C.e. scoticus* (British Isles), *C.e. atlanticus* (Norway) and *C.e. pannoniensis* (eastern mainland Europe) (Carranza et al., 2024; Geist, 1998; Mattioli et al., 2022). Subspecies delimitation will always be difficult and involve grey-area decisions, yet our genome-wide SNP analysis brings some support for the current differentiation between some of them but not for others. In particular, *C.e.atlanticus*, distributed in Norway, maintains the greatest differentiation in our dataset; likewise, Iberian red deer *C.e.hispanicus* is supported by a concordant combination of the isolated Iberian distribution range and the genomic distinction. A similar situation may be that of Scottish red deer *C.e.scoticus*, although admixture is evident in many areas. However, with respect to *hippelaphus* versus *pannoniensis* (Carranza et al., 2024), besides their clearly distinct mtDNA
phylogenetic origin (Ludt et al., 2004; Niedziałkowska et al., 2011; Skog et al., 2009), our genomic results provide much less support for a current separation of these two into distinct subspecies, and their contemporary geographical differences are probably rather a clinal gradient of variation along the continent akin to what has been shown in other widespread mammals using genomic data (e.g. McDevitt et al., 2022). This tentative taxonomic conclusion differs from recent results based on modern and ancient mtDNA by Meiri et al. (2018) who keep *pannoniensis* separate but argue for lumping *atlanticus*, *hispanicus*, *scoticus* and *hippelaphus*. Our study, thus, also showcases the taxonomic implications of the high resolution of genomic data, even in a species as detailedly studied as red deer.

## 4.6 | Conservation and management implications

All these results have management and conservation implications. On the one hand, the preservation of what can be considered sub-specific (i.e. fundamental intraspecific) distinctions should be a priority (e.g. Carranza et al., 2003; Zachos & Hartl, 2011). The Iberian red deer populations sampled show a clear differentiation in relation to the other European populations, so they represent a good example of an evolutionarily significant unit that deserves protection. It is likely that some Spanish populations may have suffered introgression based on historical reports (Carranza et al., 2003) and microsatellite analyses (Fernández-García et al., 2014; Queirós et al., 2020), but our genome-wide assessment indicated that, at least for our sampled populations, Iberian red deer are still well differentiated and little admixed, thus deserving conservation efforts. Even within Iberian red deer, our results show that we might identify significant units in the western populations (for which it may be warranted to distinguish MO with respect to SP1 and SP2) and in the central populations (SM1 and SM2), in accordance with previous mtDNA analyses (Carranza et al., 2016; Queirós et al., 2019). One main threat to the integrity of Iberian red deer is the current incursion of individuals from restocked French populations southwards through the Pyrenees, for which a recent study recommends implementing measures to prevent contacts (Pérez-González et al., 2023). Norwegian and Scottish populations also deserve specific conservation programmes, as they maintain distinctiveness mostly produced by a history of isolation. In eastern Europe, the ancestral pre-glacial differentiation of A and C mtDNA lineages (Doan et al., 2022; Meiri et al., 2013; Niedziałkowska et al., 2011; Skog et al., 2009) seems mostly clinal at the genomic level. However, the use of genome-wide approaches to identifying natural populations free from anthropogenic impact on their genetic composition may be a desirable conservation effort. Other areas of Europe, such as the Netherlands (de Jong et al., 2020), seem to be too affected by human translocations. But even so, genomic procedures may go a long way towards identifying conservation targets. It should be stressed that genomic studies can provide precise information at the level of the individual that can be difficult to obtain with previous genetic methods. Thus, for most areas of the red deer range where considerable anthropogenic

impact has already occurred, genomic studies are able to unveil genetic signatures of recent and past human-mediated effects. They can establish baseline information on effective population size, inbreeding, landscape connectivity and barriers, as well as the presence of allochthonous ancestry, and become an essential tool for decision-making in management and conservation.

## AUTHOR CONTRIBUTIONS

J.C and J.P.-G. conceived the study and wrote/edited the manuscript; G.A., M.d.J., J.P.G and P.M. performed bioinformatics and statistical analyses; C.B., E.P, G.V. and A.M. carried out laboratory work; F.Z., A.M., M.N., M.S., S.C., N.B., L.C., K.R., C.S.-A., A.B., A.G.E., J.R.O., J.M.S., S.M.M. and J.P. provided tissue samples or genotypes; J.C. obtained funding; J.C. and A.M. managed projects. All authors provided feedback on the manuscript and have approved the final manuscript for submission.

## AFFILIATIONS

## ACKNOWLEDGEMENTS

We thank people from Tarvisio (Mario De Bortoli, director of the Tarvisio Malborghetto hunting reserve) and Austria (Luigi Marcati, Seefeld in Tyrol hunting reserve) who helped with the sample collection. Armando Caballero helped in demographic estimation with GONE. Rafael Serrano helped with historical data to interpret N<sub>e</sub> variations. João Queirós and one anonymous reviewer made valuable comments that greatly contributed to improve the manuscript.

<footnote n="1">¹Wildlife Research Unit (UIRCP), University of Córdoba, Córdoba, Spain</footnote>

<footnote n="2">²Biology and Ethology Unit, Veterinary Faculty, University of Extremadura, Cáceres, Spain</footnote>

<footnote n="3">³Senckenberg Biodiversity and Climate Research Institute (SBiK-F), Frankfurt am Main, Germany</footnote>

<footnote n="4">⁴Natural History Museum Vienna, Vienna, Austria</footnote>

<footnote n="5">⁵Department of Genetics, University of the Free State, Bloemfontein, South Africa</footnote>

<footnote n="6">⁶Department of Evolutionary Biology, University of Vienna, Vienna, Austria</footnote>

<footnote n="7">⁷Research Institute for the Environment and Livelihoods, Charles Darwin University, Casuarina, Northwest Territories, Australia</footnote>

<footnote n="8">⁸Department of Natural Resources and the Environment, Atlantic Technological University, Galway, Ireland</footnote>

<footnote n="9">⁹Mammal Research Institute Polish Academy of Sciences, Białowieża, Poland</footnote>

<footnote n="10">¹⁰Department of Wildlife Biology and Management, Hungarian University of Agriculture and Life Sciences (MATE), Gödöllő, Hungary</footnote>

<footnote n="11">¹¹Department of Game Management, Ministry of Agriculture, Budapest, Hungary</footnote>

<footnote n="12">¹²Department of Preclinical Sciences and Pathology, Norwegian University of Life Sciences, Ås, Norway</footnote>

<footnote n="13">¹³French Biodiversity Agency, Vincennes, France</footnote>

<footnote n="14">¹⁴Department of Animal and Food Science, Universitat Autònoma de Barcelona, Barcelona, Barcelona, Spain</footnote>

<footnote n="15">¹⁵Dirección General d'Ecosistemes Forestals i Gestió del Medi, Barcelona, Spain</footnote>

<footnote n="16">¹⁶Department of Integrative Ecology, Estación Biológica de Doñana (CSIC), Seville, Spain</footnote>

<footnote n="17">¹⁷Faculty of Biotechnology, University of Wrocław, Wrocław, Poland</footnote>

<footnote n="18">¹⁸Instituto de Investigación en Recursos Cinegéticos, IREC (CSIC, UCLM, JCCM), Ciudad Real, Spain</footnote>

<footnote n="19">¹⁹Department of Biology, York University, Toronto, Ontario, Canada</footnote>

<footnote n="20">²⁰Institute of Ecology and Evolution, School of Biological Sciences, University of Edinburgh, Edinburgh, UK</footnote>
Projects UCO-1262830 and P18-RT-4642 (Government of Andalusia and ERDF of European Union) to JC provided financial support.

## CONFLICT OF INTEREST STATEMENT

The authors declare no conflicts of interest.

## DATA AVAILABILITY STATEMENT

Individual genotype data are available on Dryad repository: doi:10.5061/dryad.1rn8pk13b.

## BENEFIT-SHARING STATEMENT

Benefits Generated: A research collaboration was developed with scientists from the countries providing genetic samples, all collaborators are included as co-authors, the results of research have been shared with the provider communities and the broader scientific community (see above), and the research addresses a priority concern, in this case, the conservation of organisms being studied. More broadly, our group is committed to international scientific partnerships, as well as institutional capacity building. The authors declare that the research described in the publication complies with relevant national laws.

## ORCID

Juan Carranza https://orcid.org/0000-0002-0368-7173

Magdalena Niedziałkowska https://orcid.org/0000-0002-7746-4116

Jose Antonio Godoy https://orcid.org/0000-0001-7502-9471

S. Eryn McFarlane https://orcid.org/0000-0002-0706-458X

## REFERENCES

Alexander, D. H., Novembre, J., & Lange, K. (2009). Fast model-based estimation of ancestry in unrelated individuals. *Genome Research*, 19(9), 1655-1664. https://doi.org/10.1101/gr.094052.109

Allendorf, F. W., Hohenlohe, P. A., & Luikart, G. (2010). Genomics and the future of conservation genetics. *Nature Reviews Genetics*, 11(10), 697-709. https://doi.org/10.1038/nrg2844

Allendorf, F. W., Leary, R. F., Spruell, P., & Wenburg, J. K. (2001). The problems with hybrids: Setting conservation guidelines. *Trends in Ecology & Evolution*, 16, 613-622.

Andersen, B. G., & Borns, H. W. (1994). *The ice age world*. Scandinavian University Press. https://www.nature.com/articles/nrg2844.

Apollonio, M., Scandura, M., & Šprem, N. (2014). Reintroductions as a management tool for European ungulates. In R. Putman & M. Apollonio (Eds.), *Behaviour and management of European ungulates* (pp. 46-77). Whittles Publishing.

Biello, R., Zampiglia, M., Fuselli, S., Fabbri, G., Bisconti, R., Chiocchio, A., Mazzotti, S., Trucchi, E., Canestrelli, D., & Bertorelle, G. (2022). From STRs to SNPs via ddRAD-seq: Geographic assignment of confiscated tortoises at reduced costs. *Evolutionary Applications*, 15(9), 1344-1359. https://doi.org/10.1111/eva.13431

Branco, M., Monnerot, M., Ferrand, N., & Templeton, A. R. (2002). Postglacial dispersal of the European rabbiy (*Oryctolagus cuniculus*) on the Iberian Peninsula reconstructed from nested clade and mismatch analyses of mitochondrial DNA genetic variation. *Evolution*, 56(4), 792-803. https://doi.org/10.1111/j.0014-3820.2002.tb01390.x

Brauning, R., Fisher, P. J., McCullock, A. F., Smithies, R. J., Ward, J. F., Bixley, M. J., Lawley, C. T., Rowe, S. J., & McEwan, J. C. (2015).

Utilization of high throughput genome sequencing technology for large scale single nucleotide polymorphism discovery in red deer and Canadian elk. *BioRxiv*, 027318.

Brumfield, R. T., Beerli, P., Nickerson, D. A., & Edwards, S. V. (2003). The utility of single nucleotide polymorphisms in inferences of population history. *Trends in Ecology & Evolution*, 18(5), 249-256. https://doi.org/10.1016/S0169-5347(03)00018-1

Carranza, J., & Martínez, J. G. (2014). Erratum to: Mitochondrial D-loop phylogeny signals two native Iberian red deer (*Cervus elaphus*) lineages genetically different to Western and Eastern European red deer and infers human-mediated translocations. *Biodiversity and Conservation*, 23(3), 555-557. https://doi.org/10.1007/s10531-014-0623-8

Carranza, J., Martínez, J. G., Sánchez-Prieto, C. B., Fernández-García, J. L., Sánchez-Fernández, B., Álvarez-Álvarez, R., Valencia, J., & Alarcos, S. (2003). Game species: Extinction hidden by census numbers. *Animal Biodiversity and Conservation*, 26(2), 81-84.

Carranza, J., Mattioli, S., Putman, R., Pérez-Barbería, J. F., & Gordon, I. G. (2024). Red Deer *Cervus elaphus* (Linnaeus, 1758). In M. Melleti & S. Focardi (Eds.), *Deer of the world: Ecology, conservation and management. Chapter 10*. Springer. Nature.

Carranza, J., Salinas, M., Andrés, D., & Pérez-González, J. (2016). Iberian red deer: Paraphyletic nature at mtDNA but nuclear markers support its genetic identity. *Ecology and Evolution*, 6(4), 905-922. https://doi.org/10.1002/ece3.1836

Clutton-Brock, T., Guinness, F., & Albon, S. (1982). Red Deer: Behaviour and ecology of two sexes. University of Chicago Press.

de Jong, J. F., van Hooft, P., Megens, H.-J., Crooijmans, R. P. M. A., de Groot, G. A., Pemberton, J. M., Huisman, J., Bartoš, L., Iacolina, L., van Wieren, S. E., Ydenberg, R. C., & Prins, H. H. T. (2020). Fragmentation and translocation distort the genetic landscape of ungulates: Red deer in The Netherlands. *Frontiers in Ecology and Evolution*, 8, 535715. https://doi.org/10.3389/fevo.2020.535715

de Jong, M. J., Jong, J. F., Hoelzel, A. R., & Janke, A. (2021). SambaR: An R package for fast, easy and reproducible population-genetic analyses of biallelic SNP data sets. *Molecular Ecology Resources*, 21(4), 1369-1379. https://doi.org/10.1111/1755-0998.13339

de Jong, M. J., Niamir, A., Wolf, M., Kitchener, A. C., Lecomte, N., Seryodkin, I. V., Fain, S. R., Hagen, S. B., Saarma, U., & Janke, A. (2023). Range-wide whole-genome resequencing of the brown bear reveals drivers of intraspecies divergence. *Communications Biology*, 6(1), 153. https://doi.org/10.1038/s42003-023-04514-w

De Kort, H., Prunier, J. G., Ducatez, S., Honnay, O., Baguette, M., Stevens, V. M., & Blanchet, S. (2021). Life history, climate and biogeography interactively affect worldwide genetic diversity of plant and animal populations. *Nature Communications*, 12(1), 516. https://doi.org/10.1038/s41467-021-20958-2

Doan, K., Mackiewicz, P., Sandoval-Castellanos, E., Stefaniak, K., Ridush, B., Dalén, L., Węgleński, P., & Stankovic, A. (2018). The history of Crimean red deer population and *Cervus* phylogeography in Eurasia. *Zoological Journal of the Linnean Society*, 183(1), 208-225. https://doi.org/10.1093/zoolinnean/zlx065

Doan, K., Niedziałkowska, M., Stefaniak, K., Sykut, M., Jędrzejewska, B., Ratajczak-Skrzatek, U., Piotrowska, N., Ridush, B., Zachos, F. E., Popović, D., Baca, M., Mackiewicz, P., Kosintsev, P., Makowiecki, D., Charniauski, M., Boeskorov, G., Bondarev, A. A., Danila, G., Kusak, J., ... Stanković, A. (2022). Phylogenetics and phylogeography of red deer mtDNA lineages during the last 50 000 years in Eurasia. *Zoological Journal of the Linnean Society*, 194(2), 431-456. https://doi.org/10.1093/zoolinnean/zlab025

Doan, K., Zachos, F. E., Wilkens, B., Vigne, J.-D., Piotrowska, N., Stanković, A., Jędrzejewska, B., Stefaniak, K., & Niedziałkowska, M. (2017). Phylogeography of the Tyrrhenian red deer (*Cervus elaphus corsicanus*) resolved using ancient DNA of radiocarbon-dated sub-fossils. *Scientific Reports*, 7(1), 2331. https://doi.org/10.1038/s41598-017-02359-y
Elblinger, E., Bokor, J., Bokor, Á., Altbäcker, V., Nagy, J., Szabó, J., Sárdi, B., Bálteanu, A. V., Rónai, Z., Rózsa, L., Rátky, J., Anton, I., & Zsolnai, A. (2022). Parentage testing and looking for single nucleotide markers associated with antler quality in deer (*Cervus elaphus*). *Archives Animal Breeding*, 65(3), 267-274. https://doi.org/10.5194/aab-65-267-2022

Fernández-García, J. L., Carranza, J., Martínez, J. G., & Randi, E. (2014). Mitochondrial D-loop phylogeny signals two native Iberian red deer (*Cervus elaphus*) lineages genetically different to Western and Eastern European red deer and infers human-mediated translocations. *Biodiversity and Conservation*, 23(3), 537-554. https://doi.org/10.1007/s10531-013-0585-2

Fitak, R. R. (2021). OptM: Estimating the optimal number of migration edges on population trees using Treemix. *Biology Methods and Protocols*, 6(1), bpab017. https://doi.org/10.1093/biomethods/bpab017

Frank, K., Bleier, N., Tóth, B., Sugár, L., Horn, P., Barta, E., Orosz, L., & Stéger, V. (2017). The presence of Balkan and Iberian red deer (*Cervus elaphus*) mitochondrial DNA lineages in the Carpathian Basin. *Mammalian Biology*, 86, 48-55. https://doi.org/10.1016/j.mambio.2017.04.005

Frank, K., Szepesi, K., Bleier, N., Sugár, L., Kusza, S., Barta, E., Horn, P., Orosz, L., & Stéger, V. (2022). Genetic traces of dispersal and admixture in red deer (*Cervus elaphus*) populations from the Carpathian Basin. *European Journal of Wildlife Research*, 68(5), 55. https://doi.org/10.1007/s10344-022-01602-w

Frey, R., Volodin, I., Volodina, E., Carranza, J., & Torres-Porras, J. (2012). Vocal anatomy, tongue protrusion behaviour and the acoustics of rutting roars in free-ranging Iberian red deer stags (*Cervus elaphus hispanicus*). *Journal of Anatomy*, 220(3), 271-292. https://doi.org/10.1111/j.1469-7580.2011.01467.x

Frichot, E., & François, O. (2015). LEA: An R package for landscape and ecological association studies. *Methods in Ecology and Evolution*, 6(8), 925-929. https://doi.org/10.1111/2041-210X.12382

Galarza, J. A., Sanchez-Fernandez, B., Fandos, P., & Soriguer, R. (2015). The genetic landscape of the Iberian red deer (*Cervus elaphus hispanicus*) after 30 years of big-game hunting in southern Spain. *The Journal of Wildlife Management*, 79(3), 500-504. https://doi.org/10.1002/jwmg.854

Gauzere, J., Pemberton, J. M., Slate, J., Morris, A., Morris, S., Walling, C. A., & Johnston, S. E. (2023). A polygenic basis for birth weight in a wild population of red deer (*Cervus elaphus*). G3 Genes | Genomes | Genetics, 13(4), jkad018. https://doi.org/10.1093/g3journal/jkad018

Geist, V. (1998). *Deer of the world: Their evolution, behaviour, and ecology*. Stackpole Books.

Gómez, A., & Lunt, D. H. (2007). Refugia within refugia: Patterns of phylogeographic concordance in the Iberian Peninsula. In S. Weiss & N. Ferrand (Eds.), *Phylogeography of southern European refugia* (pp. 155-188). Springer Netherlands. https://doi.org/10.1007/1-4020-4904-8_5

Haanes, H., Røed, K. H., Flagstad, Ø., & Rosef, O. (2010). Genetic structure in an expanding cervid population after population reduction. *Conservation Genetics*, 11, 11-20. https://doi.org/10.1007/s10592-008-9781-0

Haanes, H., Røed, K. H., Mysterud, A., Langvatn, R., & Rosef, O. (2010). Consequences for genetic diversity and population performance of introducing continental red deer into the northern distribution range. *Conservation Genetics*, 11(5), 1653-1665. https://doi.org/10.1007/s10592-010-0048-1

Haanes, H., Røed, K. H., Perez-Espona, S., & Rosef, O. (2011). Low genetic variation support bottlenecks in Scandinavian red deer. *European Journal of Wildlife Research*, 57(6), 1137-1150. https://doi.org/10.1007/s10344-011-0527-9

Haanes, H., Rosvold, J., & Røed, K. H. (2013). Non-indigenous introgression into the Norwegian red deer population. *Conservation Genetics*, 14, 237-242. https://doi.org/10.1007/s10592-012-0431-1

Hewett, A. M., Stoffel, M. A., Peters, L., Johnston, S. E., & Pemberton, J. M. (2023). Selection, recombination and population history effects on runs of homozygosity (ROH) in wild red deer (*Cervus elaphus*). *Heredity*, 130(4), 242-250. https://doi.org/10.1038/s41437-023-00602-z

Hewitt, G. M. (1996). Some genetic consequences of ice ages, and their role in divergence and speciation. *Biological Journal of the Linnean Society*, 58(3), 247-276. https://doi.org/10.1006/bijl.1996.0035

Hewitt, G. M. (1999). Post-glacial re-colonization of European biota. *Biological Journal of the Linnean Society*, 68(1-2), 87-112. https://doi.org/10.1111/j.1095-8312.1999.tb01160.x

Hewitt, G. M. (2000). The genetic legacy of the quaternary ice ages. *Nature*, 405(6789), 907-913. https://doi.org/10.1038/35016000

Hewitt, G. M. (2004). Genetic consequences of climatic oscillations in the quaternary. *Philosophical Transactions of the Royal Society of London. Series B: Biological Sciences*, 359(1442), 183-195. https://doi.org/10.1098/rstb.2003.1388

Hoffmann, G. S., Johannesen, J., & Griebeler, E. M. (2016). Population dynamics of a natural red deer population over 200 years detected via substantial changes of genetic variation. *Ecology and Evolution*, 6(10), 3146-3153. https://doi.org/10.1002/ece3.2063

Hohenlohe, P. A., Funk, W. C., & Rajora, O. P. (2021). Population genomics for wildlife conservation and management. *Molecular Ecology*, 30(1), 62-82. https://doi.org/10.1111/mec.15720

Johnston, S. E., Huisman, J., Ellis, P. A., & Pemberton, J. M. (2017). A high-density linkage map reveals sexual dimorphism in recombination landscapes in red deer (*Cervus elaphus*). G3 Genes | Genomes | Genetics, 7(8), 2859-2870. https://doi.org/10.1534/g3.117.044198

Jonsson, L. (1995). Vertebrate fauna during the Mesolithic on the Swedish west coast. In A. Fischer (Ed.), *Man and sea in Mesolithic* (pp. 147-160). Oxbow books.

Kardos, M., Luikart, G., & Allendorf, F. W. (2015). Measuring individual inbreeding in the age of genomics: Marker-based measures are better than pedigrees. *Heredity*, 115(1), 63-72. https://doi.org/10.1038/hdy.2015.17

Kasarda, R., Moravčíková, N., Trakovická, A., Krupová, Z., & Ondrej, K. (2017). Genomic variation across cervid species in respect to the estimation of red deer diversity. *Acta Veterinaria*, 67(1), 43-56. https://doi.org/10.1515/acve-2017-0005

Khederzadeh, S., Kusza, S., Huang, C., Markov, N., Scandura, M., Babaev, E., Šprem, N., Seryodkin, I. V., Paule, L., Esmailizadeh, A., Xie, H., & Zhang, Y. (2019). Maternal genomic variability of the wild boar (*Sus scrofa*) reveals the uniqueness of East-Caucasian and Central Italian populations. *Ecology and Evolution*, 9(17), 9467-9478. https://doi.org/10.1002/ece3.5415

Klein, F. (1990). La réintroduction du cerf *Cervus elaphus*. *Revue d'Écologie. Supplément*, 45(5), 131-134. https://doi.org/10.3406/revec.1990.6340

Krojerová-Prokešová, J., Barančeková, M., & Koubek, P. (2015). Admixture of Eastern and Western European red deer lineages as a result of postglacial recolonization of The Czech Republic (Central Europe). *Journal of Heredity*, 106(4), 375-385. https://doi.org/10.1093/jhered/esv018

Ludt, C. J., Schroeder, W., Rottmann, O., & Kuehn, R. (2004). Mitochondrial DNA phylogeography of red deer (*Cervus elaphus*). *Molecular Phylogenetics and Evolution*, 31(3), 1064-1083. https://doi.org/10.1016/j.ympev.2003.10.003

Mackiewicz, P., Matosiuk, M., Świslocka, M., Zachos, F. E., Hajji, G. M., Saveljev, A. P., Seryodkin, I. V., Farahvash, T., Rezaei, H. R., Torshizi, R. V., Mattioli, S., & Ratkiewicz, M. (2022). Phylogeny and evolution of the genus *Cervus* (*Cervidae*, Mammalia) as revealed by complete mitochondrial genomes. *Scientific Reports*, 12(1), 16381. https://doi.org/10.1038/s41598-022-20763-x

Mahmut, H., Masuda, R., Onuma, M., Takahashi, M., Nagata, J., Suzuki, M., & Ohtaishi, N. (2002). Molecular phylogeography of the red deer (*Cervus elaphus*) populations in Xinjiang of China:
Comparison with other Asian, European, and north American populations. *Zoological Science*, 19(4), 485-495. https://doi.org/10.2108/zsj.19.485

Manichaikul, A., Mychaleckyj, J. C., Rich, S. S., Daly, K., Sale, M., & Chen, W.-M. (2010). Robust relationship inference in genome-wide association studies. *Bioinformatics*, 26(22), 2867-2873. https://doi.org/10.1093/bioinformatics/btq559

Markov, G. G., Kuznetsova, M. V., Danilkin, A. A., Kholodova, M. V., Sugár, L., & Heltai, M. (2015). Genetic diversity of the red deer (*Cervus elaphus* L.) in Hungary revealed by cytochrome b gene. *Acta Zoologica Bulgarica*, 67(1), 11-17.

Marková, S., Horníková, M., Lanier, H. C., Henttonen, H., Searle, J. B., Weider, L. J., & Kotlík, P. (2020). High genomic diversity in the bank vole at the northern apex of a range expansion: The role of multiple colonizations and end-glacial refugia. *Molecular Ecology*, 29(9), 1730-1744. https://doi.org/10.1111/mec.15427

Mattioli, S., Zachos, F. E., Rossi, L., Lister, A. M., & Corlatti, L. (2022). Red deer *Cervus elaphus* Linnaeus. 1758. In L. Corlatti & F. E. Zachos, (Eds.), *Terrestrial cetartiodactyla* (pp. 51-86). Springer. https://doi.org/10.1007/978-3-030-24475-0_19

McDevitt, A. D., Coscia, I., Browett, S. S., Ruiz-González, A., Statham, M. J., Ruczyńska, I., Roberts, L., Stojak, J., Frantz, A. C., Norén, K., Ågren, E. O., Learmount, J., Basto, M., Fernandes, C., Stuart, P., Tosh, D. G., Sindicic, M., Andreanszky, T., Isomursu, M., ... Wójcik, J. M. (2022). Next-generation phylogeography resolves post-glacial colonization patterns in a widespread carnivore, the red fox (*Vulpes vulpes*), in Europe. *Molecular Ecology*, 31(3), 993-1006. https://doi.org/10.1111/mec.16276

McFarlane, S. E., Hunter, D. C., Senn, H. V., Smith, S. L., Holland, R., Huisman, J., & Pemberton, J. M. (2020). Increased genetic marker density reveals high levels of admixture between red deer and introduced Japanese sika in Kintyre, Scotland. *Evolutionary Applications*, 13(2), 432-441. https://doi.org/10.1111/eva.12880

Meiri, M., Kosintsev, P., Conroy, K., Meiri, S., Barnes, I., & Lister, A. (2018). Subspecies dynamics in space and time: A study of the red deer complex using ancient and modern DNA and morphology. *Journal of Biogeography*, 45(2), 367-380. https://doi.org/10.1111/jbi.13124

Meiri, M., Lister, A. M., Higham, T. F. G., Stewart, J. R., Straus, L. G., Obermaier, H., González Morales, M. R., Marín-Arroyo, A. B., & Barnes, I. (2013). Late-glacial recolonization and phylogeography of European red deer (*Cervus elaphus* L.). *Molecular Ecology*, 22(18), 4711-4722. https://doi.org/10.1111/mec.12420

Mingarro, M., & Lobo, J. M. (2021). Connecting protected areas in the Iberian peninsula to facilitate climate change tracking. *Environmental Conservation*, 48(3), 182-191. https://doi.org/10.1017/S037689292100014X

Nei, M. (1972). Genetic distance between populations. *The American Naturalist*, 106, 283-292.

Niedziałkowska, M., Doan, K., Górny, M., Sykut, M., Stefaniak, K., Piotrowska, N., Jędrzejewska, B., Ridush, B., Pawełczyk, S., Mackiewicz, P., Schmölcke, U., Kosintsev, P., Makowiecki, D., Charniauski, M., Krasnodębski, D., Rannamäe, E., Saarma, U., Arakelyan, M., Manaseryan, N., ... Stanković, A. (2021). Winter temperature and forest cover have shaped red deer distribution in Europe and the Ural Mountains since the Late Pleistocene. *Journal of Biogeography*, 48, 147-159. https://doi.org/10.1111/jbi.13989

Niedziałkowska, M., Hundertmark, K. J., Jędrzejewska, B., Sidorovich, V. E., Zalewska, H., Veeroja, R., Solberg, E. J., Laaksonen, S., Sand, H., Solovyev, V. A., Sagaydak, A., Tiainen, J., Juškaitis, R., Done, G., Borodulin, V. A., Tulandin, E. A., & Niedziałkowski, K. (2016). The contemporary genetic pattern of European moose is shaped by postglacial recolonization, bottlenecks, and the geographical barrier of the Baltic Sea. *Biological Journal of the Linnean Society*, 117(4), 879-894. https://doi.org/10.1111/bij.12713

Niedziałkowska, M., Jędrzejewska, B., Honnen, A.-C., Otto, T., Sidorovich, V. E., Perzanowski, K., Skog, A., Hartl, G. B., Borowik, T., Bunevich, A. N., Lang, J., & Zachos, F. E. (2011). Molecular biogeography of red deer *Cervus elaphus* from eastern Europe: Insights from mitochondrial DNA sequences. *Acta Theriologica*, 56(1), 1-12. https://doi.org/10.1007/s13364-010-0002-0

Niedziałkowska, M., Jędrzejewska, B., Wójcik, J. M., & Goodman, S. J. (2012). Genetic structure of red deer population in northeastern Poland in relation to the history of human interventions. *The Journal of Wildlife Management*, 76(6), 1264-1276. https://doi.org/10.1002/jwmg.367

Niedziałkowska, M., Tarnowska, E., Babik, W., Konczal, M., Gharbi, K., Cezard, T., & Jędrzejewska, B. (2023). Different waves of postglacial recolonisation and genomic structure of bank vole populations in NE Poland. *Heredity*, 130(5), 269-277. https://doi.org/10.1038/s41437-023-00600-1

Niedziałkowska, M., Tarnowska, E., Ligmanowska, J., Jędrzejewska, B., Podgórski, T., Radziszewska, A., Ratajczyk, I., Kusza, S., Bunevich, A. N., Danila, G., Shkvyria, M., Grzybowski, T., & Woźniak, M. (2021). Clear phylogeographic pattern and genetic structure of wild boar *Sus scrofa* population in Central and Eastern Europe. *Scientific Reports*, 11(1), 9680. https://doi.org/10.1038/s41598-021-88991-1

Novo, I., Ordás, P., Moraga, N., Santiago, E., Quesada, H., & Caballero, A. (2023). Impact of population structure in the estimation of recent historical effective population size by the software GONE. *Genetics Selection Evolution*, 55(1), 86. https://doi.org/10.1186/s12711-023-00859-2

Nussey, D. H., Pemberton, J., Donald, A., & Kruuk, L. E. B. (2006). Genetic consequences of human management in an introduced island population of red deer (*Cervus elaphus*). *Heredity*, 97(1), 56-65. https://doi.org/10.1038/sj.hdy.6800838

Paradis, E., & Schliep, K. (2019). Ape 5.0: An environment for modern phylogenetics and evolutionary analyses in R. *Bioinformatics*, 35(3), 526-528. https://doi.org/10.1093/bioinformatics/bty633

Passilongo, D., Reby, D., Carranza, J., & Apollonio, M. (2013). Roaring high and low: Composition and possible functions of the Iberian Stag's vocal repertoire. *PLoS One*, 8(5), e63841. https://doi.org/10.1371/journal.pone.0063841

Pemberton, J., Johnston, S. E., & Fletcher, T. J. (2021). The genome sequence of the red deer, *Cervus elaphus* Linnaeus 1758. *Wellcome Open Research*, 6, 336. https://doi.org/10.12688/wellcomeoprenres.17493.1

Pemberton, J., Kruuk, L. E. B., & Clutton-Brock, T. (2022). The unusual value of long-term studies of individuals: The example of the Isle of Rum red deer project. *Annual Review of Ecology, Evolution, and Systematics*, 53(1), 327-351. https://doi.org/10.1146/annurev-econsys-012722-024041

Pemberton, L. W., Cogan, N. O. I., & Forster, J. W. (2013). St AMPP: An R package for calculation of genetic differentiation and structure of mixed-ploidy level populations. *Molecular Ecology Resources*, 13(5), 946-952. https://doi.org/10.1111/1755-0998.12129

Pérez-González, J., Frantz, A. C., Torres-Porras, J., Castillo, L., & Carranza, J. (2012). Population structure, habitat features and genetic structure of managed red deer populations. *European Journal of Wildlife Research*, 58(6), 933-943. https://doi.org/10.1007/s10344-012-0636-0

Pérez-González, J., Gort-Esteve, A., Ruiz-Olmo, J., Anaya, G., Broggini, C., Millán, M. F., Vedel, G., de la Peña, E., Membrillo, A., Seoane, J. M., Azorit, C., & Carranza, J. (2023). Red deer in the Pyrenees: A risky secondary contact zone for conservation genetics. *The Journal of Wildlife Management*, 87(6), e22454. https://doi.org/10.1002/jwmg.22454

Pérez-Suárez, G., Palacios, F., & Boursot, P. (1994). Speciation and phyly in western mediterranean hares (*Lepus castroviejoi*, *L. euro-paeus*, *L. granatensis*, and *L. capensis*) revealed by mitochondrial DNA phylogeny. *Biochemical Genetics*, 32(11-12), 423-436. https://doi.org/10.1007/BF00566063

Peters, L., Huisman, J., Kruuk, L. E., Pemberton, J. M., & Johnston, S. E. (2022). Genomic analysis reveals a polygenic architecture of antler